Phylogeny, distribution and potential metabolism of candidate bacterial phylum KSB1

Collection of KSB1 genomes

A total of 42 nonredundant KSB1 genomes were downloaded from GTDB, NCBI, JGI databases and published papers (Table S1) (November, 2020) and were filtered with following cutoff values: completeness score ≥ 50%, contamination rate ≤ 10% and QS > 50 (QS = completeness-5*contamination) (Almeida et al., 2019). KSB1 genomes were also binned from the metagenomes of the Mariana sediments collected 5,400 to 10,911 m depths during R/V DY37-II, TS01 and TS03 (Cui et al. 2021). Fastp (Chen et al., 2018a) (v.0.20.0) was used for the quality control of metagenome raw data. Repeated Illumina sequences were removed by Fastuniq (Xu et al., 2012) and the clean reads were assembled with SPAdes (v.3.13) (Bankevich et al., 2012). Contigs >2,000 bp were used for genome binning with MetaWRAP (v.1.2) integrated with three binning tools, followed by a treatment with bin_refinement module (Uritskiy, DiRuggiero & Taylor, 2018). The MAGs with qualified completeness (higher than 50%) and contamination (lower than 10%) were selected by CheckM (Parks et al., 2015). The MAGs affiliated with KSB1 were selected from the classification result of GTDB-tk (v.1.4.0) software (Chaumeil et al., 2019) integrated with GTDB release95 database. One KSB1 genome named vent-69 had been obtained from a Mid-Atlantic hydrothermal sediment metagenome (accession: SAMN10350645).

Calculation of KSB1 relative abundance in different niches

16S rRNA gene sequences of the KSB1 genomes were extracted to create a dataset. The 16S miTags of the Tara Ocean data were downloaded from http://ocean-microbiome.embl.de/data/16SrRNA.miTAGs.tgz. KSB1 16S miTags (metagenomic Illumina tags) were identified from the Tara Ocean miTags (Moran, 2015) by BLASTn (Gish & States, 1993) (v.2.9.0) (-evalue 1e–05) against the 16S rRNA dataset of KSB1. The mapped KSB1 16S miTags with 97% identity and at least 100 bp in length were further selected from the BLASTn result as KSB1 16S miTags for calculation of their relative abundance in the marine water samples. Raw data of public metagenomes (fastq files) were downloaded from NCBI and were subjected to quality control as described above. The relative abundance of KSB1 MAGs in the metagenomes was calculated by coverM (https://github.com/wwood/CoverM, -m relative_abundance; -min-read-aligned-length 50; -min-read-percent-identity 0.99; -min-covered-fraction 0) after mapping with bwa mem (Li & Durbin, 2009) (default parameter) and sorting with samtools (default parameter) (Li et al., 2009).

Gene annotation and metabolic reconstruction

Genomes from public databases and this study were used for gene annotation. Open reading frames (ORFs) were predicted by Prodigal (Hyatt et al., 2010) (v.2.6.3) and were searched against KEGG (Kyoto Encyclopedia of Genes and Genomes) database (release 92) by Kofamscan (Aramaki et al., 2020) (v.1.0.0; -f mapper), COG (Cluster of Orthologous Groups of proteins) database (COG_2019_v11.0) (Tatusov et al., 2000) by BLASTp (Eddy, 1995) (BLAST+ v.2.9.0) (-evalue 1e–05) and CAZy database (dbCAN-HMMdb-V7) by hmmscan (Eddy, 1995) (v.3.2.1) with default settings. Phage was predicted in http://phaster.ca/.

Phylogenetic analysis of genomes and proteins

Single copy marker proteins in the KSB1 MAGs and reference genomes were identified with GTDB-tk classify_wf (Chaumeil et al., 2019). The marker proteins were filtered by their presence in 80% of the MAGs. The selected marker proteins were concatenated and used for reconstruction of a phylogenomic tree with iqtree2 (Minh et al., 2020) (v.2.1.0; -m MFP) after multiple sequence alignment using MAFFT (Katoh & Standley, 2013) (v7.453) and alignment optimization using trimAl (Capella-Gutierrez, Silla-Martinez & Gabaldon, 2009) (v.1.4). A phylogenomic tree of the genomes from FCB superphylum (Table S2) was constructed by iqtree2 with MFP model using 43 conserved proteins selected by checkM (Parks et al., 2015).

