Introduction

In contrast to the ongoing decline in the rate of tuberculosis, Mycobacterium avium–intracellulare complex pulmonary disease (MAC-PD) is increasing in many parts of the world1,2. The standard therapy for MAC-PD involves multidrug chemotherapy that includes clarithromycin and several antituberculous drugs. Unfortunately, treatment failure is common, resulting in relapse, disease progression, and death3. Therefore, it is critical to investigate the bacterial physiological characteristics that underlie the ability of MAC-PD strains to establish and maintain chronic drug-recalcitrant infections.

Recent advances in comparative genomics have revealed major differences in the genomic features of reference type strains and circulating pathogenic clinical isolates. Extensive genomic diversity has been recognized across clinical isolates of nontuberculous mycobacteria (NTM). Mycobacterium avium has been divided into four subspecies including subspecies, avium, hominissuis, silvaticum, and paratuberculosis 4. Whereas M. intracellulare has been divided into two subgroups, typical M. intracellulare (TMI) and M. paraintracellulare–M. indicus pranii (MP-MIP), with currently two additional distinguished subspecies, yongonense, and chimaera5,6. Despite the increasing number of sequenced MAC-PD strains, the molecular genetic mechanisms that underlie their virulence and pathogenesis are still poorly understood.

The burgeoning field of functional genomics has allowed the rapid assessment of molecular aspects of organisms on a genome-wide scale. Transposon (Tn) sequencing (TnSeq) is a functional genomic approach that combines saturation-level transposon-insertion mutagenesis with next-generation sequencing (NGS) to comprehensively evaluate the fitness costs associated with gene inactivation in bacteria. Over the last decade, TnSeq has been applied to numerous bacterial species using various approaches7. Whereas initial studies applied TnSeq to reference strains, recent approaches have included the comparative assessment of representative clinical strains to better understand the fundamental aspects of pathogen biology8. We have previously used TnSeq to globally characterize the genes that are essential for growth versus those that are essential for hypoxic pellicle formation in the M. intracellulare type strain ATCC139509. However, given the recent expansion of MAC-PD and the well-recognized high levels of genomic diversity in this complex of bacteria, data obtained from the study of a decades-old reference strain is not expected to represent the diversity of extant clinical isolates. ATCC13950 was originally isolated from the abdominal lymph node of the 34-month-old female10. Even though the etiological bacterial species are the same, the clinical manifestation of infant lymphadenitis differs markedly from MAC-PD, which often occurs in middle-aged and elderly female patients with no predisposing immunological disorder. Therefore, it is essential to investigate a collection of strains that captures the diversity of recent clinical MAC-PD strains in order to define key molecular aspects of their pathogenesis and virulence.

In this study, we used TnSeq to analyze a set of recently isolated clinical M. intracellulare strains with diverse genotypes. We integrated the gene essentiality data for these strains and the type strain to identify the common essential and growth-defect-associated genes that represent the genomic diversity of this group of pathogens. We also identified the greater essentiality of genes involved in gluconeogenesis, the type VII secretion system, and cysteine desulfurase in the clinical MAC-PD strains than in the type strain. Furthermore, the requirement for these genes was confirmed in a mouse lung infection model. The profiles of gene essentiality in clinical MAC-PD strains suggest the mechanism of hypoxic adaptation in NTM in terms of promising drug targets, especially for strains that cause clinical MAC-PD.

Results

Common essential and growth-defect-associated genes representing the genomic diversity of M. intracellulare strains

To obtain comprehensive information on the genome-wide gene essentiality of M. intracellulare, we included nine representative strains of M. intracellulare with diverse genotypes in this study (Fig. 1a, b). We used TRANSIT11 to obtain 2–6 million reads of Tn insertions. The average number of Tn insertion reads was > 50 at Tn-inserted TA sites, which confirmed the TnSeq data (Supplementary Table 1). The numbers of essential genes in the clinical MAC-PD strains, calculated with the hidden Markov model (HMM; a transition probability model), tended to be smaller than that in ATCC13950 (Supplementary Table 2). In total, 131 genes were identified as essential or growth-defect-associated with the HMM analysis across all M. intracellulare strains, including the eight MAC-PD clinical strains and ATCC13950 (Fig. 2a, Supplementary Table 3). The genes identified as universal essential or growth-defect-associated were involved in fundamental functions, such as glyoxylate metabolism, purine/pyrimidine metabolism, amino acid synthesis, lipid biosynthesis, arabinogalactan/peptidoglycan biosynthesis, porphyrin metabolism, DNA replication, transcription, amino acid tRNA ligases, ribosomal proteins, general secretion proteins, components of iron–sulfur cluster assembly, some ABC transporters, and the type VII secretion system. The genes identified corresponded to the genes that encode conventional and under-development drug targets, such as those encoding gyrase, arabinosyltransferase, enoyl-(acyl-carrier-protein) reductase, alanine racemase, DNA-directed RNA polymerase, malate synthase12 and F0F1 ATP synthase (Fig. 2b). These data suggest that the set of universal essential or growth-defect-associated genes identified with the HMM analysis of M. intracellulare strains with diverse genotypes has potential utility as antibiotic drug targets.

Phylogenetic tree of the M. intracellulare strains and strategy of the experiment in this study.

a Phylogenetic tree of the M. intracellulare strains used in this study. The tree was generated based on average nucleotide identity (ANI) with the neighbor-joining method. TMI: typical M. intracellulare; MP-MIP: M. paraintracellulare–M. indicus pranii. b Strategy of the experimental procedure.

Identification of the essential and growth-defect-associated genes across all nine M. intracellulare strains used in this study.

a Functional categories of 131 genes identified as universal essential or growth-defect-associated with an HMM analysis. The genes were categorized according to the Kyoto Encyclopedia of Genes and Genomes (KEGG) database. b Universal essential or growth-defect-associated genes corresponding to the genes of existing antituberculous drug targets. c The number of the essential or growth-defect-associated genes in the accessory genome, and the number of strain-dependent essential or growth-defect-associated genes in the core genome.

