1 INTRODUCTION

Microbial communities are recognized as efficient players in the biogeochemical cycles of carbon, nitrogen, phosphorus, and sulfur (Fuhrman, 2009; Chaffron et al., 2010; Louvado et al., 2015). Investigations of microbial communities in sediments are therefore critical for understanding the marine ecosystem. A fundamental question is the bacterial community assembly, which is simultaneously regulated by stochastic processes (e.g., birth, death and dispersion) and the environmental conditions (Lennon et al., 2012; Coutinho et al., 2015). The environmental variables in a marine system are defined by the hydrological environment (Bi et al., 2011; Dixon et al., 2014).

As straits separate and connect two large bodies of water, they are influenced by strong tidal currents, transportation, and the exchange of various materials (e.g., nutritional chemicals and pollutants) (Li et al., 2015). Deep straits are characterized by stratification, oceanic currents and salinity/density differences on either side, while in shallow straits there are always large tidal currents forced by winds and oceanic inputs (Li et al., 2015). Intense bioturbation, high sedimentation rates, and low deposition of organicrich sediments can also be found in shallow straits (Behbahani et al., 2015). Coastal straits are also strongly influenced by human activities, which exacerbate the environmental complexity. Thus, coastal straits are representatives of complex hydrological environments. Investigations of the bacterial communities in coastal straits will therefore help to understand the microbial community overall in marine ecosystems.

Owing to the advancement of high-throughput DNA sequencing technologies, studies of complex marine microbial communities based on 16S rRNA genes have made great progress in recent years. In Fram Strait, located between Greenland and Svalbard at a depth greater than 3 500 m, marine bacteria in sediments exhibited distinct biogeographic patterns and depth was one of the important variables influencing bacterial diversity (Jacob et al., 2013). Moreover, interactive effects on the degradation activity of marine bacteria were induced by multiple environmental changes in the sediments in a deep strait (Piontek et al., 2015). However, few studies on the bacterial community in shallow straits have been reported.

The Bohai Strait (average depth around 20 m), which separates the Bohai Sea from the Yellow Sea, is adjacent to well-developed cities, Dalian and Yantai in the north and south, respectively, and classified as a shallow strait (Li et al., 2015). The Yellow Sea Warm Current (YSWC), which is a branch of the Tsushima Warm Current, flows into the Bohai Sea mainly in winter through the Bohai Strait. In addition, the Yellow Sea Cold Water Mass (YSCWM) that has cold and high salinity bottom waters partly presents in the northeastern Bohai Strait (Li et al., 2006). Therefore, the Bohai Strait is an ideal model to study coastal and shallow straits. Over the past few years, studies have investigated water circulation (Naimie et al., 2001; Li et al., 2015), suspended sediment transport (Bi et al., 2011), and zoobenthos communities (Mu et al., 2002). However, there are few studies on sediment microbial communities in this region. Thus, evaluating the roles of the sediment bacterial communities in the Bohai Strait is very important for understanding the microbial ecology in coastal and shallow straits overall.

So far, most studies on marine microbial communities have focused on their geographical diversity with emphasis on describing taxonomic composition, alpha/beta diversity, and the influences of environmental variables (Signori et al., 2014; Zheng et al., 2014). In addition to taxonomic composition, investigations on microbial associations should be included to understand the full assembly of bacterial communities (Faust and Raes, 2012; Coutinho et al., 2015). Those co-occurring bacteria may apply different strategies to avoid competition for resources in their environments (Macalady et al., 2008). The co-occurrence patterns among microorganisms can reveal the ecological associations and the key populations within the communities (Chaffron et al., 2010; Faust and Raes, 2012). Comprehensive analysis of the relationship between microbial association networks and related ecological functions will help to elucidate the roles of sediment microbial communities in straits.

In this study, six stations in the Bohai Strait were selected to explore microbial communities in the surface sediments in the summer of 2014. We investigated the taxonomic composition of microbial communities, their ecological associations, the key populations, and the environmental influences.

2 MATERIAL AND METHOD 2.1 Survey stations, sediment sampling, and physicochemical analysis

