{'fraction_non_alphanumeric': 0.037931498855321634, 'fraction_numerical': 0.034923912555550565, 'mean_word_length': 4.971053336906995, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
A proportion of mutations fixed in the genomes of in vitro selected isogenic drug-resistant Mycobacterium tuberculosis mutants can be detected as minority variants in the parent culture

We studied genomic variation in a previously selected collection of isogenic Mycobacterium tuberculosis laboratory strains subjected to one or two rounds of antibiotic selection. Whole genome sequencing analysis identified eleven single, unique mutations (four synonymous, six non-synonymous, one intergenic), in addition to drug resistance-conferring mutations, that were fixed in the genomes of six monoresistant strains. Eight loci, present as minority variants (five non-synonymous, three synonymous) in the genome of the susceptible parent strain, became fixed in the genomes of multiple daughter strains. None of these mutations are known to be involved with drug resistance. Our results confirm previously observed genomic stability for M. tuberculosis , although the parent strain had accumulated allelic variants at multiple locations in an antibiotic-free in vitro environment. It is therefore likely to assume that these so-called hitchhiking mutations were co-selected and fixed in multiple daughter strains during antibiotic selection. The presence of multiple allelic variations, accumulated under non-selective conditions, which become fixed during subsequent selective steps, deserves attention. The wider availability of ‘deep’ sequencing methods could help to detect multiple bacterial (sub)populations within patients with high resolution and would therefore be useful in assisting in the detailed investigation of transmission chains.


INTRODUCTION
Drug-resistant Mycobacterium tuberculosis, the causative agent of tuberculosis (TB), is an increasing problem that is severely complicating effective treatment of the otherwise curable disease. Multidrug-resistant TB (MDR-TB) is caused by M. tuberculosis bacilli that are by definition resistant to first-line drugs rifampicin and isoniazid. Resistance in M. tuberculosis is predominantly conferred by genomic mutations and MDR strains typically carry resistance-mediating mutations in rpoB, conferring resistance to rifampicin, and in katG, conferring resistance to isoniazid (Musser 1995;Ramaswamy and Musser 1998). Resistance mutations can impair bacterial fitness which may affect survival within the host and onward transmission (Cohen, Sommers and Murray 2003;Mariam et al., 2004;Gagneux et al., 2006;Luciani et al., 2009). Some evidence is available that resistant mutants can restore this fitness loss by the acquisition of adaptive mutations (Sherman et al., 1996;Casali et al., 2012Casali et al., , 2014Comas et al., 2012).
Mycobacterium tuberculosis is a monomorphic organism with no evidence of horizontal transfer of genetic elements and a low recombination rate (Sreevatsan et al., 1997;Brosch et al., 2002;Smith et al., 2006Smith et al., , 2009Gagneux and Small 2007;Achtman 2008); within the members of the M. tuberculosis complex (M. tuberculosis, M. bovis, M. africanum, M. microti and M. pinnipedii), there is a genomic resemblance of >99%. The population structure of M. tuberculosis furthermore suggests that phylogenetic lineages, which can be distinguished on the basis of few genetic markers (Coll et al., 2014), have co-evolved with and adapted to their respective hosts (Caws et al., 2008).
Typing of M. tuberculosis strains is commonly achieved by one or a combination of the following methods: RFLP, CRISPR typing, VNTR/MLVA, MLST or SNP typing. The latter has been gaining popularity, as it is particularly appropriate for the classification of a monomorphic organism such as M. tuberculosis, is less prone to homoplasies and is compatible with other sequencebased tools, such as whole genome or deep sequencing. With the cost of next-generation sequencing methods dropping rapidly and their availability increasing, many researchers have used these tools to study the micro-evolution of M. tuberculosis isolates. Genomic analysis of M. tuberculosis isolates recovered from patients has shown that, even over longer periods of time, the genome of the pathogen remains relatively conserved (Schürch et al., 2010). Likewise, after transmission or selection for antibiotic resistance, minimal genetic variability has been detected (Niemann et al., 2009;Ford et al., 2011Ford et al., , 2012Saunders et al., 2011;Casali et al., 2012Casali et al., , 2014Comas et al., 2012;Bryant et al., 2013;Clark et al., 2013;), although in some cases an (initial) increased level of genetic divergence has been observed for drug resistance mutations (Mariam et al., 2011).
The evolution of M. tuberculosis is probably most accurately elucidated by genomic analysis of sequential isolates from a single patient or a transmission chain. However, as there may be multiple mycobacterial foci within a single patient, each being subjected to severe sequential bottlenecks, multiple mutants may emerge simultaneously within the patient. Thus, sequencing single isolates may dramatically overestimate the relative prevalence and likelihood of particular mutations and underestimate the dynamics of traits under evolution or very strong selection, such as drug resistance. A recent study focusing on the emergence of quinolone resistance within the host has indeed shown that multiple resistance-associated mutations can occur inside one patient and that the distribution of these mutants can vary significantly between sequential isolates (Streicher et al., 2012).
We investigated the in vitro evolution of a collection of wellcharacterized isogenic drug-resistant mutants after antibiotic selection on solid medium by whole genome sequencing. Our in vitro findings corroborate data from previous studies where clinical isolates under evolution show a similar degree of genomic stability (Schürch et al., 2010;Ford et al., 2011;Mariam et al., 2011;Bryant et al., 2013;Clark et al., 2013). Additionally, in the major-ity of the strains studied we also found evidence of hitchhiking mutations that were present as minority variants in the parent strain, as well as fixation of additional mutations not associated with drug resistance. On the basis of these signature mutations a phylogenetic relationship could be inferred between all strains included in this study. We believe that this illustrates the great potential of whole genome sequencing analysis, in particular deep sequencing identifying minority variants, to aid investigations on patient-to-patient transmission of M. tuberculosis, particularly when clinical or epidemiological data is missing or incomplete.

Strains
Drug-resistant mutants used in this study were all derived from a cultivar of the M. tuberculosis laboratory strain MTB72 (ATCC 35801) which has been propagated in our laboratory for more than 20 years and is susceptible to all commonly used anti-TB drugs. All strains are derived from experiments that have been previously described (Anthony et al., 2005;Bergval et al., 2007;Bergval et al., 2009).  Fig. 1 and Bergval et al. (2009)] from a liquid culture from MTB72, by plating an aliquot from a log-phase MTB72 culture (ca. 3 weeks, 10 8 CFU) onto antibiotic-containing solid medium. Strains R4, R46, R181 and R190 were selected with 8 μg mL −1 rifampicin, strains R181 and R190 in a different experiment and with a different starting culture of MTB72, but in the same month. Strain B15 was selected with 8 μg mL −1 rifabutin and strain B66 with 0.8 μg mL −1 rifabutin in two different experiments using two starting cultures from MTB72 over a period of two months. The rifampicin MIC and the rifabutin MIC was previously determined for all rifampicin-resistant strains (Anthony et al., 2005); R190 showed resistant to rifampicin (32 μg mL −1 ), but susceptible to rifabutin (<0.8 μg mL −1 ), which allowed us to subject strain R190 to a second round of selection ( Fig. 1). Strains RB14, RB15, RB16 and RB19 are spontaneous, double-resistant, mutants derived from strain R190 selected by plating an aliquot from an R190 liquid culture onto solid medium containing 2 μg mL −1 rifabutin.

