Sulfur cycling, one of the key biological processes in marine sediments, is dominated by sulfate-reducing bacteria (SRB) and sulfur-oxidizing bacteria (SOB). Dissimilatory sulfate reduction, which is mediated by SRB, is considered the main process in the biomineralization of organic matter in marine sediments and might account for up to 50% of organic matter mineralization in most continental shelf sediments (Jørgensen, 1982). In addition, as much as 12–29% of the organic carbon flux on the ocean seafloor is channeled through sulfate reduction (Bowles et al., 2014). Interestingly, 80–95% of the massive amount of hydrogen sulfide formed through sulfate reduction is recycled within sediments and gradually oxidized back to sulfate (Jørgensen and Nelson, 2004). Thus, SRB and SOB control the key processes of organic matter degradation and the biogeochemical cycling of sulfur and carbon. Studying the community structure and diversity of SRB and SOB is important for revealing the roles of these bacteria in the biogeochemical cycles of carbon and sulfur and for providing insight into the biological factors driving the marine sulfur cycle.

SRB and SOB show high diversity, including both phylogenetic and metabolic diversity (Ghosh and Dam, 2009; Müller et al., 2014). To explore the environmental abundance and diversity of SRB and SOB, genes encoding key enzymes in the sulfate reduction and sulfur oxidation biochemical pathways have been used as molecular markers. For example, the dissimilatory sulfite reductase gene dsrAB, which encodes a key enzyme that catalyzes the last, essential step in the dissimilatory sulfate reduction pathway, has been frequently employed as a functional gene in studies of SRB in various environments (Foti et al., 2007; Pester et al., 2012; He et al., 2015). The known SOB have been demonstrated to use different enzymes, pathways, and electron transport and energy conservation mechanisms for the oxidation of sulfide. The sulfur-oxidizing (Sox) (Meyer et al., 2007), reverse dissimilatory sulfite reductase (Loy et al., 2009), adenosine-5-phosphosulfate reductase (Meyer and Kuever, 2007), and sulfide quinone oxidoreductase enzyme systems (Pham et al., 2008) play vital roles in sulfide oxidation, and among these, the Sox multi-enzyme system is considered a fundamental and primordial molecular mechanism for sulfur oxidation (Ghosh and Dam, 2009) that is widespread among the known SOB (Meyer et al., 2007). Moreover, the soxB gene, which encodes the SoxB subunit of the Sox enzyme system, has been widely employed to characterize the abundance and diversity of SOB in various environments (Kojima et al., 2014; Thomas et al., 2014; Tourna et al., 2014).

The diversity and community structure of SRB and SOB have, to date, been studied using denaturing gradient gel electrophoresis (DGGE) (Varon-Lopez et al., 2013; Yang et al., 2015), restriction fragment length polymorphism (Luo et al., 2011; Reed and Martiny, 2013), and clone library approaches (Yang et al., 2013; Purcell et al., 2014). However, obtaining comprehensive information about the highly diverse SRB and SOB communities, particularly SRB and SOB from the prime habitats of marine sediments, using these methods is challenging (Aoki et al., 2015; Zhang et al., 2016). High-throughput sequencing, which constitutes a powerful approach for achieving complete coverage of microbial communities, has been recently applied in the analysis of microbial community diversity and composition. For instance, high-throughput sequencing of the 16S rRNA gene revealed the presence of a number of novel bacterial groups (Zhu et al., 2013; Mahmoudi et al., 2015) and a high diversity of bacteria (Wang et al., 2015). Moreover, the community structure and diversity of functional microbes, particularly nitrogen-cycling microbes, such as ammonia-oxidizing bacteria, ammonia-oxidizing archaea and denitrifying microbes, have been investigated through high-throughput sequencing analyses of functional genes (Tago et al., 2014; Saarenheimo et al., 2015). However, only a few studies have evaluated the community structure and diversity of SRB and SOB using high-throughput sequencing (Meyer et al., 2016; Cui et al., 2017).

