INTRODUCTION
Tuberculosis (TB) remains the leading cause of single bacterial infectious disease mortality worldwide, and is responsible for at least 1 million deaths each year [1–3]. The COVID-19 pandemic has had devastating effects on global TB control efforts, and the consequences are gradually becoming apparent [1]. In the future, infections due to drug sensitive or resistant mycobacteria are expected to increase. Therefore, an urgent need exists to increase understanding of novel mechanisms of drug action and resistance.
Bedaquiline (BDQ) is a diarylquinoline-based, species-selective chemical inhibitor of the Mycobacterium tuberculosis (Mtb) ATP synthase and the first novel TB drug to be approved by the FDA in 40 years [4, 5]. BDQ exhibits bactericidal activity against Mtb and a series of mycobacteria, such as Mycobacterium bovis, Mycobacterium avium, and Mycobacterium smegmatis [6]. After the clinical introduction of BDQ in 2012, several genetic mutations or resistance-associated variants were identified. For instance, mutations in AtpE (ATP synthase c chain) confer target-based BDQ resistance [6]. Mutations in Rv0678, an MmpR transcriptional repressor of the MmpS5-MmpL5 system, cause efflux-based BDQ resistance [7]. In addition, PepQ (an aminopeptidase) has been linked with low-level BDQ resistance [8].
Mycobacterium marinum (M. marinum), the closest genetic relative of the M. tuberculosis complex, is a non-tuberculous mycobacterium (NTM) that causes a TB-like disease in fish and can infect humans when damaged skin is exposed to a contaminated aqueous environment. This pathogen spreads through the lymphatic system and causes particularly serious symptoms among immunocompromised patients. Fortunately, previous research has indicated that BDQ efficiently kills M. marinum in vitro and therefore may potentially serve as a powerful weapon in the arsenal for fighting this zoonotic pathogen [9, 10]. Here, to better understand the mechanism of action and resistance of BDQ in M. marinum, we characterized and sequenced 58 BDQ-resistant mutants isolated in vitro from an M. marinum M strain, and discovered several new mutations linked with BDQ resistance that have not previously been reported in Mtb or other mycobacteria.
MATERIALS AND METHODS
Bacterial strain and culture conditions
The Mycobacterium marinum M strain, a kind gift from Professors Qian Gao and Chuan Wang (Fudan University), was used in this research. Wild-type bacteria and BDQ-resistant isolates were cultured in BSL2 at 30°C in regular Middlebrook 7H9 broth supplemented with 0.2% glycerol and 10% ADC, or on 7H10 agar supplemented with 0.5% glycerol and 10% OADC.
MIC determination
BDQ was purchased from MedChemExpress Co for MIC testing and isolation of resistant mutants. The wild-type and BDQ-resistant bacteria were first cultured at 30°C for 3–4 days, then adjusted to OD600 of 0.001 as the starting culture to determine MIC values in a 96-well plate. All MIC plates were cultured in aerobic conditions at 30°C for 7 days. The MIC for each bacterium was determined through the regular Alamar Blue assay [11].
Laboratory-based adaptive evolution of BDQ-resistant mutants
A 5 mL volume of M. marinum culture (∼1.3×108 CFU) [12] was inoculated onto ten 7H10 agar plates (0.5 mL each) containing 0.1 μg/mL (1× MIC) BDQ, for the first round of laboratory evolution. After incubation at 30°C for 4 weeks, all colonies appearing on the plates were picked and resuspended in 20 μL 7H9 medium, then inoculated on a 7H10 plate with 2× MIC BDQ. After another 4 weeks, the surviving colonies were transferred to 7H10 plates with a higher concentration BDQ (4× MIC), for selection of highly drug-resistant mutants.
Genomic DNA extraction
The liquid-cultured M. marinum strain was harvested after a 10 min centrifugation at 12,000×g. The genomic DNA of bacteria was extracted with a Wizard® Genomic DNA Purification Kit (Promega) according to the manufacturer’s protocol. Purified genomic DNA was quantified with a TBS-380 fluorometer (Turner BioSystems Inc., Sunnyvale, CA). High-quality DNA (OD260/280 ≥1.5, ≥150 ng) was used in further experiments.
Library construction and genome sequencing
The draft genome sequence analyses of the M. marinum parental strain and resistant isolates were performed on the Illumina NovaSeq6000 sequencing platform (MajorBio Co., Shanghai, China). Briefly, DNA samples were sheared into 400–500 bp fragments with a Covaris M220 Focused Acoustic Shearer according to the manufacturer’s protocol. Illumina sequencing libraries were prepared from the sheared fragments with a NEXTflex™ Rapid DNA-Seq Kit. The 5′ primer ends were first end-repaired and phosphorylated. Next, the 3′ ends were A-tailed and ligated to sequencing adapters. Third, the adapter-ligated products were cloned with PCR. The prepared libraries then were used for paired-end Illumina sequencing (2×150 bp) on the Illumina NovaSeq6000 sequencing platform.
Genome assembly and annotation
The data generated from the Illumina platform were used for bioinformatics analysis. All of analyses were performed with the free online Majorbio Cloud Platform (www.majorbio.com) from Shanghai Majorbio Bio-pharm Technology Co., Ltd. The detailed procedures are as follows. Raw reads obtained after sequencing were filtered with fastp software (version 0.19.6) [13] and assembled with SOPA de novo version 2.04 [14]. For each isolate, 1.09–1.94 gigabase (169.6-fold to 306.04-fold genome coverage) sequences were generated after barcodes were trimmed. A total of 324 to 1414 contigs, and 3653376 to 6480401 raw paired reads were generated; the genome sizes were 4.85–6.94 megabases. The WGS data have been submitted to the Sequence Read Archive of the National Center for Biotechnology Information as fastq files (accession number PRJNA903131). Glimmer [15] was used for CDS prediction, tRNA-scan-SE was used for tRNA prediction, and Barrnap was used for rRNA prediction. The predicted CDSs were annotated from the NR, Swiss-Prot, Pfam, GO, COG and KEGG databases with sequence alignment tools such as BLASTP, Diamond and HMMER. Briefly, each set of query proteins was aligned with the databases, and annotations of best-matched hits (e-value < 10–5) were obtained for gene annotation. Mutations in proline-glutamic acid/proline–proline–glutamic acid family genes and in regions with repetitive sequences were excluded from the analysis.
RESULTS
Isolation of M. marinum mutants resistant to BDQ
First, we determined that the MIC of BDQ for the susceptible M. marinum M parent strain was 0.1 μg/mL, which was close to the original reported value. As shown in Fig 1, to conduct the laboratory evolution for generating the BDQ-resistant mutant, 1.3×108 bacteria were plated onto 7H10 plates containing a 1× MIC concentration of BDQ (0.1 μg/mL). In total, 176 isolates grew and formed visible clones after 4 weeks of antibiotic selection (Fig 1). These 176 resistant isolates were transferred to 7H10 plates containing a 2× MIC concentration of BDQ (0.2 μg/mL). Finally, 57 of 176 isolates survived in this increased stress condition and were transferred to 7H10 plates containing a 4× MIC concentration of BDQ (0.4 μg/mL). Eventually, only one BDQ resistant isolate was recovered. The MIC of this isolate was then determined to be 0.6–1.0 μg/mL with an Alamar Blue assay (Fig 2).
Liquid MIC tests of representative bedaquiline-resistant isolates.
WT: wild type. Strain IDs (mutant genes): 002 (MMAR_1007), 016 (YrbE3A-1), 036 (MMAR_1007, YrbE3A-1, MMAR_1793), 059 (MMAR_1007, YrbE3A-1), 067 (MMAR_1007, MMAR_1604), 093 (YrbE3A-1, MMAR_1873), 176 (MMAR_1007, YrbE3A-1, MMAR_1963) and 4×057 (YrbE3A-1, atpB).
Mutations identified in BDQ-resistant mutants by WGS
A total of 58 BDQ-resistant mutants as well as the parental M strain were subjected to WGS. As shown in Table 1, we identified seven different gene mutations in total of 58 mutants. These genomic deep sequencing results were confirmed by PCR amplification and DNA sequencing (Table 2). Notably, the highest level of drug resistance (6–10× MIC) was associated with a mutation in AtpB, the primary biochemical target of BDQ in Mtb, thus indicating that BDQ also primarily blocks the functional role of ATP synthase in M. marinum. In contrast, numerous mutations and insertions mapped to the gene MMAR_1007(46/58), which encodes the homolog of Rv0678 (MmpR) in Mtb. This finding suggested that the mechanism of efflux-based low-level BDQ resistance (against 2× MIC) is largely conserved between Mtb and M. marinum.
Mutations identified in bedaquiline-resistant mutants of M. marinum by WGS analysis.
| Locus tag | Gene product | Nucleotide mutation | Amino acid change | Frequency | Mtb locus tag (gene/product) |
|---|---|---|---|---|---|
| MMAR_4049 | YrbE3A-1 | G563A | G188D | 54/58 | rv1964 (YrbE3A)a |
| MMAR_1007 | Transcriptional regulator | G269C | R90P | 1/58 | rv0678 (MmpR)b |
| MMAR_1007 | Transcriptional regulator | C296T | A99V | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | C394T | R132 stop | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G122A | G41D | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | T298G | F100V | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G404A | R135H | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G370 deletion | 124 codon shift | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | C170A | A57E | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G215T | R72L | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G110A | G37D | 2/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | T131G | L44W | 3/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | T350C | L117P | 2/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | C247T | L83F | 3/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | T449G | L150R | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G310T | E104 stop | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G70T | G24C | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | 289G insertion | 97 codon shift | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | C251T | A84V | 3/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G493 deletion | 165 codon shift | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | A434 deletion | 145 codon shift | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | A202G | S68G | 2/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | 301C insertion | 101 codon shift | 3/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | T299G | F100C | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G203A | S68N | 2/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | A494 deletion | 165 codon shift | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | A488G | E163G | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G74A | G25E | 2/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | C265T | R89W | 2/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G197A | G66E | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | G361A | G121S | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | T2C | start lost | 1/58 | rv0678 (MmpR) |
| MMAR_1007 | Transcriptional regulator | 375T insertion | 126 codon shift | 1/58 | rv0678 (MmpR) |
| MMAR_4093 | ATP synthase A chain AtpB | A511C | I171L | 1/58 | rv1304 (AtpB) |
| MMAR_1963 | Metallo-beta-lactamase superfamily protein | G250T | E84 stop | 1/58 | rv2752c |
| MMAR_1837 | Conserved transmembrane protein | 214ACATCGCCG insertion | 71DIA insertion | 1/58 | rv2869c (Rip1)c |
| MMAR_1793 | Conserved hypothetical protein | G123T | W41C | 1/58 | rv2915c |
| MMAR_1604 | Conserved hypothetical membrane protein | C320T | T107I | 1/58 | - |
a rv1964 encodes integral membrane protein YrbE3A.
b rv0678 encodes the MmpS5–MmpL5 efflux pump repressor MmpR.
c rv2869c gene encodes the regulated intramembrane proteolysis (RIP) metalloprotease Rip1.
Primers used to verify gene mutations in bedaquiline-resistant mutants.
| Gene | Primer | Sequence (5′-3′) |
|---|---|---|
| MMAR_4049 | YrbE3A-1-F | ATGGTTGCTCCTGTTGCCGTGGCCA |
| YrbE3A-1-R | TCACAGCGTCACCTTGATCGCCACC | |
| MMAR_1007 | MMAR-1007-F | ATGAAGTTGTTGTAGACCTATCGGG |
| MMAR-1007-R | TCGACAACATCGGGTCGGGCAACAA | |
| atpB | atpB-F | ATGACTGAGTCGATCCTGGCCGCCC |
| atpB-R | TTAGTGGTGGTCCTCTTCTAGCTCC | |
| MMAR_1963 | MMAR-1963-F | GTGAACGAAGAACTTCCCCCACCAG |
| MMAR-1963-R | TCAAACTTCAATGACCGTGGGCACG | |
| MMAR_1837 | MMAR-1837-F | ATGATGTTCGTTGTCGGCATTGTGC |
| MMAR-1837-R | TTATTGGAACAACCTGATCGGGTTG | |
| MMAR_1793 | MMAR-1793-F | GTGCGACTGCACGTGCGGGGGCGGG |
| MMAR-1793-R | CTAGCGGTGACCTAGTTGGTTGGGG | |
| MMAR_1604 | MMAR-1604-F2 | TCAATACCAGGCAAGCGCCCCAGG |
| MMAR-1604-R2 | CGCTTGGCCAGCCGGTCCATCCCG |
Interestingly, we also discovered a panel of novel mutations that have not previously been reported in Mtb or other mycobacteria. More than 93% of the mutants (54/58) contained a single mutation (G563A) within MMAR_4049, which encodes the integral membrane protein YrbE3A-1. Both YrbE3A-1 and its Mtb homolog YrbE3A (encoded by rv1964) belong to the YrbE family and are annotated as the ABC transporter permease [16]. Unexpectedly, this enzyme is missing in the genomes of M. bovis, M. bovis BCG (Pasteur) strain and M. smegmatis. We next compared the structural similarity between YrbE3A-1 and YrbE3A and found that the mutation site (glycine to aspartate at residue 188) mapped to the linker region between two critical α helixes, thereby suggesting potential changes in structure-based protein activity (Fig 3). Another four mutations were found at a low frequency (1/58): MMAR_1963 (encoding a metallo-beta-lactamase superfamily protein), MMAR_1837 (encoding a conserved transmembrane protein), MMAR_1793 (encoding a conserved hypothetical protein) and MMAR_1604 (encoding a conserved hypothetical membrane protein). Some of these proteins have structurally conserved transmembrane domains and/or predicted enzymatic activity, but their true biological roles and how they contribute to BDQ resistance remain to be explored.
Overall structure of M. marinum YrbE3A_1.
Glycine 188, which was changed to aspartate in bedaquiline-resistant isolates, is shown in detail. Different colors in the 3D model represent the per-residue confidence score (pLDDT), ranging between 0 and 100 (dark blue: Plddt > 90; light blue: 90 > pLDDT > 70; yellow: 70 > pLDDT > 50; and orange: pLDDT < 50). Red dots and blue dots represent oxygen and nitrogen, respectively. Hydrogen bonds are represented by blue dotted lines. The 3D model of YrbE3A_1 came from the AlphaFold protein structure database (https://alphafold.com/entry/B2HQG4).
DISCUSSION
NTM infections are rapidly increasing worldwide. These bacteria are difficult to kill with many routine antibiotics, owing to their intrinsic drug resistance [17]. M. marinum, an opportunistic zoonotic pathogen, has been recognized as a model species for NTM research. Consequently, understanding of its genome, physiology, and pathogenicity has gradually improved [9]. In this research, we combined a laboratory-based adaptive evolution strategy with a WGS approach to identify the spontaneous chromosomal mutations associated with BDQ resistance in M. marinum.
A previous study in Mtb has indicated that mutations in the F0 operon (which includes the atpB, atpE, and atpF genes, whose gene products form the entire membrane-bound F0 unit of the ATP synthase) confer high-level BDQ resistance [18]. Our findings regarding atpB confirmed the above model. This specific mutation (I171L) is located at the α-helices aH5 within the a-subunit of M. marinum ATP synthase (Fig 4). This position is part of the key site for proton translocation [19]. Interestingly, the frequency of mutations in ATP synthase was relatively low with respect to that of efflux-based mutations, in contrast to previous reports in Mtb. These findings may suggest a potential bacterial fitness loss under specific environments, although we did not observe a significant growth defect of this M. marinum strain bearing the above atpB point mutation (data not shown). Future investigations on more clinical Mtb or NTM strains are required to further determine the relationship between atpB mutation and BDQ resistance.
Overall structure of M. marinum AtpB.
Isoleucine 171, which was changed to leucine in bedaquiline-resistant isolates, is shown in detail. The color labeling is as in Fig 3. The 3D model of AtpB came from the AlphaFold protein structure database (https://alphafold.com/entry/B2HQG4).
The unique findings of this research involved the characterization of a panel of efflux-based resistance mutants. These mutants were found at relatively high frequency but conferred only low-level BDQ resistance. In addition, earlier research has suggested that efflux-based mutations may trigger drug cross-resistance. For example, the Mtb carrying mutation on rv0678 shows a resistance phenotype to both BDQ and clofazimine [20]. In this study, we found that these mutations and insertions of MMAR_1007 (the homolog of rv0678) showed a diverse distribution, thus suggesting the genetic and functional vulnerability of this transcriptional repressor modulated MmpS5-MmpL5 system [7, 20]. Most notably, we identified that MMAR_4049 (encoding an ABC transporter permease) is a novel high-frequency mutation locus linked with BDQ resistance. A recent study has indicated that YrbE3A promotes the host innate immune response by targeting NF-κB/JNK signaling in Mtb [16]. However, further evidence of its functional role in bacterial energy metabolism or drug resistance is lacking. To better understand the role of YrbE3A-1 and the other four genes, genetic manipulation, such as knockout, knockdown or overexpression studies, as well as drug cross-resistance assays, will be necessary in the future.
In conclusion, we identified novel mutations associated with BDQ resistance in M. marinum. Our findings show that the mechanisms of BDQ resistance in M. marinum are quite complicated, and include both target-based and efflux-based mediators. Our research thus provides a landscape for better understanding of the molecular basis of BDQ resistance, which should aid in the development of novel potent anti-NTM (BDQ-based) drug regimens and facilitate the discovery of new diagnostic assays capable of detecting BDQ-resistance in M. marinum.