Comparative genomic analysis of Mycobacterium avium subspecies obtained from multiple host species

Validation of microarray sensitivity and specificity

The reference isolate MAP K-10 was compared to the sequenced MAH isolate 104 in order to evaluate the performance of the oligonucleotide microarray as a platform for comparative genomic hybridizations. The hybridizations and data analyses were performed as described in the Materials and Methods, and the resulting intensity ratios for each oligonucleotide probe were matched to the BLASTn score of the probe when searched against the complete genomes of MAP K-10 and MAH 104. Probe sequences with an E-value of greater than 1 × 10^-5 were identified as not present in the reference genome. The distribution of log transformed hybridization intensity ratios (MAH 104/MAP K-10) is presented in Figure 1. A majority of the probes (75%) were observed to have ratios falling between -0.5 and 0.5, as would be expected for the genomes of these closely related mycobacteria. The distribution of MAP K-10 specific probes peaked below a log ratio of -1.0, with a majority of the probes having ratios lower than this. Mirroring this trend, the MAH 104 specific probes peaked above a log ratio of 1.0 with most of the probes observed to have ratios greater than this. Only five of the MAP K-10 specific probes had positive log ratios, while only three of the MAH 104 specific probes had negative log ratios. In both cases, the log ratios of these outlying probes fell within the range of -0.5 to 0.5. Based on these reference hybridizations, we predict that probes on the MAP K-10 microarray with log ratios greater than 1.2 can be identified as absent from the MAP K-10 genome, while probes with log ratios less than -1.2 are missing from the test isolates relative to MAP K-10. Utilizing these cutoffs, both the false positive and false negative rates are less than 1%.

Sheep genomotype

A similar genomic profile was observed for MAP sheep isolates 397, 467, 7565, and 9040. These sheep isolates were distinguished from the other MAP isolates examined in this study by the presence of four clusters of ORFs homologous to sequences found in the MAH 104 genome but absent in MAP K-10. Additionally, these isolates are missing three clusters of ORFs that are present in the other MAP isolates examined. The sheep isolate 9040 also appears to be a member of this genomotype, although the observed differences in hybridization intensity compared to MAP K-10 are not as pronounced as the other isolates (Figure 2).

Several new or only partially described LSPs common to MAP sheep and MAH isolates were identified in the present study. One of the genomic regions included ORFs MAH2772 – 2791 (Figure 2 and Table 1). This gene cluster was partially identified as the 548 bp sequence PIG-RDA 30 (AY266302) and was mapped to a 27 Kbp region on the MAH 104 genome [9]. This cluster of genes encodes several proteins that may contribute to pathogenesis. MAH2777 is predicted to encode a cytochrome P450 enzyme. Members of this enzyme family have been implicated in basic cellular processes as well as virulence (reviewed in [10]). MAH2782 is homologous to an arylsulfatase characterized in Pseudomonas aeruginosa and like other sulfatases may modulate mycobacterial interactions with host cells (reviewed in [11]). This gene cluster is also characterized by the presence of several genes encoding putative proteins involved in lipid and energy metabolism. Seven of the identified genes encode proteins with no predicted function.

Another novel LSP found in MAP sheep and MAH isolates but absent from MAP cattle isolates was comprised of ORFs MAH3041 – 3052 (3,399,659 – 3,411,697) which are predicted to encode proteins involved in the biosynthesis of glycopeptidolipids [12, 13] (Table 1). Sequencing of this region in the sheep isolate MAP 397 revealed the presence of additional ORFs (hyp, hlpA, gdhgA and gtfC) with homology to glycopeptidolipid biosynthesis genes immediately downstream of MAH3052 (gtfB in Figure 3). These additional ORFs were not homologous to any MAH 104 sequences (and thus were not detected by the microarray analysis) but rather were homologous to MAA TMC724 (ATCC 25291), a serovar 2 isolate of MAA. A direct comparison of the glycopeptidolipid biosynthesis gene clusters from MAP 397 and MAA TMC724 (accession number AF125999) revealed that they are approximately 99% identical with the exception of a deletion at the 5' end of mtfA as well as the absence of the insertion sequence IS2534 from the MAP 397 locus (Figure 3). Despite the presence of these genes, no glycopeptidolipids have been observed to be produced by sheep isolates of MAP under standard in vitro culture conditions (T. Eckstein, personal communication); therefore, the practical contribution of glycopeptidolipids to the pathogenesis of MAP sheep isolates remains unclear at this point.