The East China Sea (ECS) is the largest marginal sea of the Northwest Pacific Ocean, with a vast continental area of 0.5 × 10¹² m². The Changjiang River (Yangtze River), which transports great quantities of freshwater, nutrients, organic matter, and other chemical elements into the ECS, has a marked effect on the sea (Milliman and Meade, 1983; Beardsley et al., 1985). The ECS is also influenced by the warm and oligotrophic Taiwan Warm Currents in the south and the Kuroshio Current (from the Western Equatorial Pacific) in the east (Jiao et al., 2005). More interestingly, the ECS has also experienced serious environmental problems, including eutrophication, harmful algal blooms, and chemical pollution (Shen et al., 2011; Cheng et al., 2012; Yu et al., 2013). Because of the strong influence of its riverine inputs, currents and environmental problems, the ECS exhibits complex physical and chemical conditions and has become an ideal study area for ecological investigations of the temporal and spatial dynamics of biota (Dang et al., 2008; Xiong et al., 2014; Zhang et al., 2014). Increasing research efforts have focused on microbial ecology in the ECS, particularly the microbial community structure and diversity of organisms associated with the sulfur cycle, such as SRB and SOB (Nie et al., 2009; Wu et al., 2009; He et al., 2015); however, because of the limitations of the applied methods (DGGE or clone libraries), only a few dominant species have been identified. In the present study, the 16S rRNA gene and two functional genes (dsrB and soxB) were combined in a high-throughput sequencing analysis, and this combination allowed for an in-depth analysis of the community structure and diversity of SRB and SOB in sediment cores from the ECS. In addition, quantitative PCR (qPCR) was performed to reveal the vertical distribution of SRB and SOB populations and provide insight into the SRB and SOB community characteristics.

In this study, a comparative analysis of the sulfur cycle-related microbial communities (SRB and SOB) in sediment cores was performed using an approach combining 16S rRNA, dsrB, and soxB gene high-throughput sequencing with qPCR analyses. The combination of these molecular techniques provided a detailed description of the SRB and SOB communities. The traditional approaches (e.g., DGGE and clone library), which have been widely applied for the detection of SRB and SOB in various environments, might notably underestimate their diversity. High-throughput sequencing based on functional genes (dsrB and soxB) can result in a greater number of sequences and almost complete coverage (Good’s coverage values were all greater than 99%), providing more comprehensive information regarding the diversity and community composition of the SRB and SOB communities. Thus, the diversity, richness and OTU numbers observed in our high-throughput sequencing analyses were notably greater than those obtained in previous studies that employed other molecular methods (Jiang et al., 2009; Yang et al., 2013; Yousuf et al., 2014; He et al., 2015). Further analysis indicated that rare species make important contributions to the greater richness and diversity of the SRB and SOB communities. For instance, as shown in Table 3A, the number of dominant OTUs in each sample was similar, ranging from 8 to 13 and representing 3.08% of OTUs (on average); however, these OTUs comprised more than 77% of the sequence abundance for the dsrB gene. Although the rare OTUs (relative abundance <0.01%) (Galand et al., 2009) accounted for 54.03% of OTUs (on average), they only comprised 0.84% (on average) of the sequence abundance. Similar results were obtained for the soxB gene (Table 3B).

The 16S rRNA sequencing analysis identified 235 and 191 OTUs as belonging to potential SRB and SOB, respectively, based on family species classification. Furthermore, these SRB OTUs were affiliated with 10 families, accounting for 7.63–19.62% (14.12% on average) of the total 16S rRNA gene sequences (Table 4A), and the SOB OTUs were affiliated with 12 families, accounting for 5.99–18.78% (12.53% on average) of the total 16S rRNA sequences (Table 4B). This finding is consistent with the results of a previous study in this region that found that numerous sequences were affiliated with potential SRB and SOB (Ye et al., 2016). Overall, these results confirm the high abundance of potential SOB and SRB groups in the ECS sediments and suggest the existence of an active sulfur cycle in this area, and these sulfur cycle-related bacterial communities have considerable potential for further exploration.