By comparative genomic analysis, we have revealed 3,153 core genes among 55 M. intracellulare strains including those used in this study5. Of these core genes, 439 genes were classified as strain-dependent essential or growth-defect-associated with the HMM analysis (Supplementary Table 4). And among 5,824 accessory genes, 222 genes were classified as essential or growth-defect-associated with the HMM analysis. The proportion of these strain-dependent and accessory essential or growth-defect-associated genes was comparable to the similar analysis of Streptococcus pneumoniae13. Of note, 17 (3.8%) of strain-dependent essential or growth-defect-associated genes and 14 (6.3%) of accessory essential or growth-defect-associated genes were also included in the 175 genes required for hypoxic pellicle formation in ATCC13950 reported in our previous study9 (Supplementary Table 4).

Genes with increased gene essentiality in clinical MAC-PD strains are also required for hypoxic pellicle formation in the type strain

The sharing of strain-dependent and accessory essential or growth-defect-associated genes with genes required for hypoxic pellicle formation in ATCC13950 prompts us to consider the possibility that the clinical MAC-PD strains increase the gene essentiality of the genes involved in hypoxic growth compared to ATCC13950. To clarify this, we compared the gene essentiality of the clinical MAC-PD strains with that in ATCC13950 by integrating the HMM analysis with the resampling analysis, a method of calculating the number of significantly altered Tn insertion reads (a gene-based permutation model). With this integrating analysis, 121 genes in the clinical MAC-PD strains were shown to have significantly fewer Tn insertions than ATCC13950 (Fig. 3, Supplementary Tables 5, 6). Of these, nine genes were identified in 3–5 clinical MAC-PD strains. These included genes encoding phosphatases (OCU_RS30185, OCU_RS35795), cysteine desulfurase (csd [OCU_RS40305]), methyl transferase (OCU_RS41540), a hypothetical protein (OCU_RS47290), type VII secretion components (mycP5 [OCU_RS38275], eccC5 [OCU_RS38345]), a glycine cleavage system aminomethyltransferase (gcvT [OCU_RS35955]), and phosphoenolpyruvate carboxykinase (pckA [OCU_RS48660]) (Figs 3, 4). The genes involved in gluconeogenesis (such as pckA), the glycine cleavage system, and the type VII secretion system were also identified as required for hypoxic pellicle formation in ATCC13950 in our previous study9. The Tn insertion reads were significantly reduced in the pathways of gluconeogenesis, including in fructose-1,6-bisphosphatase (glpX), phosphoenolpyruvate carboxykinase (pckA), and pyruvate carboxylase (pca), but were significantly increased in those of glycolysis, including in pyruvate kinase (pykA) and pyruvate dehydrogenase (aceE, lpdC), suggesting that carbohydrate metabolism was metabolically remodeled in the clinical MAC-PD strains, which mimics the metabolism observed during hypoxic pellicle formation in ATCC13950 (Fig. 5).

Detection of genes showing increased or reduced gene essentiality in the clinical M. intracellulare strains.

Left panel shows the genes identified as having fewer Tn insertion reads than the type strain ATCC13950. The fold changes in the number of Tn insertion reads calculated by a resampling analysis are represented by the color scale. Red squares indicate genes required for hypoxic pellicle formation in ATCC139509. *Genes identified with a combination of resampling and HMM analyses. †Genes identified only with resampling analysis. Right panel shows the genes identified as having more Tn insertion reads than type strain ATCC13950. NH: no homolog with ATCC13950, NS: no significant increase in Tn insertion reads by a resampling analysis (adjusted P < 0.05).

Overview of the differences in gene essentiality between the clinical M. intracellulare strains and ATCC13950, drawn with Cytoscape55.

Each central node represents the functional category of genes, assigned with a KEGG pathway analysis. Each peripheral node represents the genes showing significant changes in the number of Tn insertion reads. The size of each node represents the number of strains identified. Each edge represents a strain identified as showing significant changes in the number of Tn insertion reads. The thickness of each edge represents the adjusted P value. The color scale of each edge represents the fold change in the number of Tn insertion reads calculated by a resampling analysis. The graphical organization was referenced to the previous publication by Carey, et al8.

Preferential gene requirements for gluconeogenesis in the clinical M. intracellulare strains inferred from the TnSeq results.

a Carbohydrate pathway showing the changes in Tn insertion reads. Genes identified with the resampling analysis are framed by squares. The essentiality of gluconeogenesis-related genes (pckA, glpX) was higher and that of glycolysis-related genes (aceE, lpdC) was lower in the clinical M. intracellulare strains compared with those in ATCC13950. b Data on gene essentiality calculated with a resampling analysis in the clinical M. intracellulare strains compared with ATCC13950. Red indicates a reduction in Tn insertion reads in the clinical M. intracellulare strains. Blue indicates an increase of Tn insertion reads in the clinical M. intracellulare strains. Asterisks indicate statistical significance of log2 fold changes (log2FCs). The graphical organization was referenced to the previous publication by Carey, et al8.

On the other hand, 162 genes were detected with significantly more Tn insertions in the clinical MAC-PD strains than in ATCC13950 (Fig. 3). Of these, 12 genes were identified in more than seven strains. The increased Tn insertion reads in the genes encoding superoxide dismutase (OCU_RS25945) and pyruvate dehydrogenase (aceE [OCU_RS35720]) in the clinical MAC-PD strains suggest the inessentiality of the detoxification of oxygen radicals and glycolysis metabolism in these strains (Figs 35). These genes also included those encoding amidohydrolases (OCU_RS34400, OCU_RS34425, OCU_RS34445), some transcription factors (OCU_RS28440, OCU_RS30230), a CoA transferase (OCU_RS34390), and some predicted transporters (OCU_RS28435, OCU_RS39950, OCU_RS39955, OCU_RS44685). Furthermore, the increase in Tn insertion reads in predicted multidrug-resistant genes, such as those encoding the DrrAB-ABC transporter complex (OCU_RS39950, OCU_RS39955)14 and a putative lipoprotein LpqB (OCU_RS44685), encoded by the two-component system mtrAB operon1516, in the clinical MAC-PD strains were found.