Culture and selection of mutants
Strains H15 (March 2006), H71 (August 2006) and H103 (November 2006) are selected by 20 μg mL −1 (H15 and H71) and 0.4 μg mL −1 (H103) isoniazid; in the case of H71, the parent culture was treated with hydrogen peroxide before an aliquot was plated on selective medium (Bergval et al., 2009). The strains were selected in three different experiments, using three different starting cultures of MTB72.

Genetic analysis
For all mutants, the mutation responsible for resistance to drug(s) in question was determined via Sanger sequencing (Anthony et al., 2005;Bergval et al., 2009); two regions of rpoB were analysed for strains resistant to rifampicin and/or rifabutin for the isoniazid-resistant strains parts of the mabA-inhA intergenic region and of katG, known to contain mutants in clinical isolates, were analysed. Before sequencing the associated genes, the isoniazid-resistant strains were first assessed by an adapted multiplex ligation-dependent probe amplification  assay, to determine if a deletion in or of katG was present (Bergval et al., 2009). The genomes of strains MTB72 and R190 were previously analysed by comparative genome hybridization (Bergval et al., 2009).

DNA isolation and sequencing
For each strain, DNA was extracted via the cetyltrimethylammonium bromide procedure (Soolingen et al., 1991) from 1 mL of a logarithmic phase liquid culture (ca. 3 weeks growth) (Anthony et al., 2005;Bergval et al., 2009). The DNA concentration of the samples used in this study was between 53 and 142 ng mL DNA, as determined with the NanoDrop (Thermo Fisher Scientific Inc., Wilmington, Delaware, USA). The sequencing was carried out at the KAUST Genomics Facility, using Illumina HiSeq 2000 technology with paired-end reads of length 100 base pairs.