The amino acid sequences encoded by narG, nrfA and nosZ identified in KSB1 genomes were searched against the NCBI_nr database. The most similar homologous sequences from different phyla were retrieved for phylogenetic tree construction as described above but with MFP+LM model. The protein sequences of NosZ were obtained from FunGene database (Fish et al., 2013) and clustered by CD-HIT (Li & Godzik, 2006) (-c 0.76 -A 0.8). An unroot phylogenetic tree was built for KSB1 as mentioned above.

MAG availability

The MAGs of KSB1 binned from the Mariana sediments and the Mid-Atlantic Ridge hydrothermal sediment were submitted to the National Omics Data Encyclopedia (NODE) with OEP002159 as project accession number and OER184014 as run ID.

Phylogenomics of KSB1

A total of 89 genomes of KSB1 were collected from public databases, including GTDB, JGI and NCBI. In addition, 14 KSB1 genomes provided by published papers were recruited manually as well (Table S1). Fifteen KSB1 MAGs retrieved from metagenomes of the Mariana sediments with depths ranging from 5,400 to 10,953 m and the Mid-Atlantic hydrothermal vent by this study were added into the KSB1 dataset. After dereplication and quality control for the 118 genomes, 44 of them were retained for further study (Table S1), which included 11 MAGs from this study in size of 2.87~5.37 Mbp (Table 1).

The phylogenomic tree constructed by using 39 conserved marker proteins displayed four phylogenetic clades of KSB1, which were then named as ‘clades I–IV’ (Fig. 1A). The 11 KSB1 MAGs from this study were distributed into two clades (clade II and clade IV). Particularly, most KSB1 MAGs of the Mariana sediments were grouped into clade IV. According to the microbial phylogenetic tree of hydrothermal sediments (Dombrowski, Teske & Baker, 2018) and a previous phylogenetic inference of KSB1 (Youssef et al., 2019), KSB1 might be affiliated with the Fibrobacteres-Chlorobi-Bacteroidetes (FCB) superphylum or be a sister phylum of the superphylum. To examine this hypothesis, we constructed a phylogenomic tree using some high-quality genomes of the FCB superphylum with Proteobacteria and Therrabacteria serving as the outgroup (Table S2). The topological structure of the tree showed that KSB1 was placed into one monophyletic branch within the FCB superphylum and was adjacent to SAR406 (Huang & Wang, 2020) and ‘Candidatus Tianyabacteria’ (Cui et al., 2021) (Fig. S1).

The genome size and GC content of the 44 KSB1 MAGs were calculated. The mean genome size of clade I with six MAGs was around 7 Mbp, which was significantly larger than other three clades (t-test; p < 0.001) (Fig. 1B). The median GC content of clade I was 52.50%, significantly higher when compared to other clades with a t-test (p < 0.001 for clade II with 11 MAGs; p = 0.004 for clade III with 12 MAGs; p < 0.001 for clade IV with 15 MAGs) (Fig. 1C). The isolation sources of the 44 MAGs were summarized to demonstrate distribution specificity of KSB1 clades. These MAGs of clade I were uniquely from bioreactors, while clade IV was exclusively from marine sediments (Fig. 1D). In contrast, MAGs of clade III were isolated from different environments including freshwater sediment, wetland sediment and tailing pond (Fig. 1D), indicating that clade III is broadly distributed. For the KSB1 MAGs of clade II, 18.18% genomes were obtained from marine sediments and 81.82% genomes were identified in hydrothermal vent environments (Fig. 1D).