Requirement of genes with increased gene essentiality in the clinical MAC-PD strains for mouse lung infection

Hypoxia is the major environmental condition in natural water17, biofilms18, and infections in humans and animals characterized by caseous granuloma lesions19. The impact of hypoxia on mycobacteria under various ecological circumstances implies that the genes with increased essentiality in the clinical MAC-PD strains compared to ATCC13950 and the genes required for hypoxic pellicle formation in ATCC13950 may be required for infection in vivo. To examine whether the set of genes with increased essentiality in the clinical MAC-PD strains also shows increased essentiality during infection in vivo, we performed TnSeq in mice infected intratracheally with the Tn mutant library strains of M.i.198 or M.i.27, which we had demonstrated to be virulent in mice (Supplementary Fig. 1a). The time course of the changes in the bacterial burden showed a pattern similar to those of the wild-type strains M.i.198 and M.i.27, respectively, except that it was not possible to harvest sufficient colonies (as few as 104/mouse) in the few mice infected with the M.i.27 Tn mutant strain in week 8 and week 16 (Supplementary Fig. 1b, Supplementary Table 7). With a resampling analysis, 629 and 512 genes were shown to be required for in vivo infection with transposon mutant library strains M.i.27 and M.i.198, respectively, between day 1 and week 16 of infection (Fig. 6a, Supplementary Fig. 1c, d, Supplementary Tables 8, 9). Of these genes, 72 (11.9%) and 50 (9.76%) in M.i.27 and M.i.198, respectively, were also required for hypoxic pellicle formation in ATCC13950. Forty-one genes also required for hypoxic pellicle formation in ATCC13950 were identified in both M.i.27 and M.i.198, and included genes associated with succinate production (multifunctional 2-oxoglutarate metabolism enzyme) and the glyoxylate cycle (isocitrate lyase), cysteine and methionine synthesis (thiosulfate sulfurtransferase [CysA], 5-methyltetrahydrofolate-homocysteine methyltransferase [MetH]), the serine/threonine-protein kinase PknG system, proteasome subunits (PrcA, PrcB), DNA repair UvrC, secretion protein EspR, cell-wall synthesis proteins (such as phosphatidylinositol [polyprenol-phosphate-mannose-dependent alpha-[1-2]-phosphatidylinositol pentamannoside mannosyltransferase], peptidoglycan (penicillin-binding protein A pbpA), and L,D-transpeptidase 2), putative membrane proteins and transporters (such as the MmpL family, transport accessory protein MmpS3, transport protein MmpL11, vitamin B12 transport ATP-binding protein BacA, ABC transporter substrate-binding protein, and ABC transporter permease), and mammalian cell entry (Mce) proteins (OCU_RS48855–48875) (Supplementary Table 10).

Significance of increased gene essentiality profiles in vitro on the pathogenesis of mouse lung infection.

a Consistency between the genes required for infection in mouse lungs and the 175 genes required for hypoxic pellicle formation in ATCC13950. b Changes in the requirements for genes with increased essentiality observed in the clinical M. intracellulare strains after the infection of mouse lungs. Highlighted genes are also required for hypoxic pellicle formation by ATCC13950. Gene essentiality of the strain itself: gene essentiality compared with that of ATCC13950. Gene requirements in mouse lungs: changes in the gene requirements in the infected mouse lungs compared with those before infection. Gene set enrichment analysis: The genes listed as core enrichment are shown as “Yes” and the genes not listed as core enrichment are shown as “No.”

The genes only identified in M.i.27 that were also required for hypoxic pellicle formation in ATCC13950 included genes associated with gluconeogenesis (pckA, glpX), the glycine cleavage system (gcvT), the type VII secretion system (eccD5, espG5, eccC5 [OCU_RS38280-OCU_RS38285-OCU_RS38345], eccC3 [OCU_RS48145]), and D-inositol 3-phosphate glycosyltransferase (Supplementary Table 10). These three systems, gluconeogenesis, the glycine cleavage system, and the type VII secretion system, also showed increased gene essentiality in 3–5 clinical MAC-PD strains than in ATCC13950 (Figs 3, 6). Although M.i.27 was excluded as a relevant strain showing increased gene essentiality of these three metabolic systems, the increase in gene requirements for infection in vivo suggests a positive relationship between hypoxic growth in vitro and bacterial growth in vivo. In contrast, M.i.198 already showed increased gene essentiality in these three metabolic systems under aerobic conditions in vitro, so the genes associated with these three systems were not identified as having fewer Tn insertion reads during infection in vivo. Genes only identified in M.i.198 and required for hypoxic pellicle formation in ATCC13950 included genes associated with succinate synthesis (2-oxoglutarate oxidoreductase subunit KorA) or encoding the HTH-type transcriptional repressor KstR and menaquinone reductase (Supplementary Table 10). As increased gene essentiality in 3–5 clinical strains (excluding M.i.198) compared with ATCC13950, an increased requirement for the RNA methyltransferase OCU_RS41540 was identified in M.i.198 during infection in vivo. Despite the differences in the profiles of the genes required for in vivo infection between strains, these data suggest that increased gene essentiality for hypoxic growth confers advantages for pathogenesis in vivo.

We have performed a statistical enrichment analysis of gene sets by GSEA20 to investigate whether the genes required for growth of M.i.27 and M.i.198 in mouse lungs overlap the gene sets required for hypoxic pellicle formation in ATCC13950 as well as the gene sets showing increased gene essentiality observed in the clinical MAC-PD strains listed in Fig 3. Of these gene sets including a total of 180 genes, 54% (98 genes) and 40% (73 genes) were listed as enriched in genes required for mouse infection in M.i.27 and M.i.198, respectively (P < 0.01 in M.i.27 and P = 0.016 in M.i.198) (Supplementary Fig. 2, Supplementary Table 11). The enriched genes both in M.i.27 and M.i.198 included genes encoding fructose -1,6-bisphosphatase (glpX), a type VII secretion protein (eccC5), glycine cleavage system aminomethyltransferase (gcvT), some phosphatase (OCU_RS30185) and a hypothetical protein (OCU_RS47290) (Fig. 6). The enriched genes only in M.i.27 included genes of phosphoenolpyruvate carboxylinase (pckA), cysteine desulfurase (csd) and the ESX-5 type VII secretion component serine protease (mycP5). The enriched genes only in M.i.198 included genes encoding some phosphatase (OCU_RS35975) and class I SAM-dependent RNA methyltransferase (OCU_RS41540). These data suggest the sharing of the hypoxia-adaptation pathways between in vitro hypoxic pellicle formation and in vivo bacterial growth.

