Genome Organization of Sphingobium indicum B90A: An Archetypal Hexachlorocyclohexane (HCH) Degrading Genotype.

Among sphingomonads, Sphingobium indicum B90A is widely investigated for its ability to degrade a manmade pesticide, γ-hexachlorocyclohexane (γ-HCH) and its isomers (α-, β-, δ-, and ε-HCH). In this study, complete genome of strain B90A was constructed using Single Molecule Real Time Sequencing (SMRT) and Illumina platform. The complete genome revealed that strain B90A harbors four replicons: one chromosome (3,654,322 bp) and three plasmids designated as pSRL1 (139,218 bp), pSRL2 (108,430 bp) and pSRL3 (43,761 bp). The study determined the precise location of lin genes (genes associated with the degradation of HCH isomers), for example, linA2, linB, linDER, linF, linGHIJ, and linKLMN on the chromosome; linA1, linC, and linF on pSRL1 and linDEbR on pSRL3. Strain B90A contained 26 copies of IS6100 element and most of them (15 copies) was found to be associated with lin genes. Duplication of several lin genes including linA, linDER, linGHIJ, and linF along with two variants of linE, that is, linEa (hydroquinone 1,2-dioxygenase) and linEb (chlorohydroquinone/hydroquinone 1,2-dioxygenase) were identified. This suggests that strain B90A not only possess efficient machinery for upper and lower HCH degradation pathways but it can also act on both hydroquinone and chlorohydroquinone metabolites produced during γ-HCH degradation. Synteny analysis revealed the duplication and transposition of linA gene (HCH dehydrochlorinase) between the chromosome and pSRL1, possibly through homologous recombination between adjacent IS6100 elements. Further, in silico analysis and laboratory experiments revealed that incomplete tyrosine metabolism was responsible for the production of extracellular brown pigment which distinguished strain B90A from other HCH degrading sphingomonads. The precise localization of lin genes, and transposable elements (IS6100) on different replicons now opens up several experimental avenues to elucidate the functions and regulatory mechanism of lin genes acquisition and transfer that were not completely known among the bacterial population inhabiting the HCH contaminated environment.

Introduction

Sphingobium indicum B90A has been extensively studied for the degradation of hexachlorocyclohexane (HCH) isomers (Sahu et al., 1990; Kumari et al., 2002; Lal et al., 2010). Strain B90A was isolated from HCH treated sugarcane rhizosphere soil, Cuttack, India and found to degrade the most recalcitrant HCH isomer, β-HCH (Sahu et al., 1990). Till date many sphingomonads have been isolated that harbor lin genes (Nagata et al., 1999; Boltner etal. 2005; Kaur et al., 2013; Tabata et al., 2013). Among these, only strain B90A is reported to harbor two copies of linA (HCH dehydrochlorinase) with different enantioselectivity towards (+) and (−) enantiomers of α-HCH (Suar et al., 2005). In addition, linB gene (haloalkane dehalogenase) of strain B90A mediates an additional step of β-HCH dehalogenation (Sharma et al., 2006). Genetic studies on this strain also provided the first evidence of close association between IS6100 and lin genes (Dogra et al., 2004).

Although, genetics and biochemistry of HCH degradation by strain B90A had been studied in detail (Lal et al., 2010) but, complete genetic repertoire especially the localization of lin genes and IS6100 could not be ascertained. Here, we report the complete genome of S. indicum B90A by combining the data obtained from Single Molecule Real Time Sequencing (SMRT) and Illumina platform. As a result, genome was assembled into one chromosome and three plasmids (pSRL1, pSRL2, and pSRL3) and also resolved the controversy of copy number and precise location of lin genes and IS6100 on the replicons. For instance, the analysis revealed the presence of 26 copies of IS6100 elements instead of 11 (Dogra et al., 2004) or seven (Anand etal. 2012) as reported previously. Further, the data analysis and subsequent laboratory experiments revealed the genetic basis of extracellular brown pigment production by strain B90A which involves accumulation, oxidation, and polymerization of homogentisic acid due to absence of hmgA gene (homogentisate dioxygenase) in the complete genome reported in this study.

Materials and Methods

Genome Sequencing and Assembly of the Strain B90A

Genome of S. indicum B90A was re-sequenced using third generation SMRT sequencing from Genome Quebec, Canada. For this, a 20 kb library was constructed with genome coverage of ∼390×. The genomic reads were assembled using Hierarchical Genome Assembly Process (HGAP) at Genome Quebec, Canada. In silico gap filling was performed using the combined data from SMRT and Illumina Hiseq platforms (Anand etal. 2012).