A previous study has indicated that genome size might correlate with GC content (Wu et al., 2012). The large genome size of clade I may be ascribed to a large number of genes and mobile genetic elements in the genomes. It has been reported that mobile DNA density increased when the genome size was enlarged (Newton & Bordenstein, 2011). Prediction of transposases against KEGG and COG databases revealed that these genes were more frequently present in clades I and IV (Table S3). Therefore, transposases might be one of the factors driving the genome expansion of clade I (Fig. 1B). Phage or CRISPR/Cas arrays might be an effective mobile element (Al-Shayeb et al., 2020; Ali et al., 2020; Sanderson et al., 2020) for the expansion as there were more than one copy of cas2 gene in the clade I MAGs (Table S3) and phage was predicted in 66.67% of the MAGs of clade I (Table S3). In addition, multicopy genes could be predicted in MAGs of clade I (Table S3). Although studies indicated that a larger genome of microbes endows greater versatility (Nielsen et al., 2021), genome size as a trait was nearly independent from cell size and growth rate of various bacterial ecological types (Nielsen et al., 2021; Westoby et al., 2021). Considering the high diversity of isolation sources, the roles of KSB1 bacteria might be differentiated notably.

Relative abundance of KSB1 in different niches

To evaluate the distribution of KSB1 in different environments, relative abundance of KSB1 was assessed using percentage of 16S miTags and coverage of KSB1 MAGs in metagenomes. The metagenome raw data for the calculation were downloaded from NCBI (Table S4). Our results showed that the relative abundance of MAGs belonging to clade I was most abundant in bioreactor (Fig. 2; Table S5), indicating that KSB1 of clade I may be likely one of the representative bacterial groups in wastewater treatment. Clade IV was more prevalent in marine sediment than other environments. All the four clades were present in sea water and were more abundant than in groundwater (Fig. 2).

To examine vertical distribution of KSB1 in marine waters, relative abundance of KSB1 16S miTags in those of Tara Ocean project (Moran, 2015) was calculated. Clade II was relatively more abundant in the oceans, compared to the other clades (Fig. S2). Clade II was the most abundant at 5-m surface layer (0.65%) than other zones (Fig. S2; Table S5), whereas clade IV was detectable in marine waters below 200 m (Fig. S2; Table S5). Nevertheless, the low abundance of KSB1 in the Tara Ocean data suggests limited distribution of the KSB1 bacteria in oxic marine waters.

Carbohydrate-active enzymes (CAZYmes) in KSB1 genomes

KSB1 may be broadly distributed due to their potential functions on degrading complex carbohydrates. An analysis of CAZymes in the 44 KSB1 MAGs revealed diverse GH and GT classes (Fig. 3A). Particularly, GH109, GH23, GT4, GT2, CE10 (esterase) and CBM50 (LysM) could be identified in all MAGs (Table S6). GH23 includes lytic transglycosylases with helice D, F and beta-sheet, acting on peptidoglycan to cleave the glycosidic linkage between N-acetylglucosaminyl and N-actetylmuramoyl residues to produce cyclic 1,6-anhydro-N-acetylmuramic acid (anhMurNAc) (Harding et al., 2020; Scheurwater, Reid & Clarke, 2008). CBM50 can interact with chitin and peptidoglycan (Bertucci et al., 2019). GT2 and GT4, including α-glucosyltransferase and chitin synthase, were dominant families of GTs (Bohra, Dafale & Purohit, 2019). These results indicate that the five types of CAZyme subfamilies were found across KSB1 and potentially allow them to obtain nutrients from various organic substrates such as chitin, peptidoglycan and other components of cell wall.

Clades I and III contained more CAZYmes of different subfamilies, compared to the other clades (Fig. 3B; Table S6), which indicated that clades I and III of KSB1 encode a broader repertoire of CAZymes than clades II and IV apart from CAZYmes involved in auxiliary activities (AA). Furthermore, CBM, GH and PL were absent in many MAGs of clades II and IV (Fig. 3C; Table S6). Notably, CAZymes with a potential secretion signal (Dombrowski, Teske & Baker, 2018) were different in distribution among four clades. GH28, GH33 and PL9 were absent in clades II and IV (Fig. 3C) and this seems to be a result of low availability of labile carbohydrates in deep-sea hydrothermal vent or hadal sediments (Richardson et al., 1995). As a contrast, such CAZYmes were more common in KSB1 MAGs of clade I that were identified in bioreactor sludge enriched with abundant and complex organic carbon sources (Mata et al., 2020) (Fig. 3C). The high variety of CAZymes in clade III agrees with their diverse isolation sources as shown in Fig. 1D.