Preferential hypoxic adaptation of clinical MAC-PD strains evaluated with bacterial growth kinetics

The metabolic remodeling, such as the increased gene essentiality of gluconeogenesis and the type VII secretion system, and the overlap of the genes required for mouse lung infection and those required for hypoxic pellicle formation involved by confer these metabolic pathways suggests the preferential adaptation of the clinical MAC-PD strains to hypoxic conditions. To confirm this, we investigated the growth kinetics of each strain under aerobic and hypoxic conditions. In terms of growth kinetics, most of the clinical MAC-PD strains (except M001) entered logarithmic (log) phase earlier than ATCC13950 under 5% oxygen (Fig. 7a, Supplementary Fig. 3). However, under aerobic conditions, the clinical MAC-PD strains entered log phase at a similar time to ATCC13950. The growth rate at mid-log point (midpoint) was significantly reduced under hypoxic conditions in a total of 5 clinical MAC-PD strains (M.i.198, M.i.27, M003, M019 and M021) compared to aerobic conditions (Fig. 7b, Supplementary Fig. 3). The growth rate at midpoint under hypoxic conditions was significantly lower in these 5 clinical MAC-PD strains than in ATCC13950. By contrast, the growth rate at midpoint was significantly increased under hypoxic conditions compared to aerobic conditions in ATCC13950. And there was no significant change in the growth rate between aerobic and hypoxic conditions in the rest 3 clinical MAC-PD strains (M018, M001 and MOTT64). The data of more rapid adaptation to hypoxia followed by reduced growth rate once adapted indicate the greater metabolic advantage on continuous logarithmic bacterial growth under hypoxic conditions in the clinical MAC-PD strains than in ATCC13950.

Comparison of the timing of entry into logarithmic growth and logarithmic growth rate by the clinical M. intracellulare strains and ATCC13950.

a Time at the inflection point (midpoint) on the sigmoid growth curve. *significantly earlier than an aerobic culture of ATCC13950; #significantly earlier than a hypoxic culture of ATCC13950; †significantly later than an aerobic culture of ATCC13950; ‡significantly later than a hypoxic culture of ATCC13950. b Growth rate at midpoint of the growth curve in each strain. #significantly slower than a hypoxic culture of ATCC13950; †significantly slower than an aerobic culture of the corresponding strain; ‡significantly faster than an aerobic culture of ATCC13950. Open bars: aerobic; closed bars: 5% O2. Data are shown as the means ± SD of triplicate experiments. Data from one experiment representative of three independent experiments (N = 3) are shown.

Effects of knockdown of universal and accessory/strain-dependent essential or growth-defect-associate genes in clinical MAC-PD strains

To facilitate the functional analysis of the TnSeq-hit genes, it is necessary to construct the genetically-engineered strains. Still now, there are very few reports of gene engineering of MAC strains and none of them has been accompanied by the follow-up studies2126. Moreover, no attempts have been reported to construct the knockdown strains in MAC strains. With an intention to break through this current situation, we used the pRH2052/pRH2521 clustered regularly interspaced short palindromic repeats interference (CRISPR-i) system to construct knockdown strains where the expression of the target gene is suppressed27. Of the nine M. intracellulare strains analyzed in this study, seven could be transformed with the vectors for single guide RNA (sgRNA) expression and dCas9 expression. However, the other two strains, M001 and MOTT64, could not be transformed with these CRISPR-i vectors. To confirm the gene essentiality detected with the HMM analysis, we evaluated the consequent growth inhibition in the knockdown strains of representative universal essential or growth-defect-associated genes, including glcB, inhA, gyrB, and embB. The knockdown strains showed growth suppression with comparative growth rates of knockdown strains to vector control strains that do not express single sgRNA reaching 10-1 to 10-4 (Fig. 8a).

Validation of the inhibition of bacterial growth caused by the suppression of the expression of TnSeq-hit genes using the CRISPR-i system.

Open bar: day 3; closed bar: day 7. a Comparative growth rates of the knockdown strains relative to those of the vector control strains in the representative universal essential or growth-defect-associated genes: glcB, inhA, gyrB and embB. b Comparative growth rates of the knockdown strains relative to those of the vector control strains in the representative accessory and strain-dependent essential or growth-defect-associated genes: pckA, glpX, csd and ESX-5 type VII secretion components. Data are shown as the means ± SD of triplicate experiments. Data from one experiment representative of two independent experiments (N = 2) are shown.

As representative accessory and strain-dependent essential or growth-defect-associated genes, we evaluated the 9 genes showing reduced Tn insertion reads in the 3– 5 clinical MAC-PD strains as listed in Fig. 3, as well as glpX, a gene of gluconeogenesis where four MAC-PD strains showed reduced Tn insertion reads in a resampling analysis (Fig. 5b). The growth was suppressed in the pckA-knockdown strains in M.i.198 and M003, in the glpX-knockdown strain in M003, and in the cysteine desulfurase gene (csd)-knockdown strains in M.i.198 and M003, with comparative growth rates of knockdown strains to vector control strains reaching 10-2 to 10-4 (Fig. 8b). The growth tended to be suppressed in the type VII secretion protein gene (eccC5)-knockdown strain in M003 and in the type VII secretion-associated serine protease mycosin gene (mycP5)-knockdown strain in M003, with comparative growth rates of knockdown strains to vector control strains reaching less than 10-1. By contrast, the growth was not suppressed in these knockdowns trains in ATCC13950, with comparative growth rates of knockdown strains to vector control strains not reaching 10-1, which was consistent to the inessentiality of these genes under aerobic conditions in ATCC13950. As for M.i.198 and M003, the consistency between TnSeq data and CRISPR-i data was supported by a correlation between the growth suppression of the knockdown strains and the fitness costs evaluated by a resampling analysis (Supplementary Table 12).