The dsrB gene sequencing analysis performed in this study identified two phyla (Proteobacteria and Firmicutes), and Deltaproteobacteria (within Proteobacteria) represented the dominant SRB. Our data were consistent with previous observations obtained from various environments, such as marsh sediments (Jemaneh et al., 2013), the Great Salt Lake (Kjeldsen et al., 2007), mangrove sediments (Varon-Lopez et al., 2013), wastewater treatment plants (Biswas et al., 2014), and marine and estuary sediments (Jiang et al., 2009; Wu et al., 2009; He et al., 2015), suggesting that Deltaproteobacteria have a wide adaptation range and play major roles in global sulfur cycling. The soxB gene sequencing results revealed that the SOB community was dominated by Alphaproteobacteria and unclassified groups, followed by Gammaproteobacteria and Betaproteobacteria. For comparison, Gammaproteobacteria and Epsilonproteobacteria are considered dominant Sox organisms in marine sediments (Sievert et al., 2008; Lenk et al., 2011; Akerman et al., 2013; Dyksma et al., 2016); however, no Epsilonproteobacteria groups were observed in the present study. This discrepancy could be due to primer bias; the primers used for soxB gene amplification amplify the soxB genes of Epsilonproteobacteria and Chloroflexi very poorly (Meyer and Kuever, 2007). In addition, Bradyrhizobiaceae within Alphaproteobacteria and Burkholderiaceae within Betaproteobacteria were the dominant families in our study area; however, previous research revealed that these families might adapt to oligotrophic sulfur environments, such as soil (Jung et al., 2005; Masuda et al., 2010) and rhizosphere soils (Anandham et al., 2008; Masuda et al., 2016). Therefore, further study is required to determine whether the families described here are the result of contamination or actually exist. More interestingly, most of the domain and core OTUs obtained through dsrB and soxB gene high-throughput sequencing were difficult to clearly assign at the family or genus level, partially due to the limited number of reference sequences. For instance, 15 OTUs for the dsrB gene, with a relative abundance ranging from 40.91 to 73.36%, and 238 OTUs for the soxB gene, with a relative abundance ranging from 20.05 to 45.82%, could not be assigned at the family or genus level. Overall, these findings suggest that the ECS environment contains a high abundance of yet unknown SRB and SOB lineages, which should be studied in more detail. Indeed, this phenomenon is consistent with the following recent findings: 62% of dsrB sequences in surface sediments from the ECS exhibit low similarity with previously cultured SRB (Zhang et al., 2016); 57.6% of dsrAB clones from polluted harbor sediments belong to unclassified groups of SRB (Zhang et al., 2008); and 60% of sequences obtained through a dsrA gene analysis show no clear affiliation with a known SRB (Quillet et al., 2012). Furthermore, similar results have been reported by Yousuf et al. (2014) for soxB. In addition, we noticed that these unclassified OTUs exhibit high similarity with environmental dsrB/soxB sequences, as demonstrated through BLAST identity searches. For example, most of the unclassified dsrB OTUs displayed high similarities (>90%) with uncultured SRB recovered from Japanese fish farm sediments (Kondo et al., 2012) and ECS sediments (He et al., 2015; Zhang et al., 2016). Similarly, most of the unclassified soxB sequences exhibited high similarity (>80%) with sequences recovered from coastal soil ecosystems (Yousuf et al., 2014). These observations suggest that these unclassified sequences are widely distributed in marine and terrestrial environments and play a vital role in the biogeochemical cycling of sulfur and carbon.

A comparison of the SRB community composition identified based on the 16S rRNA and dsrB genes showed that 16S rRNA gene sequencing detected 10 SRB families, whereas dsrB gene sequencing only detected four SRB families (Tables 4A,B). Moreover, it was difficult to obtain more detailed taxonomic information for the SRB community based on the 16S rRNA and dsrB genes. Overall, the dsrB gene proved to be less efficient for the proper annotation of taxonomic information for SRB groups. This phenomenon might be attributed to the limited number of reference sequences used to annotate taxonomic information for the dsrB gene. In contrast, the results obtained from the high-throughput sequencing of the soxB gene revealed 15 SOB families, whereas the high-throughput sequencing of the 16S rRNA gene revealed 12 SOB families. In addition, the soxB gene can provide more detailed taxonomic information (at the genus or even the species level) for the SOB community (data not shown) than the 16S rRNA gene. Taken together, the results indicate that the soxB gene might better define SOB groups.

The 16S rRNA, dsrB and soxB gene high-throughput sequencing data obtained in the present study revealed a high diversity of SRB and SOB in sediments from the ECS, and the classification of SRB and SOB communities based on 16S rRNA gene and functional gene was analyzed and compared according to the results of a previous study (Meyer et al., 2016). Overall, the present study provided more comprehensive information regarding the SRB and SOB community characteristics and suggested the existence of an active sulfur cycle in the study area. Significantly, the sulfate amount never decreases below approximately 24 mM in the core samples, which appears to contradict the conclusion that SRB and SOB communities are responsible for active sulfur cycling. In fact, this phenomenon will be observed based on the fact that the sulfate consumed by sulfate reduction is recycled within sediments by oxidation of reduced inorganic sulfur compounds to sulfate. As previously reported, 80–95% of the massive amount of hydrogen sulfide formed through sulfate reduction is recycled within sediments and gradually oxidized back to sulfate (Jørgensen et al., 1990; Jørgensen and Nelson, 2004). Jørgensen et al. (1990) reported that sulfate reduction rates increased in the top 10 cm of core in the Belt Sea; however, the concentration of sulfate varied little (approximately 25 mM). Similarly, Leloup et al. (2009) revealed that substantial sulfate reduction rates were measured in the top 25 cm of core, and the sulfate concentration profile was approximately 20 mM. However, the diversity of sulfur cycle-related bacteria and the potential sulfur cycle activity might still be overestimated or underestimated because certain named SRB/SOB might did not actually take part in the sulfur cycle in sediments, such as Pelotomaculum and Sporotomaculum, which are not able to grow with sulfite and/or sulfate as electron acceptors even though they are closely affiliated with SRB (Brauman et al., 1998; Imachi et al., 2006). Instead, certain unnamed SRB/SOB, such as archaeal anaerobic methanotrophs, might be involved in the sulfur cycle in sediments (Treude et al., 2014). Fortunately, these groups are only found in special habitats and account for a small portion of the identified communities (Brauman et al., 1998; Imachi et al., 2006; Miyashita et al., 2009; Treude et al., 2014). Thus, the use of 16S rRNA and functional genes (dsrB and soxB) as proxies for determining the microbial groups that play active roles in sulfur cycling and for studying the diversity of SOB and SRB is a powerful approach.