Genome Annotation and Pathways Mapping

Functional annotation of strain B90A was performed with Rapid Annotation Subsystem Technology (RAST) server v 4.0 (Aziz etal. 2008) using Glimmer (Delcher et al., 2007). Single copy essential genes were identified by AMPHORA2 (Wu and Scott 2012). Insertion sequences were identified using IS finder (Siguier et al., 2006), and Island Viewer3 (Dhillon etal. 2015) was used for the prediction of genomic islands. PHAST online server (Zhou 2011) was used to identify integrated phages. Gene specific to lin system and xenobiotic degradation were mapped on B90A replicons using BRIG (Alikhan et al., 2011).

Genome and lin System Synteny of HCH Degraders

Pairwise Average Nucleotide Identity (ANI) (Konstantinidis and Tiedje, 2004) between all the replicons of three HCH degraders namely S. indicum B90A, S. japonicum UT26S and Sphingomonas sp. MM1 was calculated. Strain UT26S consists of two chromosomes; Chr-1 and Chr-2 and three plasmids; pCHQ1, pUT1, and pUT2. Strain MM1 possesses only one chromosome and five plasmids; pISP0, pISP1, pISP2, pISP3, and pISP4 (supplementary table 2, Supplementary Material online). These replicons were clustered hierarchically using Pearson correlation method and heat map was constructed on MeV4.9.0 (Saeed et al., 2003). Replicons identity was visualized using CIRCOS (Krzywinski et al., 2009). Synteny mapping of lin gene clusters among strains and highly identical replicons was done using annotation from BLASTX (Altschul et al., 1990) and Kablammo server (Wintersinger et al., 2015). Recruitment plot between strain B90A plasmids and Illumina platform reads from HCH dumpsite metagenome was constructed using MUMmer3.23 (Kurtz et al., 2004). Syntenic region between chromosome and pSRL1 were visualized using Microscope server (Vallenet et al., 2009).

Results and Discussion

General Genome Features of Sphingobium indicum B90A Replicons

The data obtained from SMRT (390X) and Illumina Hiseq (2000X) (Anand etal. 2012) were used to assemble the genome into four replicons: one chromosome (Size: 3,654,322 bp, CDS: 3533, coding density: 88.40%) and three plasmids namely, pSRL1 (Length: 139,218 bp, CDS: 128, coding density: 85.70%), pSRL2 (Length: 108,430 bp, CDS: 128, coding density: 89.90%) and pSRL3 (Length: 43,761 bp, CDS: 46, coding density: 88.97%) (fig. 1). The detection of 31 essential single copy genes indicated that the genome is complete (Wu and Scott 2012). Previously, the draft genome of strain B90A was reported to consist 149 contigs (Anand etal. 2012), with genome size of 4.08 Mb and 65.0% GC content. However, the present study revealed a complete genome size of 3.9 Mb and GC content of 65.0%. These attributes are within the range of sphingomonad genomes (Genome size: 3.4–5.9 Mb, GC content: 62–68%) (Aylward et al., 2013, Verma et al., 2014; Kumar et al., 2017). GC content of plasmids pSRL1, pSRL2, and pSRL3 were 64.9%, 66.7%, and 63.32%, respectively. Almost nearly equal GC content of plasmids and chromosome of strain B90A suggested that they have been acquired either from sphingomonads or closely related species of sphingomonads and should be stable (Nishida 2012). A total of 3602 coding genes were annotated in the complete genome in which 33% codes for hypothetical proteins. Thus, the function of nearly one-third genes is yet to be ascertained.

. 1.—

Circular plot of genetic elements of Sphingobium indicum B90A (A) Chromosome, (B) plasmid pSRL1, (C) plasmid pSRL2, and (D) plasmid pSRL3. From inside to outside; first, second, and third circle show GC content, GC skew, and replicon, respectively. All lin pathway genes and IS6100 elements were mapped over replicons along with Type II/IV secretion pathway and two component system.

In strain B90A, 26 copies of IS6100 element were identified, of which 20 were present on chromosome, four on pSRL1 and one each on pSRL2 and pSRL3. So far this number is maximum among the reported range from five in S. japonicum UT26S (Nagata et al., 2011) to 24 in S. quisquilarum P25 (Verma et al., 2014). In the previous studies based on southern blot hybridization and draft genome, only 11 (Dogra et al., 2004) and seven (Anand etal. 2012) copies of IS6100 were reported, respectively. This discrepancy might be due to the low resolution of the DNA bands obtained from southern hybridization. Moreover, in the draft genome, these insertion sequences might have been missed as repeat regions are more likely to be mis-assembled in de-novo assembly (Salzberg and Yorke 2005). Nevertheless, the role of IS6100 in genome shuffling is well known (Sangwan et al., 2014) and high copy number reflected ongoing shuffling and rearrangements in B90A genome especially in relation to lin genes.