In contrast to the cases of M.i.198 and M003, the growth of knockdown strains of the accessory and strain-dependent essential genes was not effectively suppressed in the other MAC-PD strains despite positively detected in TnSeq (Fig. 5b, Supplementary Fig. 4). These data suggest the difficulties in the effective operation of genetically modified system using foreign genes from different bacterial species in MAC-PD strains. Although the suppression of gene expression was confirmed in the knockdown strains of these accessory and strain-dependent essential or growth-defect-associate genes, the gene expression levels were not suppressed to < 1% as reported in M. tuberculosis (Mtb)28 (Supplementary Fig. 5). The differences between strains as well as between genes in the degree of suppression of gene expression required to inhibit bacterial growth may explain such inconsistencies between the TnSeq and CRISPR-i data.

To be summarized, we could validate the universal essential or growth-defect-associated genes by the CRISPR-i system in the M. intracellulare strains used in this study. And, although not successful with all strains, we could validate the accessory and strain-dependent essential or growth-defect-associated genes, especially genes of gluconeogenesis, type VII secretion components and iron-sulfur cluster assembly in several MAC-PD strains like M.i.198 and M003.

Discussion

In this study, using TnSeq of nine genetically diverse M. intracellulare strains, we have identified 131 shared essential and growth-defect-associated genes with an HMM analysis. The genes identified are promising drug targets, representing the genomic diversity of the clinical MAC-PD strains and reflecting the actual clinical situation. The strategy of narrowing the essential genes of bacteria is a useful approach for refining the database of drug targets. We have also demonstrated an overview of the differences in gene essentiality required for hypoxic growth between the ATCC type strain and the clinical M. intracellulare strains responsible for MAC-PD, with the evidence of earlier entry to log phase followed by slower growth rate observed in most of these clinical strains compared to ATCC13950. This is the first study to demonstrate the positive relationship between the profiles of increased gene essentiality required for hypoxic metabolic remodeling and the phenotypes advantageous for bacterial survival under hypoxic conditions.

In this study, 121 genes were detected as significantly reduced Tn insertion in the clinical MAC-PD strains compared to ATCC13950, which was higher than the case of Mtb showing 10–50 genes with altered Tn insertion reads detected by a resampling analysis8. This may reflect the genomic diversity of M. intracellulare compared with that of Mtb. The proportion of accessory genes is larger in NTM than in Mtb29. In a clinical trial of phage therapy for M. abscessus infection, the number of accessory genes was shown to decrease with the loss of virulence in M. abscessus during the course of therapy30, suggesting some role for accessory genes for its clinical pathogenesis. Our data of the overlap of the accessory essential or growth-defect-associated genes such as pckA and mycP5 (Supplementary Table 5) with genes required for mouse lung infection (Fig. 6b, Supplementary Fig. 2) supports the idea that accessory genes (non-core genes) may encode functions required for bacterial survival, especially under nonoptimal conditions, such as hypoxia and during the infection of human and animal bodies.

Our TnSeq data are consistent with metabolomic data that have shown the importance of the gluconeogenic carbon flow of TCA cycle intermediates for bacterial growth under hypoxic conditions in terms of the essentiality of PckA and the reverse TCA cycle3134. Gluconeogenesis is essential for bacterial growth when fatty acids are used as a major nutrient33. PckA catalyzes the conversion of the metabolites derived from the TCA cycle to the metabolites of gluconeogenesis32. The reaction catalyzed by PckA is essential for bacterial growth when the available nutrients are predominantly fatty acids35,36. The main nutrients for mycobacteria inside hypoxic granulomas in an infected human or mouse body are lipids, such as fatty acids and cholesterol, rather than carbohydrates33,37. In the present study, M. intracellulare grew under 5% oxygen, whereas Mtb stops growing under 1% oxygen and enters a dormant state32. The depletion of phosphoenolpyruvate in the pckA deletion mutant of Mtb reduced the phosphoenolpyruvate–carbon flux essential for bacterial growth and drug sensitivity32. The present study supports the idea that PckA is essential for the microaerobic growth of mycobacteria. Our TnSeq data imply the preferential adaptation of the clinical MAC-PD strains for circumstances in vivo, which are characterized by hypoxia.

On the other hand, such metabolic shift towards gluconeogenesis and reverse TCA cycle can also be observed under reactive nitrogen specie and in vivo infection as well. This suggests that the metabolic shift observed under hypoxia is a general stress response rather than being specific to hypoxia due to the nature of the provided/available carbon source38. Several factors other than hypoxia can be estimated to be related to the change of the gene essentiality in the clinical MAC-PD strains such as host-related factors and general stress conditions (i.e. nutrition of abundant fatty acids and less carbohydrates, low pH circumstances)33,3537. Except for our data of the profile of genes required for hypoxic pellicle formation in M. intracellulare9, no fundamental genome-wide data are available on the change of the gene essentiality in response to each stress factors in non-tuberculous mycobacteria. Further studies are needed to discriminate which factors are specifically involved in the increased essentiality of each gene. Although of interest, it is beyond the scope of the current study.

About 10% of the genes required for mouse infection are also required for hypoxic pellicle formation. The genes showing increased gene essentiality in the clinical MAC-PD strains, such as pckA and genes encoding the ESX-5 type VII secretion system, including eccC5 and mycP5, also showed increased gene requirements in mouse infections. In a TnSeq study of mice infected intravenously with M. avium, the genes of the ESX-5 type VII secretion system were required for infection of the liver and spleen39. The ESX-5 type VII secretion system is involved in the transport of fatty acids and the virulence of mycobacteria40,41. These similarities in the profiles of gene essentiality prompt us to propose some role of biofilm formation in mycobacterial pathogenesis in animals and humans. In the present study, we shed light on the impact of hypoxia on biofilm formation, the acquisition of virulence, and the in vivo pathogenesis of NTM in terms of global gene essentiality. Although the morphologies of the pellicle and macroscopic lesions (including cavitation and nodular bronchiectasis) are not similar, monitoring the expression of the genes identified in this study may provide clues to the mechanism of biofilm formation in vivo.