Data analysis
Trimmomatic software (Bolger, Lohse and Usadel 2014) was used to remove low-quality reads and trim low-quality 3 ends of reads. Nucleotide positions in the raw reads with a quality score lower than Q20 were removed. High-quality reads were then mapped to the H37Rv reference genome (Camus et al., 2002) (Genbank accession: AL123456.3) using the BWA-mem algorithm (Li 2014). SAMtools (Li et al., 2009) andGATK (DePristo et al., 2011) were used to call SNPs. SNP sites were then selected as the intersection dataset between those obtained from both programs. Polymorphic sites at non-unique regions of the genome (Derrien et al., 2012) were filtered out. The depth of coverage was calculated from alignment files using a custom Perl script. For each sample, SNP genotypes were called if supported by at least 20-fold total coverage, otherwise they were classified as missing. The read depth and coverage of all variants identified in all strains included in this study are depicted in Fig. S1 (Supporting Information). Multiple genotypes were allowed by adding alleles supported by at least 5% of total number of reads. We excluded positions in the genome with more than 20% missing genotypes across all samples. In the final dataset, we retained 912 SNP sites, 83% (912/1096) of the original SNP sites derived using SAMtools and GATK. Artemis (Rutherford et al., 2000) was used to verify these variants in the BAM files and to determine the read depth. All sequences have been deposited to Sequence Read Archive (http://trace.ncbi.nlm.nih.gov/Traces/sra/), NCBI's repository for high-throughput sequencing data and are publicly available under accession number SRP049364.

RESULTS
Our collection of isogenic mutants included a pansusceptible parent strain (MTB72), nine first-generation (single-step) drugresistant daughter strains and four second-generation (twostep) daughter strains (Fig. 1). The nine monoresistant strains that were directly derived from MTB72 by a single round of antibiotic selection (R-, B-and H-series, Fig. 1) all acquired at least one additional mutation in a gene not directly associated with drug resistance (Table 1).
A total of 21 non-drug resistance mutations were detected, eight of which were identified in multiple daughter strains. Interestingly, these eight alleles present in multiple daughter strains could already be detected in the genome of the parent strain (indicated as mixed alleles in Table 1). The mutations that could be detected in the parent strain show various linkage patterns: mutations in ctpV, Rv1215c, ceoB and Rv2693c are shared by strains R4, R46, R181, R190 (and offspring) and B66; mutations in dnaB, ispE and gcvB are shared by strains R46 and B66; the mutation in ilvA is shared by H15, H71 and H103 (Table 1).
Eleven mutations were uniquely found in six strains and could not be detected in the genome of the parent strain: R4 (hycQ-T138T, zwf1-R136C, murF-A501V), R181 (lppY-G97R), R190 (dnaA-V270I, Rv0466-S93S, shared with offspring: zmp1-D637D, Rv1215c-M243I), B15 (Rv2974c-L251L), H103 Table 1. Mutations and codon changes identified in strains analysed in this study. Depicted are only those loci for which at least one of the strains included had acquired a mutation. The genomes of all strains were mapped against the H37Rv reference genome and the genomic location of each mutation therefore refers to this genome. Drug resistance-conferring mutations are depicted in Fig. 1 and are not included in this figure. Sequences for strain MTB72, the susceptible parent strain, are depicted in bold. If multiple alleles were observed above the expected error rate, the sequence of both alleles is given, with the dominant allele first. The sequence in the daughter strains was regarded as a mutant if the sequence differed from the sequence of the dominant allele.
: non-synonymous mutation, a : the sequence of the dominant allele is different from the H37Rv sequence, whereas the sequence of the minority variant is similar to that of H37Rv, b : the mutation was acquired in the parent strain. (wag22-Rv1760-A(-)63G) and H71 (Rv3446c-V241I) ( Table 1). These singletons were fixed in the genomes of the strains carrying them; no other alleles could be detected [see Table S1 (Supporting Information) for coverage and read depth of all variants]. Six out of eleven mutations were non-synonymous (highlighted in blue, Table 1), four were synonymous and one was intergenic. The mutations in crp and pks2 observed in H15 (Table 1) were accompanied by the wild-type allele which indicates that these variants were acquired in H15 post antibiotic selection. Both mutations are non-synonymous and the percentages of both mutant alleles were higher than the corresponding wild-type allele; the allelic distributions detected in H15 were 14/0/0/52 and 0/40/0/9 versus 111/0/0/0 and 0/0/0/59 in MTB72 for crp and pks2, respectively (order of the bases A/C/G/T).
One of the first-generation mutants, R190, was subjected to a second round of selection with rifabutin and four of its daughter strains (RB14, RB15, RB16 and RB19), with double rpoB mutations, were also included in this study (Fig. 1). Strain R190 acquired, in addition to the rpoB-S522L mutation associated with rifampicin resistance, four unique mutations (in dnaA, zmp1, Rv1215c and Rv0466, Table 1), two of which were absent in its offspring (dnaA-V270I and Rv0466-S93S), illustrating that these mutations were acquired and/or selected for after the daughter strains were selected. No additional mutations apart from those associated with rifabutin resistance were identified in the four daughter strains of R190.
Based on the data acquired by whole genome analysis, we inferred a phylogenetic tree (Fig. 2), using the data on drug resistance mutations to resolve the branches in the cases where this was necessary (e.g. the rpoB-S531L mutation which was acquired separately by four different strains). This tree is highly similar to the evolution assumed on the basis of our experimental procedure (Fig. 1), with the exception of strain B66 that showed to be more related to the R strains (namely R46, Fig. 2) than we anticipated. Most branches are represented by two or more mutational events, since it was not possible to determine the exact order and timing for each mutation separately.

DISCUSSION
The population structure of M. tuberculosis is highly clonal, with all phylogenetic members sharing more than 99% of their genome sequence. Genomic studies have found that epidemiologically linked strains, extracted from (patient-to-patient) transmission chains (Schürch et al., 2010;Bryant et al., 2013;Clark et al., 2013;Torok et al., 2013) or within the patient (Saunders et al., 2011;Sun et al., 2012), are genetically stable with estimated mutation rates of 0.3-0.5 SNPs/year in low-endemic areas Walker et al., 2013).
We have analysed the genomes of an isogenic collection of in vitro selected drug-resistant strains to determine if mutations, other than those conferring drug resistance, had accumulated during selection or cultivation. Our results show that during selection, mutations not associated with drug resistance were co-selected with the resistance trait. Some of these hitchhiking mutations were already detectable in the genome of the parent strain as minority variants and were likely accumulated over time in only a discrete amount of subpopulations [Tables 1 and 2 and Fig. S2 (Supporting Information)] (Barrick and Lenski 2013).
At the loci that showed mixed alleles in the MTB72 genome (Table 1), we assumed that the dominant allele was the original sequence and the minority variant represented the mutation (Table 2). However, for ctpV, Rv1674c, ceoB and Rv2693c, the dominant allele was different to that of the reference genome from H37Rv; instead, the minority allele corresponded to this reference genome. It is possible that the mutant allele had become dominant in the MTB72 culture in these specific locations. Conversely, the mutations at these loci in the MTB72 may have been historical and the minority variants would therefore represent 'back mutations'. From the data we have available, it was not possible to determine which scenario had occurred.
None of the strains acquired a mutation in rpoA or rpoC, some of which have been suggested to compensate for the fitness loss inflicted by rpoB mutations, particularly rpoB-S531L (Casali et al., 2012;Comas et al., 2012). However, compared to H37Rv, MTB72, and its entire offspring, carry an rpoC-G594E mutation. This mutation is not associated with compensation, but is a genotypic marker, specific for a Euro-American sub-lineage containing X and Haarlem spoligotype families (lineage 4.1 as described in Coll et al. (2014). Although an epistatic effect on the fitness cost of the rifampicin resistance mutations acquired of this rpoC-G594E mutation cannot be ruled out.
We cannot be sure if the singletons observed in the genomes of R4, R181, B15, H103, H71 and R190 (and its offspring) (Table 1) were either co-selected and therefore pre-existing as rare variants in the MTB72 culture, or acquired post-selection. However, the short cultivation period for these strains between selection and sequence analysis makes it probable that these mutations were co-selected.
The two unique mutations that were observed in strain R190 could only have been acquired after selection with rifabutin, since none of the R190 daughter strains carried either of these mutations. The fixation of these two post-selection mutations Table 2. Read depth and coverage of loci identified as mixed alleles in the parent strain. Depicted are only those loci where mutations were acquired by at least one daughter strain and mixed alleles were observed in the parent strain MTB72. The order of the bases is A/C/G/T. Data highlighted in bold indicates that the base was detected in >5% of the reads. dnaB  ctpV  ispE  ilvA  Rv1674c  gcvB  ceoB  Rv2693c  Genomic position  62690  1078856  1130414  1763680  1899754  2077437  3009506  raises the possibility that one or both confer a benefit and were under selective pressure, possibly restoring a fitness deficit resulting from the rpoB-S522L mutation. This appears most likely for the mutation in dnaA (Table 1), which is non-synonymous and thus suggestive of a phenotypic effect. Our findings demonstrate that there can be a detectable amount of minority variants present in a bacterial population, possibly reflecting the genetic diversity in (cultured) patient isolates that are otherwise assumed to be largely homogenous. Studies with Streptococcus pneumoniae also show that bacteria are much more likely than previously anticipated to acquire mutations in addition to those that confer the specific trait selected for (Stevens and Sebert 2011). Furthermore, this concept may have consequences for the researcher aiming to investigate single isogenic mutants to assess their associated (in vitro) characteristics, as we show here that even during a single selection step, unintended genetic changes may be introduced. Especially against a stable genetic background such as that of M. tuberculosis these hitchhiking mutations may be of greater influence on the phenotype than anticipated.

Gene/CDS
Sequencing technologies that produce greater numbers of higher quality reads of longer length increase the possibility to detect minority variants with low copy numbers, allowing the true genetic dynamics of a bacterial population under evolution to be studied. During patient-to-patient transmission, rigorous bottlenecks can dramatically reduce the effective, net mutation rate by selecting out most of the newly acquired mutations within the bacterial populations transmitted (Sun et al., 2012). Deep sequencing (> 500-fold coverage) therefore has the potential to be particularly appropriate to study the dynamics of M. tuberculosis strain(s) within a single patient.
Here, we have shown that a single selection event can result in the one-step fixation of a potentially informative number of pre-existing mutations. Minority variants in clinical isolates should thus not be ignored, as the potentially characteristic signatures of these hitchhiking or co-selected mutations in the genomes of infected strains can be used to help resolve the details of patient-to-patient transmission.