Interestingly, along with linA1 and linA2 (will be discussed in later section), strain B90A also contained two variants of linE, that is, linEa (chromosome) and linEb (pSRL3) which differed by 53% and encode for hydroquinone 1,2-dioxygenase and chlorohydroquinone/hydroquinone 1,2-dioxygenase, respectively (supplementary table 1, Supplementary Material online). This demarcates its ability to metabolize both chlorohydroquinone and hydroquinone formed during degradation of γ-HCH (Nagata et al., 1999). Also, duplication of linD, linR, linF, and linGHIJ (supplementary table 2, Supplementary Material online) with linE variants suggested an efficient lower HCH degradation pathway.

Genomic Islands (GIs) and Integrated Phage

In strain B90A and other HCH degrading sphingomonads (Kaur et al., 2013; Negi etal. 2014) accessory genome, especially genomic islands were not described precisely. Genomic islands were predicted using SIGIHMM, IslandPick and IslandPath-DIMOB, and GIs identified by atleast two of the methods were considered. A total of 14 genomic islands (GIs) were predicted and designated as GI1-GI14 (fig. 1). Nitrite-sensitive transcriptional repressor NsrR and nitric oxide dioxygenase (NOD) were identified in GI3. NsrR and NOD are essentially needed for endurance against toxic nitric oxide produced by denitrifying bacteria in the rhizosphere by converting it into nitrous oxide (N2O) or nitrate () (Tucker et al., 2008). In GI5 and GI9, genes for type IV secretion (T4SS) system were identified that are known to mediate the extracellular transport of virulence factors (Yeo and Waksman 2004). Out of 11 genes of virB-virD system (Yeo and Waksman 2004), all except virB2 and virB7 were found in strain B90A However, traC and traI, homologs of virB2 and virB7 were present in GI6 which suggest complete set of virulence factors are present highlighting its role in combating the stress in rhizospheric environment.

A complete integrated prophage in chromosome (Location: 1255173–1303366 bp) was identified. It contained 67 CDSs in which nearly 50% codes for phage-associated proteins while the other codes for hypothetical proteins. Integration and maintenance of a phage by strain B90A might be an adaptation against environmental stress (Burrus etal. 2004). But the precise role of this phage in the genome of strain B90A is yet to be explored due to the hypothetical proteins.

Evidence of linA (HCH Dehydrochlorinase) Duplication: Origin of linA1 from linA2

The linA gene that encodes for HCH dehydrochlorinase holds its importance in the degradation of HCH isomers as it mediate the first step of dehydrochlorination of α-, γ-, δ-, and ɛ-HCH (Lal et al., 2010). But, its origin and acquisition is still not known. For instance, unlike HCH degrading sphingomonads, that have been reported to contain a single copy of linA gene, strain B90A harbor two variants, linA1 and linA2 with enantioselective preference towards (+) and (−) α-HCH (Suar et al., 2005), respectively. The linA gene of strain UT26S (Nagata etal. 2001) is 100% identical to linA2 of strain B90A (Lal et al., 2010; Verma et al., 2014). Previous study by Sangwan et al., (2014) revealed that linA is environment specific and was even absent in the last common ancestor of Sphingobium spp. Exceptionally, linA1 has only been reported in strain B90A (Kumari et al., 2002; Suar etal. 2004; Boltner etal. 2005; Ceremonie etal. 2006; Malhotra et al., 2007). linA1 and linA2 were identified on pSRL1 and chromosome, respectively, which were aligned to fetch the syntenic regions between the replicons. The analysis showed that flanking sequences of linA1 (upstream) and linA2 (downstream) possess sequence similarity of ∼15 kb (with >80% identity) (supplementary fig. 1, Supplementary Material online). In addition, 22 bp at C-terminus of linA1 and N-terminus of IS6100 were identical, as reported previously (Dogra et al., 2004). Along with the oppositely oriented IS6100 elements, this 15 kb region reflects the possibility of homologous recombination event resulting into duplication of linA in strain B90A. This is further supported by the previous study in which under stress condition, gene that confer adaptation against stress was duplicated between replicons of proteobacterial strains (Slater etal. 2009). Further, genes on the plasmid were found more prone to the sequence variation as compared to the chromosome (Slater etal. 2009). Hence, this suggests that linA2 was perhaps acquired initially by sphingomonads under HCH stress followed by its duplication and sequence divergence as linA1. This gene might have been selected during evolution as it not only enhanced the degradation ability of strain B90A but also possess enantioselectivity towards (+) enantiomer of α-HCH as ancient known linA or linA2 prefers (−)α-HCH as the substrate (Suar et al., 2005). Hence, sequence identity between the flanking regions of both the linA and multiple copies of IS6100 provides an evidence of linA duplication followed by sequence variation in linA1 which led to the current scenario of strain B90A genome. Further, this data and metagenome of HCH dumpsite (Sangwan etal. 2014, Shrivastava etal. 2015) suggest that linA is very recent and evolving at a rate beyond scientific imagination.