In this study, we revealed the pattern of hypoxic adaptation characteristic to the clinical MAC-PD strains by the growth curve analysis: an early entry to log phase followed by a reduced growth rate (Fig. 7). This pattern implicates the continuous impact on the infected hosts during logarithmic bacterial growth, which may be involved in the persistent and steadily progressive illness for years in humans3. The advantage of slow growth in mycobacteria may be inferred by the role of a histone-like nucleoid protein, Mycobacterial DNA-binding protein 1 (MDP1), also designated HupB4244. The expression of MDP1 is enhanced in the stationary phase, resulting in reduced growth speed and long-term survival with increased resistance to reactive oxygen species (ROS)4244. The reduced growth rate after entry to log phase may allow long-term survival of mycobacteria under ROS-rich circumstances inside host cells. The pattern of hypoxic adaptation as observed in the clinical MAC-PD strains may enhance the pathogenesis by extending the period of the repeated process of bacterial invasion and replication in infected macrophages.

Contrary to our expectation, not all clinical strains showed rapid adaptation from aerobic to hypoxic conditions relative to ATCC13950. These data suggest that some, but not a major proportion of, clinical M. intracellulare strains do not have an advantageous phenotype for continuous logarithmic growth under hypoxia. In Mtb, the adaptation to hypoxia is considered to be necessary for clinical pathogenesis because the conditions in vivo, such as inside macrophages or granulomas, are hypoxic38,45. However, our strain M001 entered log phase later than ATCC13950 under 5% oxygen. Our inability to construct knockdown strains in M001 prevented us from clarifying its clinical pathogenesis without preferential hypoxic adaptation. However, because we have shown that most of the M. intracellulare strains with different genotypes showed preferential adaptation to hypoxia in this study, our results largely reflect the clinical situation of the etiological agents of MAC-PD. A possible strategy for identifying factors that confer clinical pathogenesis will be to focus on the major clinical strains that show preferential hypoxic adaption.

The data of validation of the universal essential or growth-defect associated genes by the CRISPR-i system was overall acceptable in the clinically MAC-PD strains whose knockdown strains could be constructed (Fig. 8a). However, we have encountered the difficulties of validating the TnSeq-hit genes by CRISPR-i experiment, especially accessory and strain-dependent essential or growth-defect-associated genes in the MAC-PD strains other than M.i.198 or M003 (Fig. 8b, Supplementary Fig. 4). As a possible reason, the low efficiency of the CRISPR-i system may obscure the increase in gene essentiality suggested by the TnSeq data (Supplementary Fig. 5). The efficiency of the suppression of gene expression by CRISPR-i in accessory and strain-dependent essential or growth-defect-associated genes was lower in our M. intracellulare strains than in Mtb28. In Mtb, the expression of mmpL3, a gene encoding a mycolic acid transporter, is reported to be suppressed to 6% under 10 ng/ml anhydrotetracycline (aTc) and by < 1% under 100 ng/ml aTc in the knockdown strains, compared to strains expressing scrambled single guide RNAs (sgRNA)28. The gene expression levels were suppressed in our knockdown strains by around 10%–40% under 50–200 ng/ml aTc relative to the levels found in strains not expressing sgRNA (Supplementary Fig. 5). The growth inhibition effect is theoretically dependent upon the degree to which gene expression is suppressed. On the other hand, the inconsistent growth phenotypes of knockdown strains of accessory and strain-dependent essential or growth-defect-associated genes between MAC-PD strains implies some differences in the degree of suppression of gene expression required to inhibit bacterial growth between genes as well as between bacterial strains. Despite the issue of the lower efficiency of the CRISPR-I system in M. intracellulare than in Mtb, this study has opened the avenue for functional analysis of NTM that leads to the development of novel drugs that target universal and accessory/strain-dependent essential metabolic pathways that has significant impact on bacterial growth, especially by making use of M.i.198 and M003 as optimal NTM strains for gene manipulation.

As in the glpX-knockdown strains, csd knockdown caused growth suppression in several MAC-PD strains. This type II Csd is involved in the biosynthesis of the iron– sulfur clusters of the SUF machinery. Compared with the type I cysteine desulfurase IscS, the SUF machinery is active under oxidative-stress and iron-limiting conditions46. Representatives of iron–sulfur cluster proteins include components of the electron-transport system, such as NADH dehydrogenase I (Nuo) and succinate dehydrogenase (Sdh), which are components of respiratory complexes I and II, respectively. Under hypoxic conditions, the non-proton-translocating type II NADH dehydrogenase Ndh is upregulated, and cytochrome bd oxidase (CydAB) and the proton-pumping type I NADH dehydrogenase (Nuo) are downregulated38. In the reverse TCA cycle on the same reaction pathway as Sdh, fumarate reductase (Fnr) is also an iron–sulfur cluster protein that produces succinate, and is a proposed alternative electron acceptor under hypoxic conditions34,47. The increased essentiality of the csd gene observed in the clinical MAC-PD strains seems to reflect the switching of the electron-transport system to less-energy-efficient respiration to support preferential hypoxic survival, including in vivo.

We tried to use inhibitors of phosphoenolpyruvate carboxykinase to confirm the increased essentiality of the pckA gene. However, we detected no positive effect because millimolar 3-mercaptopicolinic acid is required to show the effect of inhibition, even in eukaryotic cells, which have no cell wall48. Moreover, 3-alkyl-1,8-dibenzylxanthine49 also only effectively inhibits PckA in eukaryotic cells. However, PckA has been proposed as a drug target for Mtb50, so our TnSeq data seem to be consistent with studies directed towards antimycobacterial drug discovery.