The sediment environments of the ECS are highly heterogeneous (Dang et al., 2008; Fang et al., 2013); this characteristic has a strong influence on the spatial heterogeneity of the sediment microbial communities in the ECS, as observed in previous studies (Dang et al., 2008; Feng et al., 2009; Chen et al., 2014; Ye et al., 2016). Based on the sediment grain size combined with the overlying water masses, the ECS can be divided into three domains: the inner shelf mud area, the outer shelf sand area, and the slope plus Okinawa Trough mud area (Liu et al., 2006; Xu et al., 2016). S31 is located within the inner shelf mud deposits along the Zhejiang coast, whereas S33 is located at the periphery of the Zhejiang coastal mud area (Supplementary Figure S1). Thus, most of the environmental parameters showed substantial variation among the different samples and among the various sediment depths (Supplementary Table S1 and Figure S2). The concentration of DIN, for instance, ranged from 91.44 to 449.21 μM in station S31 and from 18.81 to 120.73 μM in station S33 (Supplementary Table S2). ²¹⁰Pb has been proven to serve as an effective proxy for sedimentary dynamics processes, such as erosion–transportation–deposition, boundary scavenging, and resuspension process, in coastal and shelf regions (Huh and Su, 1999; Su and Huh, 2002). The excess ²¹⁰Pb activity varied markedly between the two sediments (Supplementary Figure S2). At station S31, ²¹⁰Pb_ex activity was stable from the surface sediment to a depth of 8 cm; however, at a depth below 12 cm, ²¹⁰Pb_ex activity increased gradually and displayed characteristics of reverse accumulation, indicating that this sedimentary layer has experienced sedimentary dynamic events. At station S33, ²¹⁰Pb_ex activity remained stable until a depth of 5 cm and then gradually declined. This finding indicates that the sedimentary dynamic environment at S33 was more stable than that at S31, as observed in our previous study (He et al., 2015). Our results showed that the abundances of 16S rRNA, dsrB and soxB genes show substantial fluctuations at S31 (Figure 9), and none of the gene abundance values were significantly correlated with sediment depth. In contrast, all of the genes (16S rRNA, dsrB and soxB) gradually decreased with increasing depth at S33 (Figure 9) and were significantly negatively correlated with sediment depth. These results suggest that the sedimentary dynamic environment is likely an important factor in controlling the vertical distribution of the abundances of total bacterial 16S rRNA gene and functional genes (dsrB and soxB). Further analysis revealed similar vertical fluctuations in the abundance of the 16S rRNA, dsrB and soxB genes at S31, but the range of fluctuation obtained for the abundance of the 16S rRNA and soxB genes was smaller than that observed for the dsrB gene. In addition, the ratios of the functional genes (dsrB and soxB) to the total bacterial 16S rRNA gene were calculated: for the dsrB gene, the ratio ranged from 0.22 to 40.56% in S31 and from 0.06 to 4.2% in S33, and the ratios calculated for the soxB gene were less than 1% at both stations. These results indicate that the abundance of SRB might be more sensitive to the sedimentary dynamic environment than that of SOB and total bacteria.

The UPGMA clustering analysis results demonstrated different community structures of total bacteria, SRB, and SOB at the two sampling sites and among the different sediment depths, which might be due to the differences in the sediment depths and sedimentary dynamic environments between the two sampling stations. Overall, the sedimentary dynamic environment is most likely related to unmeasured variables, such as salinity, organic matter, dissolved oxygen, temperature, and H₂S, which are widely accepted as important factors in controlling the community structures of total bacteria and sulfur-cycling bacteria (SRB and SOB) (Headd and Engel, 2013; Ye et al., 2016; Zhang et al., 2016). In addition, an RDA analysis revealed that sediment depth and DIN have a significant influence on the community structure of total bacteria, SRB, and SOB. In summary, our results suggest that environmental parameters, such as sediment depth and DIN, and the sedimentary dynamic environment, which showed differences between the two sampling stations, can significantly influence the community structure of total bacteria, SRB, and SOB.

YZhe, TM, and ZY conceived and designed the experiments; YZha, XW, and HH performed the experiments and analyzed the data; YZha and YZhe wrote the paper.

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.