Core metabolic genes detected in KSB1 genomes

The functional genes responsible for utilization of carbohydrates, associated to nitrogen, sulfur and selenate metabolisms, were detected in the KSB1 MAGs (Fig. 4). The core genes of gltA and sucD related to citric acid cycle (TCA) were present in four clades of KSB1 (Fig. 4). The core genes of rpe, rpiA, rpiB and tktA involved in pentose phosphate pathway (PPP) were also present in four clades of KSB1 (Fig. 4). These results indicate that central carbon metabolism is nearly complete in KSB1. About half of MAGs of clade II encoded genes involved in hydrocarbon metabolism (Fig. 4), although the MAGs mainly came from hydrothermal environment (Fig. 1D). There has been reported abiogenic hydrocarbon in hydrothermal system (McDermott et al., 2015; Proskurowski et al., 2008), which may be the carbon source for growth of KSB1 that encoded genes associated to hydrocarbon metabolism. KSB1 MAGs harbored more than one copy of paaF, pccB, epi, mcmA1 and mcmA2 genes that are involved in short-chain hydrocarbon transformation (Fig. 4). They were abundantly present in clades I and IV of KSB1. It had been reported that an operon consisting of 14 paa genes encode enzymes to degrade phenylacetate (Teufel et al., 2010). PaaABCDE catalyze phenylacetyl-CoA to ring 1,2-epoxyphenylacetyl-CoA (Teufel et al., 2010). Paa converts 3-hydroxyadipyl-CoA to 3-oxoadipyl-CoA with NADH as a byproduct (Teufel et al., 2010). PaaABH coding genes were detected in clade I of KSB1 (Table S7), which indicated that the clade I might catabolize phenylacetic acid or act on the intermediates of the whole pathway to obtain energy. However, most of the paa genes (paaZCDEGIJKXY) were not identified in clade I; catalytic experiment of Paa complex of clade I in phenylacetic acid utilization is needed in future work. mcmA1 and mcmA2 encode α-subunit of methylmalonyl-CoA mutase that takes part in propanoate pathway (Han et al., 2013). PccB participates in the conversion of propionyl-CoA to methylmalonyl-CoA (Wongkittichote, Ah Mew & Chapman, 2017), which might be subsequently converted by McmA1 or McmA2 to enter TCA cycle (Bobik & Rasche, 2001). EPI (MCEE) is a methylmalonyl-CoA epimerase responsible for degradation of odd chain-length fatty acids and branched-chain amino acids (Dobson et al., 2006). The interconversion of D- and L-methylmalonyl-CoA might be performed by EPI, which is a key step in propanoyl-CoA to succinyl-CoA for TCA cycle (Dobson et al., 2006). Since all these genes have been identified in the KSB1 genomes, the metabolism of propionyl-CoA to succinyl-CoA might be employed by KSB1 for short hydrocarbons degradation into TCA cycle (Fig. 4), which was similar to previous report that KSB1 was involved in anaerobic degradation of hydrocarbon in hydrothermal sediments (Dombrowski et al., 2017). HADH that was present in clades I and IV is a 3-hydroxyacyl-CoA dehydrogenase gene involved in fatty acids metabolism as described in a previous study that KSB1 has the capacity for beta-oxidation (Baker et al., 2015). Coupled with the other core genes such as pccB, mcmA1, mcmA2 and epi, the result suggests that the clades I and IV might break down fatty acids and hydrocarbons for energy.