Six stations were sampled to cover the Bohai Strait (Fig. 1). The R and E sections represented two sides adjacent to the Bohai and Yellow seas, respectively. Between 27 August and 6 September 2014, composite surface sediment samples (0-5 cm) were collected in triplicate with a box grab at each station. The samples (500 g each) were pooled in sterile polythene bags and stored at -80℃ in the dark until processing. The physical parameters, including temperature, were determined in bottom water on-site using a Sea-Bird 911 plus CTD probe (Sea-Bird Electronics, Bellevue, WA, USA). Total carbon (TC), total nitrogen (TN), and mineral nitrogen (ammonium [NH₄⁺] and nitrate [NO₃^-]) in the samples were determined after freezedrying, as described by Guan et al. (2014).

2.2 DNA extraction, PCR, Illumina-Miseq sequencing, and sequence data processing

To extract genomic DNA, the MoBio PowerSoil^TM kit (MoBio, Carlsbad, CA, USA) was used according to the manufacturer's instructions by adding 0.5 g sediment to the kit. For each station, three replicates were extracted and then mixed. Six genomic DNA samples were obtained in total and then sent to Chengdu Institute of Biology, Chinese Academy of Sciences for sequencing analysis. The V4-V5 region of the bacterial 16S rRNA gene using the universal primers 515F (5'-GTGCCAGCMGCCGCGG-3') and 907R (5'-CCGTCAATTCMTTTRAGTTT-3') were sequenced on the Illumina-Miseq platform (Illumina Inc., San Diego, CA, USA). PEAR software (Zhang et al., 2014) was used to merge the paired reads with default settings and QIIME software with default parameters (Caporaso et al., 2010) was applied to process the merged sequences. Reads were filtered according to the quality scores and then split into samples according to the barcode information (split_ libraries.py). Chimeric sequences were identified and removed using UCHIME (identify_chimeric_seqs.py and filter_fasta.py) (Edgar et al., 2011). All samples without chimeric sequences were normalized to the same sequencing depth level at 18 000 sequences. Operational taxonomic units (OTUs) were clustered with a 90% similarity cutoff (Barberán et al., 2012) by using the strategy "pick_closed_reference" and UCLUST (Edgar, 2010). The phylogenetic affiliation of each 16S rRNA gene sequence was analyzed using Greengene (version 13.8, http://www.greengene.com/). Those OTUs unknown at the kingdom level were removed.

2.3 Statistical analysis