In summary, we have identified 131 universal essential and growth-defect-associated genes that covers the genomic diversity of M. intracellulare strains responsible for MAC-PD, which narrowed down the drug targets against MAC-PD. Furthermore, we have demonstrated a preferential adaptation to hypoxia in clinical NTM strains compared with the type strain in terms of the increased essentiality of genes required for hypoxic growth, such as those encoding phosphoenolpyruvate carboxykinase, the ESX-5 type VII secretion system, and cysteine desulfurase. These genes were shown to have a significant impact on infection in vivo by increasing their gene requirements. Our data highlight the importance of using clinical strains and hypoxic conditions in experimental research, and provide fundamental resources for developing novel drugs directed towards nontuberculous mycobacterial strains that cause severe clinical diseases.

Methods

Bacterial strains and phages used to prepare Tn mutants

A total of nine strains were used in this study, including seven clinical strains isolated from patients with MAC-PD in Japan (M.i.198, M.i.27, M018, M003, M019, M021, and M001)51, M. intracellulare genovar paraintracellulare type strain MOTT64 isolated from MAC-PD in South Korea52, and the M. intracellulare type strain ATCC13950 (Fig. 1). The Tn mutant libraries were constructed using glucose as a carbon source under aerobic conditions according to the transduction method used for ATCC139509, except that a higher concentration of kanamycin (100 μg/ml) was used to select the transposon mutants than was used for ATCC13950 (25 μg/ml). In brief, a hundred milliliter of the bacteria were cultured aerobically in 7H9/ADC broth until mid-log phase (optical density at a wavelength of 600 nm [OD600] of approximately 0.4–0.6). After harvesting the bacteria, 10 mL of high-titer mycobacteriophage phAE18053 was transduced into the bacteria 1-day at 37°C for inserting Tn and plated on 7H10/OADC agar plates containing kanamycin. After 2 or 3 week-cultivation at 37°C aerobically, the colonies were harvested en masse and used as Tn mutant libraries. The kanamycin-resistant colonies (> 1 × 105) were harvested and evenly resuspended in 7H9/ADC medium containing 20% glycerol, aliquoted, and stored at −80 °C until further use.

Next-generation sequencing and TnSeq analysis

DNA libraries were prepared as reported previously54. The TnSeq libraries were sequenced with the NextSeq 500 System, 150-bp single-end run (Illumina, San Diego, CA, USA). The sequence reads were analyzed as previously described9. Gene essentiality was calculated with resampling and HMM analyses using TRANSIT11. First, an HMM analysis was conducted to calculate the fluctuation of the number of Tn insertion reads individually for each strain. Next, a resampling analysis was conducted to calculate the difference of the number of Tn insertion reads between each clinical strain and ATCC13950 (Supplementary Table 6). By comparing the data obtained by HMM and resampling analyses, we summarized the difference of gene essentiality between the clinical MAC-PD strains and ATCC13950. To compare the gene essentiality of the clinical MAC-PD strains and ATCC13950, the annotation for the clinical MAC-PD strains was adapted from that of ATCC13950 by adjusting the START and END coordinates of each open reading frame (ORF) in the clinical MAC-PD strains according to their alignment with the corresponding ORFs of ATCC13950. The function of TnSeq-hit genes was categorized according to the Kyoto Encyclopedia of Genes and Genomes (KEGG) database.

Comparison of growth curves under aerobic and hypoxic conditions

First, we precultured the bacteria in 7H9/ADC broth until the bacteria reached early log-phase (OD600 approximately 0.1–0.3). Next, the bacterial cultures were diluted with new 7H9/ADC broth (final OD600 0.003) for the aerobic and hypoxic growth assay. Bacterial cells were grown statically in 96-well plates containing 250 μl of 7H9/ADC broth in each well under aerobic conditions or 5% oxygen at 37 °C with humidification (APM-30D incubator, Astec, Tokyo, Japan). The culture was sampled every 1–3 days. At sampling, the culture was mixed by pipetting 50 times. An aliquot (4 μl) of the culture and its diluents were streaked onto 7H10/OADC agar and the CFUs counted. The growth curve was fitted to the logistic 4P model with the JMP software (SAS Institute Inc., Cary, NC, USA) to calculate the inflection point (midpoint) and growth rate at midpoint in the log growth phase. Entry into the log phase was evaluated from the timing of the inflection point on the growth curve.

Animals

Female C57BL/6JJcl mice were purchased from CLEA Japan (Tokyo, Japan). Six-week-old mice were used for the infection experiment (Supplementary Fig. 1). All the mice were kept under specific-pathogen-free conditions in the animal facility of Niigata University Graduate School of Medicine, according to the institutional guidelines for animal experiments. The mice were housed (four per cage) under a 12-h/12-h light/dark cycle. Food and water were available ad libitum. Kanamycin was added to the drinking water to a final concentration of 100 μg/ml to maintain the mice’s resistance to the Tn mutant strains. Euthanasia was conducted with an intraperitoneal injection of a combined anesthetic (15 μg of medetomidine, 80 μg of midazolam, and 100 μg of butorphanol), followed by cervical dislocation.

Mouse TnSeq experiment

The protocol was based on our previous mouse lung infection experiment51. Twenty mice each were used for the Tn mutant strains M.i.198 and M.i.27. The mice were injected intratracheally with 3.25 × 108 and 5.93 × 107 CFUs of the Tn mutant strains, respectively. The infected lungs were expected to be harvested from five mice at day 1, week 4, week 8, and week 16 after infection. However, the lungs were harvested from only four mice infected with the M.i.27 Tn mutant strains in week 8 because we experienced difficulties with sampling, and the lungs were harvested from only three, one, and two mice infected with M.i.198 Tn mutant strains in week 8, week 12, and week 16, respectively, because morbidity was high in the mice infected with M.i.198. The harvested lungs were homogenized in 4.5 ml of distilled water with the gentleMACS™ Tissue Dissociator (Miltenyi Biotec, Bergish Gladbach, Germany) and the homogenates were plated on Middlebrook 7H10–oleic acid ADC (7H10-OADC) agar plates and cultured for 3 weeks. The colonies were harvested en masse for TnSeq.