narG, nrfA and nosZ for denitrification were only present in clade I of KSB1 (Fig. 4). Nitrate could be reduced to nitrite by nitrite oxidoreductase encoded by narG and the following reduction of nitrite to ammonia could be finished with the function of nrfA (Giblin et al., 2013). This indicates that MAGs affiliated with KSB1 clade I may be involved in dissimilatory nitrate reduction in bioreactor, which was not reported previously (Youssef et al., 2019). nosZ gene involved in the last step of denitrification to produce nitrogen gas (Giblin et al., 2013), was revealed in five out of six MAGs of clade I (Fig. 4). However, the other genes related to denitrification (nirK or norBC) could not be found in MAGs of clade I (Table S7). The genes only identified in the clade I might be the results of lateral gene transfer as indicated by the genome expansion (Fig. 1B). Phylogenetic trees were built to examine the origin of NarG and NrfA predicted in KSB1 clade I. Our results showed that the NarG sequences of KSB1 clade I were grouped with the homologs from Acidobacteria, Rokubacteria and NC10, while the NrfA sequences were adjacent to those derived from Anaerolineae, ‘Candidatus Jettenia’ and ‘Candidatus Brocadia’ (Fig. S3). In addition, the NosZ sequences of KSB1 clade I were approximate to the homologs from Ignavibacteria and were associated with sec-type signal peptide (Fig. S4; Table S8) (Jones et al., 2013). There were some genes responsible for sulfur metabolism such as sat (sulfate adenylyltransferase), phsA (polysulfide reductase chain A) and cysK (cysteine synthase) in the KSB1 MAGs (Fig. 4). sat and cysK were identified in almost all KSB1 except clade II. Particularly, sat, encoding a protein responsible for activation of inorganic sulfate, was used for catalysis of sulfate to adenylyl sulfate (Fauque & Barton, 2012). PhsA, mostly present in clade I, converts thiosulfate to sulfide for synthesis of cysteine by CysK (Chen et al., 2018b). This suggests that KSB1 might take part in the assimilatory sulfate reduction in bioreactor. selA, selB and selD genes were predicted in MAGs of KSB1 except clade II (Fig. 4). SelA (selenocysteine synthase) and SelD (selenophosphate synthase) were required for selenocysteine synthesis (Leinfelder et al., 1990). SelB is a Sec-specific elongation factor for the incorporation of selenocysteine into proteins (Sheppard et al., 2008). With the presence of the three genes in KSB1, the biosynthesis of selenocysteine may take place in KSB1. Selenocysteine can be incorporated into proteins as well when selC (selenocysteyl-tRNAsec) is mutant (Zorn et al., 2013). The protein containing selenocysteine (selenoproteins) functions in antioxidant system and redox regulation of signal pathway (Zhang et al., 2020). The prevalence of selABD genes and absence of selC gene (Table S7) in KSB1 genomes suggest that selenoproteins might be biosynthesized by KSB1 to resist stress in diverse niches.

Metabolism reconstruction of KSB1