Replicon Based Identity among HCH Degraders

Due to the absence of complete genome, comparative studies on HCH degrading sphingomonads were limited to the draft genomes (Verma et al., 2014; Pearce etal. 2015; Kumar et al., 2017) hence, could not infer replicon based phylogeny. Here, replicons of S. indicum B90A (fig. 1) and two other HCH degraders with complete genomes, S. japonicum UT26S (Nagata et al., 2011), Sphingomonas sp. MM1 (Tabata et al., 2013) (supplementary table 2, Supplementary Material online) were clustered hierarchically with Pearson Correlation method using average nucleotide identiy (ANI) (fig. 2, Supplementary Material online). The replicons were clustered into two major clades of which one clade consists of all chromosomes with pISP0 and pSRL2, while the other clustered the remaining nine plasmids. Interestingly, both pSRL1 and pSRL2 shared high ANI (≥90%) with chromosome 2 of strain UT26S (fig. 2). Previously, many species belonging to the phylum proteobacteria have been reported to mediate inter-genomic transfers for the construction of secondary chromosome from primary chromosome and plasmids (Slater etal. 2009). This possibility cannot be ruled out as plasmids of strain B90A (pSRL1 and pSRL2) and chromosome 2 of strain UT26S showed high identity which provides a sign of probable inter-genomic transfer between replicons of B90A and UT26S.

. 2.—

Comparison between genetic elements of Sphingobium indicum B90A, Sphingobium japonicum UT26S and Sphingomonas sp. MM1. (A) Similarity plot with B90A genetic elements: B90A-CHR (chromosome), pSRL1, pSRL2, and pSRL3; UT26 genetic elements: UT26-CHR1 (chromosome 1), UT26-CHR2 (chromosome 2), pCHQ1, pUT1, and pUT2; MM1 genetic elements: MM1-CHR (chromosome), pISP0, pISP1, pISP2, PISP3, and pISP4. (B) Heat map on pairwise Average Nucleotide Identity between the replicons (chromosomes and plasmids) of the strains.

Further, pairwise ANI values of whole genome (chromosome/s and plasmids) showed that both strains B90A and UT26S shared ∼83% identity with strain MM1, 5% more than their chromosomal identity, that is, ∼78%. Generally, plasmids due to their role in horizontal gene transfer (HGT) are considered as the vehicles of genome variability (Bergstrom et al. 2000; Gogarten et al. 2002). But here, inclusion of plasmid sequences has increased the identity between the whole genome of HCH degraders which highlights their acquisition by the last common ancestor (LCA). High identity (∼90%) between pISP3 of strain MM1 and pSRL3 of strain B90A which harbored linD (reductive dechlorinase), linE (ring cleavage oxygenase), linR (transcriptional regulator) and an IS6100 element (supplementary fig. 2A, Supplementary Material online) further supports acquisition of pSRL3 or pISP3 by LCA of HCH degraders. Metagenomic recruitment plot of pSRL3 using the data obtained from the HCH dumpsite, India, also determined the co-occurrence of this plasmid at the site from where different HCH degrading sphingomonads have been isolated including Sphingomonas sp. MM1 (Tabata et al. 2013). Also, this suggest the active role of pSRL3 in acquisition of lower HCH degradation pathway in sphingomonads (supplementary fig. 2B, Supplementary Material online).

lin Gene Clusters and IS Elements: Southern Hybridization versus SMRT Sequence Analysis