A putative transcriptional regular labeled as MAH315 in MAH 104 was identified by microarray hybridization as present in MAP sheep isolates. This ORF was subsequently confirmed as present in the MAP 397 genome by sequencing. The protein encoded by MAH315 has homology to the GntR-family of transcriptional regulators which are widely distributed across bacterial species and regulate a variety of cellular processes [14, 15].

A novel LSP comprised of the ORFs MAP1432 – 1438c was identified in the present study as absent from sheep MAP isolates (Table 2). This gene cluster is predicted to encode four energy metabolism enzymes as well as a lipase (MAP1438c). MAP1432 encodes a hypothetical protein with homology to the REP13E12 family of repetitive elements that were originally described in Mycobacterium tuberculosis and have been shown to be targets of phage integration [16]. MAP2656 was initially identified as absent via microarray analysis but sequencing of MAP 397 identified a homologue with 100% identity. This represents the only observed discrepancy between the microarray and sequencing results.

Avium-silvaticum genomotype

The genome profile of MAP sheep isolate 5001 was unlike the other isolates obtained from sheep that were examined in this study. This isolate grouped together with four independent laboratory stocks of what has historically been referred to as M. avium strain 18 [17–22] as well as five isolates of MAS (Figure 2). Additionally, two isolates of MAA (ATCC25291 and ATCC35719) were observed to have genomic profiles similar to this group. The hybridization pattern of these isolates is collectively referred to here as the avium-silvaticum genomotype.

Similar to the MAP sheep isolates, the avium-silvaticum genomotype isolates contained LSPs that are absent from MAP cattle isolates but present in MAH 104. They possessed the MAH3041-3052 and MAH4893-4910 gene clusters also found in MAP sheep isolates; however, they did not possess the MAH1834-1858 or MAH2772-2791 clusters of MAH 104 ORFs (Figure 2 and Table 1). These isolates contained additional LSPs homologous to the MAH 104 genome. A cluster of MAH 104 ORFs including MAH4657-4665 (5,122,587–5,131,447) was present in the avium-silvaticum genomotype isolates. The ORFs in this LSP are predicted to encode proteins that function as transcriptional regulators and dehydrogenases, while several have no known function. A second LSP with homology to the MAH 104 genome included ORFs MAH4704-4709 (5,173,587–5,179,816) (Table 1). With the exception of two ORFs predicted to encode components of a restriction modification system, the remaining coding sequences present in this LSP have no predicted function.

All of the avium-silvaticum genomotype isolates were also observed to have greater amounts of probe hybridized to the oligonucleotide target representing MAP3814c relative to MAP K-10 (Figure 2). This ORF is 31% identical to transposase sequences and may represent part of an insertion sequence that is present in multiple copies within MAS and MAA genomes. The MAP K-10 ORFs immediately downstream of MAP3814c (MAP3815-3818) are missing in both the avium-silvaticum and hominissuis genomotype isolates based on the microarray hybridization results. Notably, ORFs MAP3815-3817 lack homology to any other sequences deposited in Genbank. This appears to provide an additional example of mycobacterial genome diversity mediated by a mobile genetic element.

The avium-silvaticum genomotype isolates could be further separated into two subgroups based on their genomic hybridization profiles. One group included the isolate 5001 and the isolates originally identified as MAS, while the other consisted of the MAA and strain 18 isolates, which the exception of isolate 6071. While both groups of isolates displayed very similar hybridization profiles, they could be distinguished by the polymorphic regions MAP1424c-1464 and MAH2829-2830. The consistent genomic profile amongst the strain 18 isolates obtained from independent laboratories suggests that these isolates have not experienced significant genomic polymorphisms despite several decades of laboratory use. Additionally, these results indicate that the microarray platform used in this study provides reproducible results when comparing isolates.

Hominissuis genomotype

Three of the MAP isolates examined (7265, 7301, and 9034) had hybridization profiles similar to MAH 104 and are collectively referred to as the hominissuis genomotype. Also included with these isolates were two mycobacterial isolates (5835cc and 5836cc) from a duck, which were distinguished from each other upon culture by pigmentation. This group is characterized by the absence of several genes and LSPs that are present in MAP, as well as the presence of additional LSPs not found in MAP. Many of these LSPs are shared with the MAS and MAA isolates examined in this study. The sequenced MAP K-10 and MAH 104 genomes have been previously compared with both bioinformatic and experimental approaches, including DNA microarrays [4, 25–28]. The LSP that includes ORFs MAP1230-1237c appears to be a distinguishing factor of the hominissuis genotype as these ORFs were absent from this genomotype but present in all of the other isolates examined in this study (Table 1). This LSP corresponds to the region previously identified by Tizard and coworkers as a low-GC genetic island present only in MAP and MAS [29]. Isolates 5835cc and 5836cc were distinguished from the other hominissuis genomotype isolates by several polymorphic regions. MAP0072c-0076 is absent from both isolates and includes several genes that are predicted to encode membrane proteins, while the deletion of MAP2523c-2532 includes several enzymes involved in energy metabolism (Figure 2). The hybridization profiles of the remaining hominissuis genomotype isolates closely mirrored MAH 104.