Gene set enrichment analysis

To investigate the overlap of the gene sets required for growth in mouse lungs with those required for hypoxic pellicle formation observed in ATCC13950 and those showing increased gene essentiality observed in the clinical MAC-PD strains, a statistical enrichment of gene sets was performed by GSEA 4.3.220. Parameters set for GSEA were as follows: permutations = 1,000, permutation type: gene_set, enrichment statistic: weighted, metric for ranking genes: Signal2Noise, max size: 500, min size: 15.

Construction and growth assay of knockdown strains

Knockdown strains were constructed with the CRISPR-dCas9-mediated RNA interference system, as reported previously27. Gene-specific 20-nt sequences were selected as protospacers for incorporation into the sgRNAs (Supplementary Table 13). pRH2502 (a vector expressing an inactivated version of the cas9 DNA sequence) and pRH2521 (an sgRNA-expressing vector whose expression is regulated by the TetR-regulated smyc promotor Pmyc1tetO) were introduced into the bacteria to obtain colonies of the knockdown strains. For selection, 25 μg/ml kanamycin and 50 μg/ml hygromycin were used for ATCC13950, and 100 μg/ml kanamycin and 50 μg/ml hygromycin were used for the clinical strains. The growth assay of the knockdown strains was performed in 7H9/ADC broth containing the above-mentioned concentrations of kanamycin and hygromycin under aerobic conditions at 37 °C with humidification (APM-30D incubator, Astec, Tokyo, Japan). During the assay, the cultures were supplemented with anhydrotetracycline (aTc) to a final concentration of 50 ng/ml (as a concentration that is not toxic to ATCC13950 and M.i.27) or 200 ng/ml (for the other clinical strains). aTc was added repeatedly every 48 h to maintain the induction of dCas9 and sgRNAs in experiments that extended beyond 48 h27. An aliquot (4 μl) of the culture and its diluents were streaked onto 7H10/OADC agar plates containing kanamycin and hygromycin, and the CFUs were counted. Growth inhibition was evaluated as the ratio of the CFUs of knockdown strains compared with those of the vector control strains lacking the guide RNA sequence. To confirm the downregulation of gene expression in the knockdown strains, quantitative reverse transcription PCR (qRT-PCR) was performed. Bacteria were cultured in 10 ml of Middlebrook 7H9 medium (Difco Laboratories, Detroit, MI, USA) containing 0.2% glycerol, 0.1% Tween 80 (MP Biochemicals, Illkrich, France), and 10% Middlebrook albumin–dextrose–catalase (7H9/ADC/Tween 80), kanamycin, and hygromycin until mid-log phase (OD600 approximately 0.4–0.6). The cultures were supplemented with aTc every 48 h, as in the growth assay for the knockdown strains. One day after the second addition of aTc, RNA was extracted by bead-beating six times at 5000 rpm for 30 s with cooling (Micro Smash MS 100R Cell Disruptor; Tomy Seiko, Tokyo, Japan) and cDNA was synthesized with the Direct-zol™ RNA Miniprep (Zymo Research, Orange, CA, USA) and ReverTra Ace qPCR RT Master Mix with gDNA Remover (Toyobo, Osaka, Japan). qRT-PCR was performed with the CFX Connect Real-time PCR Detection System (Bio-Rad, Hercules, CA, USA), according to the manufacturer’s instructions. The data were normalized to the expression of sigA.

Statistical analyses

Pearson’s correlation test was performed between the fitness costs (log2FC calculated by the resampling assay), efficiency of knockdown (assessed by qRT-PCR) and growth suppression in the knockdown strains. Statistical analyses of the growth curves were performed using a 4-parameter logistic model (logistic 4P model), Wilcoxon two-sample test, and a one-way analysis of variance (ANOVA) followed by Steel’s multiple comparisons test using the JMP software program (SAS Institute Inc., Cary, NC, USA). Data are presented as the mean ± SD. Statistical significance was set as P values less than 0.05.

Acknowledgements

The authors thank Dr. Torin Weisbrod and Dr. William R. Jacobs, Jr. (Albert Einstein College of Medicine, Bronx, NY, USA) for kindly providing mycobacterial phage phAE180, and Dr. Norikazu Hara and Dr. Takeshi Ikeuchi (Department of Molecular Genetics, Bioresource Science, Brain Center for Bioresources, Brain Research Institute, Niigata University) for next-generation sequencing with the NextSeq 500 System. This work was supported by Grants-in-Aid for Scientific Research (grant numbers 20KK0216, 21K19637, and 22H03115 to Yoshitaka Tateishi) from the Ministry of Health, Labour, and Welfare, of Japan This work was also supported by grants from the Japanese Ministry of Education, Culture, Sports, Science, and Technology, the Ministry of Health, and the Research Program on Emerging Infectious Disease from the Japan Agency for Medical Research and Development (23fk0108673h0501).

Author contributions

Y.T. conceived, designed and performed the experiments. Y.O., A.N., and S.M. contributed reagents or materials. Y. Morishige., Y. Minato, and A.D.B. contributed analytical tools. Y. Morishige., Y. Minato, A.D.B., and S.M. provided critical suggestions throughout the preparation of the manuscript. Y.T. wrote the manuscript. All authors read and approved the final manuscript.

Availability of data and materials

The datasets used or analyzed in this study are available from the corresponding author upon reasonable request.

Declarations

Ethics approval and consent to participate

All experimental protocols were approved by the Institutional Review Board of Niigata University (approval numbers: 2019-0020 for the use of human tissue samples and SA00506 for experiments on live vertebrates). All experiments on humans and human samples were conducted in accordance with relevant guidelines and regulations. Informed consent was obtained from all patients. All experiments on animals were conducted in accordance with the relevant guidelines and regulations, and the study complied with the ARRIVE guidelines.

Consent for publication

According to the protocols approved by the Institutional Review Board of Niigata University mentioned above, consent for publication was obtained from the patients enrolled in this study.

Competing interests

The authors declare that they have no competing interests.