Previously, copy number of lin genes and IS6100 in strain B90A, were predicted by southern hybridization which reported the presence of at least 11 copies of IS6100 element (Dogra et al. 2004) while 26 copies of IS6100 were identified in the complete genome. To resolve this ambiguity, in silico restriction mapping was done using BamHI and compared with the previous gel electrophoresis profile (Dogra et al. 2004) (supplementary fig. 3, Supplementary Material online). In silico profile predicted 21 fragments (containing one or more IS6100 elements) instead of 11 fragments obtained previously on the gel electrophoresis. A similar comparison was drawn for the lin genes using previous studies (Kumari et al. 2002; Dogra et al. 2004; Malhotra et al. 2007) and in silico restriction digestion with BclI, HindIII, BamHI, and BglII (supplementary fig. 4 and table 3, Supplementary Material online). Here, in silico and agarose gel profiles of lin genes were similar which might corresponds to less ambiguity due to the presence of one or two copies of lin genes as compared to 26 copies of IS6100.

Pigment Biosynthesis by Strain B90A

Production of extracellular brown pigment is a characteristic feature of strain B90A (Sahu et al., 1990) that was not observed in other HCH degrading sphingomonad. Genome analysis of strain B90A revealed the presence of gene cluster encoding incomplete pathway for tyrosine and phenylalanine metabolism. It was also found from the complete genome analysis of strain B90A that this strain lacked homogentisate dioxygenase gene, hmgA of the pathway thereby making it incapable of utilizing homogentisic acid (HGA) as substrate. HGA is a stable metabolite formed during metabolism of tyrosine and phenylalanine. Recently, Wang et al., (2015) reported that in Aeromonas spp. due to absence of hmgA gene, HGA got accumulated and further converted to pyomelanin like brown pigment. Genome analysis of strain B90A revealed the presence of four essential genes for biosynthesis of HGA such as, phhA (phenylalanine-4-hydroxylase), phhB (4-alpha-carbinolamine dehydratase), hppD (4-hydroxyphenyl pyruvate dioxygenase), and aspC (aspartate aminotransferase). In general, this pathway further extends with hmgABC that catalyzes the conversion of homogentisic acid into fumarate and acetoacetate. But due to the absence of hmgA gene in B90A, HGA accumulated and underwent auto-oxidation and polymerization resulting in pyomelanin production (supplementary fig. 5, Supplementary Material online). Further, comparative in silico annotation using the functional genes involved in HGA based pigment synthesis reported previously from P. aeruginosa UCBPPPA14 (Hunter and Newman, 2010) as reference; nearly 16 functional genes were identified in the genome of strain B90A (supplementary table 4, Supplementary Material online). The analysis also revealed the absence of tyrosinase gene in B90A genome indicating that the synthesis of pigment is not through the oxidation of tyrosine to L-3,4-dihydroxyphenylalanine (L-DOPA) and dopaquinone (Plonka and Grabacka 2006). This suggested that strain B90A does not synthesize eumelanin or phaeomelanin like pigment (Donlon et al., 1983). Also, genes encoding laccases and polyketide synthases which may contribute to melanin synthesis were not identified in the genome of strain B90A. This confirmed that HGA biosynthesis in strain B90A is responsible for the production of extracellular brown pigment. Moreover, comparative analysis of tyrosine metabolism has also revealed that phylogenetically related strains of Sphingobium that harbored hmgA gene do not produce brown pigment.

Conclusion

The study revealed the relevance of SMRT technique in constructing the complete genome of Sphingobium indicum B90A and has resolved several issues and ambiguity that existed previously. The number of plasmids and chromosomes in strain B90A was determined. Strain B90A consisted of four replicons, one chromosome and three plasmids; pSRL1, pSRL2, and pSRL3. The analysis determined the precise location of lin genes on the replicons, in order, such as linA2, linB, linDER, linF, linGHIJ, and linKLMN on the chromosome; linA1, linC, and linF on pSRL1 and linDEbR on pSRL3. Hence, along with linA, two variant/copy of linDER, linGHIJ, and linF were also identified. Complete genome analysis also led to filling of the gaps in genomic information from the previous studies that resulted in identification of 26 copies of IS6100; maximum among HCH degraders reported till date. The sequence similarity between adjacent regions of linA variants on chromosome and pSRL1 provided the strong evidence of its duplication via. homologous recombination events rather than their independent acquisition. At last, characteristic brown pigment production by strain B90A was identified as pyomelanin which formed by the accumulation and oxidation of homogentisic acid due to the absence of hmgA gene.

Supplementary Material

Supplementary data are available at Genome Biology and Evolution online.

Click here for additional data file.

Click here for additional data file.