To learn more about characteristics of KSB1 lifestyle, metabolism reconstruction was performed. Almost all genes related to flagellar biosynthetic proteins, flagellar assembly proteins and other flagellar structure proteins were identified in KSB1 genomes affiliated with clades II and IV (Fig. 5; Table S7). In addition, the stator element of the flagellar motor complex MotA, methyl-accepting chemotaxis protein MCP and the proteins of two-component system controlling chemotaxis (Table S7) were encoded by clades II and IV. The flagellar and chemotaxis responsible for motility might be useful for KSB1 clades II and IV to explore nutrients in oligotrophic deep-sea environment. The KSB1 genomes of clade I contain the genes coding for some ABC transporters for uptake of biotin, ferrous and ferric ion, branched-chain amino acids, molybdate and maltooligosaccharide (Fig. 5; Table S7). These transporters may be employed by KSB1 clade I to import organic and inorganic nutrients from bioreactor used for waste treatment. The (thio)sulfate transporter related gene was only present in clade I MAGs (Fig. 5). Thiosulfate can be oxidized to sulfite by thiosulfate/3-mercaptopyruvate sulfurtransferase (TST) or reduced to sulfide by thiosulfate reductase/polysulfide reductase chain A (PHSA) in anoxic conditions (Jorgensen, 1990). This indicates that KSB1 of clade I might take part in “thiosulfate shunt” in the bioreactor sludge. These pathways associated to sulfur might occur in the KSB1 inhabiting in sulfur-rich black mud (Tanner et al., 2000), anoxic zone with a low-H₂S rather than high-H₂S concentration (Ley et al., 2006). The genes named cydA and cydB were also predicted in clade I (Table S7). cydAB encode subunits of cytochrome bd that is expressed under low oxygen condition (Borisov et al., 2011; Kranz et al., 1983), which indicates that KSB1 could live in O₂-limited environments (Borisov et al., 2011). Nitrogen regulation genes (ntrY and ntrX) were all detected in KSB1 genomes of clade I (Table S7), which is likely required for controlling the level of cytoplasmic nitrogen (Carrica et al., 2012). In addition, the key genes (algD, algR and algZ) of alginate biosynthesis (Leech et al., 2008; Wu et al., 2015) were detected in KSB1 genomes of clade I (Fig. 5; Table S7), suggesting that they might be involved in alginate biosynthesis as one strategy for keeping carbon storage.

The molybdate transporter genes such as modA (coding for molybdate transport system substrate-binding protein) and modB (coding for molybdate transport system permease protein) were identified in KSB1 genomes of clades I and IV (Fig. 5; Table S7). In addition, the genes involved in Moco (molybdenum cofactor) biosynthesis (Nichols & Rajagopalan, 2005) were also identified (Fig. 5; Table S7). The function of Moco in KSB1 genomes is not clear yet, although it has been reported that the impairment of Moco biosynthesis affected mobility, anaerobic respiration and biofilm formation (Andreae, Titball & Butler, 2014). Moco required for molybdoenzymes was probably important for bacteria to adapt in harsh or dramatically changing redox condition (Leimkuhler & Iobbi-Nivol, 2016), and therefore these genes may encode proteins for survival of KSB1 in bioreactor or deep-sea sediments.

When the lipopolysaccharide (LPS) biosynthesis pathway was examined, the genes (lpxABCDKL and kdtA) involved in this pathway were identified in almost all KSB1 genomes except clade IV (Fig. 5; Table S7). The encoded proteins might use UDP-N-acetyl-alpha-D-glucosamine as substrate for biosynthesis of lauroyl-KDO2-lipid IV(A) and LPS (Wang, Quinn & Yan, 2015). The potential capacity of LPS biosynthesis suggests that KSB1 might be gram-negative bacteria. The pellicle (PEL) polysaccharide-dependent biofilm formation related genes (pelADEFG) (Whitfield et al., 2020a; Whitfield et al., 2020b) were present in KSB1 genomes of clade IV (Fig. 5; Table S7). Biofilm is a special colony that contains microbial cells and extracellular matrix (Davies et al., 1998). The presence of these genes in the clade IV KSB1 genomes suggests that biofilm may be important for their survival in nutrient-poor deep-sea sediments.

A set of genes involved in heme biosynthesis and ferrous iron transportation were identified in KSB1 genomes of clades I and IV (Fig. 5; Table S7). The Fe-coproporphyrin III biosynthesis might be finished in KSB1 genomes of clades I and IV due to the presence of hemABCDELHY (Fig. 5; Table S7). The presence of hemN in clade I indicated that they may biosynthesize heme A, a component of cytochrome oxidases (Hederstedt, 2012). It had been reported that hemA and hemL were necessary for heme biosynthesis and electron transfer in anaerobic respiratory metabolism (Frankenberg, Moser & Jahn, 2003; Zumft, 1997). Collectively, this might be a strategy of KSB1 clade I inhabiting bioreactor sludge to obtain energy under anaerobic conditions by denitrification using heme nitrite reductase with iron regulating heme biosynthesis.