Cattle genomotype

The remaining MAP isolates examined in this study displayed hybridization profiles similar to the sequenced cattle isolate MAP K-10 and are collectively referred to as the cattle genomotype. The distinguishing feature of this genomotype is the large variety of host species represented. This group includes isolates obtained from cattle (n = 5), birds (n = 4), goats (n = 3), bison (n = 3), humans (n = 9), a cat, and an armadillo (Figure 2). Relative to the MAP K-10 genome, no LSPs were observed in any of the cattle genomotype isolates.

SSR Genotyping

PCR amplification of two of the most discriminatory short sequence repeat (SSR) loci – mononucleotide G and tri-nucleotide GGT repeats [40] was carried out as described (2, 30) using primer sets (i) 5'-TCA GAC TGT GCG GTA TGG AA-3' and 5'-GTG TTC GGC AAA GTC GTT GT-3', and (ii) 5'-AGA TGT CGA CCA TCC TGA CC-3' and 5'-AAG TAG GCG TAA CCC CGT TC-3', respectively. The PCR products were purified with a QIAquick PCR purification kit (Qiagen Inc., Valencia, CA) and sequenced by using standard dye terminator chemistry, and the sequences were analyzed on an automated DNA sequencer (3700 DNA Analyzer, Applied Biosystems, Foster City, CA). The alleles were assigned a number congruent to the number of G and GGT residues.

Polymerase chain reaction

All primers were designed with Primer3 software [41]. PCR reactions used genomic DNA as a starting template and used standard conditions reported in detail elsewhere [4].

Microarray Design

The software program ArrayOligoSelector [42] was used to identify 70 mer oligonucleotides specific for every predicted open reading frame (ORF) in the MAP K-10 genome [39]. One 70 mer was designed for each MAP K-10 ORF with a total length of less than 4000 bp, while longer ORFs were split in half and one 70 mer was designed for each half. One 70 mer was also designed for every MAP intergenic region greater than 500 bp. Additionally, an automated annotation of the unannotated MAH 104 genome sequence (downloaded on 03/25/05) was performed using the methods described by McHardy and coworkers [43] and one 70 mer was designed for every predicted ORF that was less then 30% identical to MAP K-10 sequences as determined by BLAST analysis [44]. The annotated MAH 104 ORFs were numbered sequentially starting at the origin of replication, although for reference the nucleotide start and stop positions for each ORF or cluster are also reported. The arbitrarily assigned ORF designations were subsequently matched to the ORF designations used when the completed MAH 104 genome sequence was deposited in Genbank (accession number NC_008595) [see Additional file 1]. A detailed description of this microarray platform has been deposited in the National Center for Biotechnology Information (NCBI) Gene Expression Omnibus under the accession number GPL3433.

Bacterial growth and DNA extraction

Microarray hybridization and data analysis

Microarray hybridization and preliminary data filtering was performed as described in the National Center for Biotechnology Information (NCBI) Gene Expression Omnibus (GEO) database entries (see below). The software program Cluster 3.0 [47] was used to filter out ORFs that were lacking values for greater than 80% of the MAP isolates examined and median center the results for each isolate. Hierarchical clustering of the results was performed and visualized with Cluster 3.0 and Java TreeView [48] or MultiExperiment Viewer [49]. Putative large sequence polymorphisms (LSPs) were identified by a visual examination of the hybridization results and selected LSPs were subsequently confirmed by genome sequencing and/or PCR. All microarray data has been deposited in the NCBI Gene Expression Omnibus database as series number GSE7622 which includes sample accession numbers GSM172202 – GSM172344.

Genome sequencing

Purified genomic DNA from MAP isolate 397 was prepared as described above for microarray hybridizations. Initial shotgun genome sequencing was performed by 454 Life Sciences (Branford, CT). Additional sequences were generated by directly sequencing PCR products and genomic DNA using an ABI3100 capillary sequencer (Applied Biosystems). Sequence data was generated as part of a genome sequencing project and will be published upon completion (M. Paustian, unpublished data).