The coverage and alpha diversity of the samples were calculated in QIIME on the filtered OTUs. Alpha diversity was estimated using richness (Chao1) and diversity (Shannon) indices to represent community richness and diversity. The "vegan" packages (Dixon, 2003) in R were applied for representational difference analysis (RDA). To identify the metabolic or other functional capabilities of bacteria, a phylogenetic investigation of communities by reconstruction of unobserved states (PICRUSt) (Langille et al., 2013) was used. This predicted the functional composition of a metagenome by using marker gene data and a database of reference genomes. Following an online tutorial (http://picrust.github.io/picrust/index.html), metagenome predictions against KEGG Ortholog (Ogata et al., 1999; Okuda et al., 2008; Langille et al., 2013) were carried out based on the 16S rRNA sequence dataset. Then, thousands of predicted functions (KOs) were collapsed into higher categories (KEGG Pathways). The category "Energy Metabolism" was selected for further discussion of the ecological functions of microbial communities. Variance analysis was executed in R, with P < 0.05 considered to be statistically significant.

Network analysis was applied to visualize interactions within the sedimentary microbial community (Faust and Raes, 2012; Lupatini et al., 2014). Taxon assignment on lower taxonomic levels is not always accurate, which can confuse ecological signals and artifacts (Faust and Raes, 2012), so we only considered links and patterns above a certain taxonomic rank (Freilich et al., 2010). In this study, 647 (21.8%) and 1881 (63.5%) OTUs could be classified at the genus and family levels, respectively. Therefore, we considered ecological associations at family levels. In order to focus on the common bacteria in the survey region, those OTUs classified at family level and present in more than three samples were selected for the network analysis. Spearman correlations among the OTUs (R > 0.9, P < 0.05, adjusted with "Benjamini-Hochberg") were calculated in R to investigate the associations. After calculation, the associations were summarized at family level. The number of OTU associations (positive or negative) between two families was defined and termed NC, which is an estimate of the strength of the associations at the family level. Then, the positive and negative networks with weighted connections (NC values) were constructed using Cytoscape 3.2.1 (Smoot et al., 2011).

3 RESULT 3.1 Environmental conditions, bacterial communities and prediction metagenome

This study covered six stations in two regions (R and E) of the Bohai Strait in northern China (Fig. 1). The environmental variables are summarized in Table 1. The depths ranged from 15 m (station R7) to 50.8 m (E2), with temperatures from 7℃ (E2) to 24.94℃ (R7). The TC content in E2, R3, E6 and R7 were higher than in E4 and R5. A total of 41382 16S rRNA gene sequences were obtained across all samples after quality and OTU filtering. Good's coverage in all samples were greater than 90%, indicating that sequencing captured the majority of the bacterial diversity in the samples (Table 1). The average Shannon and Chao1 values in the R section (R3, R5 and R7, Mean_shannon=5.98±0.04, Mean_Chao1=2 385.37±351.90) were higher than in the E section (E2, E4 and E6, Mean_shannon=5.57±0.31, Mean_chao1=1 925.10±181.85). The highest diversity (Shannon=6.01) and richness (Chao1=2 650.77) values were found in R3 (Table 2).

In total, 2962 OTUs affiliated to 140 classes with 204 families were identified. At the class level, the four most abundant classes were Gammaproteobacteria (19.34%±5.02%), Deltaproteobacteria (16.72%± 5.06%), Bacilli (12.07%±2.05%), and Flavobacteriia (11.08%±5.07%). However, the relative abundance of the Flavobacteriia was lower in R3 (3.14%) than in other samples (Fig. 2a). At the family level, Lactobacillaceae (10.55%±1.89%) and Flavobacteriaceae (10.46%±4.94%) were the most abundant families. Psychromonadaceae represented a higher relative abundance in E2 (9.02%) than in other samples (Fig. 2b).

RDA was applied to investigate the major environmental influences in the Bohai Strait. In total, only 66.6% of the variation could be explained by the first two dimensionalities (Fig. 3). For the first dimensionality, temperature (R=-0.667 4, P < 0.05), TC (R=0.641, P < 0.05) and nitrate (R=-0.918 8, P < 0.05) exhibited significant correlations with RDA1. For RDA2, depth (R=0.749 6, P < 0.05) and TC (R=0.621 6, P < 0.05) were significantly correlated. Twenty-seven summary functions (e.g., energy metabolism, xenobiotic biodegradation, and metabolism), including 234 KEGG pathways (e.g., carbon fixation in prokaryotes), were predicted with PICRUSt. According to the variation of the KEGG pathways (Fig. 4a), Paths 2 (carbon fixation in prokaryotes), 3 (methane metabolism), 5 (oxidative phosphorylation) and 9 (sulfur metabolism) were relatively varied among samples. Nitrate (R=-0.948 8, P < 0.05) and TC (R=0.77, P < 0.05) were significantly correlated with RDA2 (Fig. 4b).

3.2 Family level microbial associations

The positive (P-network) and negative (N-network) networks were constructed to demonstrate microbial associations at the family level. The P-network consisted of 115 nodes and 962 edges (Fig. 5a) while 108 nodes and 750 edges were found in the N-network (Fig. 6a). The heterogeneity of both networks were severe. To reduce the complexity of the networks, a strong association was defined as a connection with a NC value greater than 5 and then two new networks with strong associations were extracted (Figs. 5b, 6b). The extracted positive network (EP-network) was comprised of 51 nodes and 118 edges (Fig. 5b). In the EP-network, the five most connected families were Flavobacteriaceae (27 connections), Piscirickettsiaceae (23 connections), Pirellulaceae (20 connections), Desulfobulbaceae (17 connections), and OM60 (16 connections). The four strongest associations were found between Piscirickettsiaceae and Desulfobulbaceae (NC=42), Flavobacteriaceae and Desulfobulbaceae (NC=35), OM60 and Flavobacteriaceae (NC=34), and Desulfobacteraceae and Flavobacteriaceae (NC=33). In the extracted negative network (EN-network), Flavobacteriaceae (32 connections), Piscirickettsiaceae (13 connections), Chitinophagaceae (6 connections), Desulfobacteraceae (6 connections), and Desulfobulbaceae (6 connections) were the five most connected families (Fig. 6b).

4 DISCUSSION 4.1 Influences on the microbial communities in the survey region

One fundamental question was how the environmental variables shape the structure of the microbial communities. The average Chao1 and Shannon values were both higher in the R region of the Bohai Strait sampled in this study (Table 2). Additionally, both of the Shannon and Chao1 indices were higher in the Bohai Sea than reported from the northern Yellow Sea (Liu et al., 2015). This suggests that the adjacent seas influence the microbial communities in the Bohai Strait. The niche theory postulates that the abundance and composition of microbial communities are determined entirely by the environmental conditions (Nemergut et al., 2013; Coutinho et al., 2015). In previous studies, different environmental factors affecting microbial communities were observed in different survey regions in the Bohai Sea. The median grain size and dissolved oxygen were the most important in intertidal sediments (Wang et al., 2013), while TN and total organic carbon were the most important in coastal sediments (Guan et al., 2014; Wang et al., 2015). However, temperature was always important (Wang et al., 2013, 2015; Liu et al., 2015). In this study, temperature, TC, depth and NO₃^- were the major factors responsible for the dissimilarity of the sediment microbial communities in the survey region (Fig. 3). Similarly, comparisons of Bohai Bay (Guan et al., 2014), Laizhou Bay (Wang et al., 2014), and Liaodong Bay (Zheng et al., 2014) indicated that temperature and nitrate were important parameters affecting sediment bacterial abundance and diversity in the Bohai Sea.

Another question was about ecological function and the related microbes. According to the KEGG pathways (Kanehisa et al., 2017), Path 2 (carbon fixation in prokaryotes) was known to exist in autotrophic bacteria and archaea including: (1) microaerophiles and anaerobes such as green sulfur bacteria, (2) strictly anaerobic bacteria and archaea (Proteobacteria, Planctomycetes, Spirochaetes, and Euryarchaeota), (3) some green non-sulfur bacteria of the family Chloroflexaceae, and (4) aerobic Crenarchaeota, Acidianus, Metallosphaera, and Sulfolobales. Path 3 is principally implemented by methanotrophs and methanogens. In this study, the class of Ignavibacteria with the highest relative abundance (0.54%, data not shown) in R3 was linked to Path 2. Path 3-related microbes might consist of Methanobacteria and Methanomicrobia and the family Methanomassiliicoccaceae. Nitrate and TC content can explain the variations in the nine KEGG pathways of energy metabolism (Fig. 4b). Moreover, terrigenous nutrient input (mainly nitrogen from the wide use of nitrogen-based fertilizer) and material exchange between the Bohai and Yellow seas may lead to different nutrient distributions (Zhang et al., 2004), indicating that variations in environmental variables induced by the complex hydrological conditions will exert influences on the microbial communities and their ecological functions in the Bohai Strait.

Microbes with high abundance are known to exert great influence on the structure and ecological function of the overall microbial community (Campbell and Kirchman, 2013, Zhang et al., 2013). Therefore, the shifts in the relative abundance of predominant microbes under different environmental conditions can reflect the adaptations of microbial communities to habitat filtering. The most abundant OTUs (8.33%±1.51%, data not shown) in the family Lactobacillaceae were affiliated with Lactococcus piscium, first isolation from salmonid fish (Sakala et al., 2002). Piscirickettsiaceae are generally considered to be pathogens causing rickettsial disease in fish and have high nutritional needs (Fryer and Hedrick, 2003); these were particularly abundant at the R7 stations (8.51%, Fig. 2b) where a fish farm is located. This implies that fishery breeding affects the microbial communities in surface sediments in the Bohai Strait. Psychromonadaceae only contains the genus Psychromonas and is usually considered a psychrophile with widespread distribution in deep seas, such as the Japan Trench, and Antarctic environments (Maruyama et al., 2000). The relative abundance of this family was higher in the E2 stations (9.02%, Fig. 2b) that have much lower temperatures (Table 1), indicating that the YSCWM influences the bacterial communities in the Bohai Strait. Desulfobacteraceae and Desulfobulbaceae were highly abundant (3.33%±1.21% and 4.85%±0.51%, respectively) across all samples, except at station E2 (Fig. 2b). These families have been previously found in marine sediments and mainly consist of sulfate reducers that completely oxidize organic matter (Elsabé et al., 2012; Na et al., 2015). The lowest relative abundances of the two families in E2 (1.72% and 1.59%, Fig. 2b) suggest that they are intolerant of cold water (Korneeva et al., 2015).

4.2 Key microbial families

We studied the interactions among different taxa within sediment microbial communities. Bacteria form complex ecological webs (Faust and Raes, 2012), and the associations are important elements for the structure of microbial communities. The interaction types encompass positive, negative, and neutral impacts. For instance, syntrophy between two species results in benefits to both, a positive impact, while competition between two species results in a negative impact (Hibbing et al., 2010; Molloy, 2014). These interactions can be reflected in their abundance correlations (Faust and Raes, 2012). Numerous studies on microbial biogeography (Campbell and Kirchman, 2013; Jacob et al., 2013; Buttigieg and Ramette, 2015; Coveley et al., 2015; Piontek et al., 2015) and the niche theory (Pontarp et al., 2012) based on habitat filtering have indicated that the members with positive correlations respond similarly to environmental conditions in their habitats. That is, if two species experience a positive correlation in their abundance across multiple environments or samples, this can be interpreted as evidence that they occupy similar ecological niches (Coutinho et al., 2015). Based on these observations, significant correlations between two species were often related to the direct or indirect interactions between them (Barberán et al., 2012; Faust and Raes, 2012; Lupatini et al., 2014; Coutinho et al., 2015). Therefore, connections with high NC values (NC > 5) could explain the strong interactions between two families, for example, Piscirickettsiaceae-Desulfobulbaceae (in the EP-network: NC=42), FlavobacteriaceaeDesulfobulbaceae (NC=35) and OM60- Flavobacteriaceae (NC=34). The key families include Flavobacteriaceae, Pirellulaceae and Piscirickettsiaceae (Figs. 5, 6).

Meanwhile, both the positive and negative networks are scale-free networks preferred by most biological systems (Boccaletti et al., 2006). One distinctive characteristic of the scale-free network is high fault tolerance (Albert et al., 2000). In other words, the forming of an ecological network is crucial for the stability of a microbial ecosystem (Faust and Raes, 2012). On the one hand, the number of connections to one node is the node's degree, and nodes with high degree values are known as hub nodes, which reduce the robustness of the network (Albert et al., 2000; Boccaletti et al., 2006) because the networks are destroyed immediately when the hub nodes are removed. On the other hand, the microbial ecosystem rapidly becomes stable when exogenetic pollutants are introduced into their habitat. For example, our previous study showed that the relative abundances of most microbes under the stress of hydrocarbons with different structures become stable after 30-day incubation (Wang et al., 2016). In this study, the families with high degree values were Flavobacteriaceae, Piscirickettsiaceae and Rhodobacteraceae in both the EP-network and ENnetwork (Figs. 5b, 6b). Flavobacteriaceae represent a major branch of the gram-negative phylum and has undergone rapid expansion incorporating many new species and genera (Bowman, 2006). There are seven phylogenetically related genera in Piscirickettsiaceae, but these genera have very different phenotypic characteristics (Marshall et al., 2014). Rhodobacteraceae exhibits global distribution and about 100 genera are currently recognized (Pujalte et al., 2014). The three key families have high phylogenetic diversity and wide distributions. The high phylogenetic diversity may imply high overlap of the ecological functions within the community. If some microbes in the three families that have certain ecological functions were removed, other microbes could perform similar functions accordingly. Our findings suggest that an important strategy to maintain the stability of the microbial ecosystem is to select families with high phylogenetic diversity as key populations in the microbial association network.

5 CONCLUSION

In this study, surface sediment bacterial communities at six sites across the Bohai Strait were analyzed in August-September 2014. Investigations into the relationship between microbial communities and environmental variables showed that temperature, total carbon, depth and nitrate were the major influences on the communities. Based on network theory and analysis, the families that are crucial for the stability and functionality of the microbial ecosystem were Flavobacteriaceae, Pirellulaceae and Piscirickettsiaceae.

6 ACKNOWLEDGEMENT

The manuscript was kindly revised by Dr. Shiliang Anthony Liu from Louisiana State University.