60
IV. Effects of bottom-water hypoxia on sediment bacterial
61
influences on respiratory metabolism, and related gene expressions of the bacterial community (Jørgensen 2006).
Previous studies on spatiotemporal variation of planktonic bacterial community composition (BCC) in relation to the progression of coastal hypoxia have revealed that succession of BCC during a transition from oxic to anoxic conditions in water column was represented by taxon replacement (Crump et al. 2007; Zaikova et al. 2010; Spietz et al. 2015), while there are several bacterial lineages, such as members of SUP05, ARCTIC96BD-19 and SAR324, which are consistently found in various types of oxygen-deficient water (Zaikova et al. 2010; Wright et al. 2012; Ulloa et al. 2012; Forth et al. 2015).
In addition to the shift in BCC during hypoxia development, the diversity of the bacterial population is also influenced by the degree of oxygen depletion. Beman and Carolan (2013) reported a non-linear relationship between the dissolved oxygen (DO) concentration and bacterial richness in oxygen minimum zones (OMZ) of the Eastern tropical North Pacific Ocean. In contrast, Spietz et al. (2015) reported a strong negative association between DO and bacterial richness in the seasonally hypoxic estuary of Hood Canal in the US. Although the question of how deoxygenation of seawater would affect bacterial community dynamics has not yet been fully addressed, it is certain that
62
DO plays a fundamental role in shaping planktonic BCC in coastal waters.
In contrast to the growing concerns on the impacts of deoxygenation for planktonic BCC, there is currently limited knowledge about how sediment BCC would respond to changing levels of DO in the water overlying seafloor. Mahmoudi et al. (2015) found sediment BCC of Caspian Sea in permanently hypoxic regions was different from that in oxic regions. Similarly, Devereux et al. (2015) reported significant differences in sediment BCC between normoxic and hypoxic periods in the northern Gulf of Mexico.
Jessen et al. (2017) also found changes in relative abundance of some phylogenetically distinct groups of bacteria (e.g., Flavobacteriia, Gammaproteobacteria, and Deltaproteobacteria) in response to decreasing oxygenation in the Crimean shelf break
of the Black Sea. Although these pioneering works demonstrated distinct shifts in sediment BCC between oxic and hypoxic conditions, their findings are based on snapshot observations and the details of how oxygen availability affects sediment BCC are poorly constrained. In order to understand the responses of sediment BCC to spreading hypoxia in the coastal sea bottom more precisely, it is necessary to examine the BCC in a defined location for a longer period of time.
In light of this, a shallow enclosed sea that has been consistently impacted by seasonal hypoxia is likely to provide a suitable research environment to conduct
long-63
term monitoring of sediment bacterial community, as it would provide a wide range of different DO conditions. Omura Bay in Nagasaki prefecture, western Kyushu, Japan is a strictly enclosed bay that represents an ideal site for such long-term monitoring (Fig.
1). The bay is connected to the open sea (East China Sea) by only two narrow channels through which cold, dense oxygenated water flows into the bottom of the bay in winter, whereas warmer and less dense water intrudes into the middle layer of the bay in early summer, leaving the bottom water below the intruding depth stagnant, and consequently leading to hypoxia (Nogami et al. 2000; Takahashi et al. 2009). Determination how the seasonal hypoxia develops and affects secondary or higher productivity within the bay is a relatively long-standing topic (Mori et al. 1973; Yokoyama 1995a; b); however, studies of the sediment microbial population during the formation and maintenance of hypoxia have just begun (Wada et al. 2012; Mori et al. 2015). Wada et al. (2012) revealed that upper sediment in the central bay was characterized by greater community respiration and diversity of bacterial components compared with the non-hypoxic sediment of the bay fringe during bottom water hypoxia. Mori et al. (2015) further revealed that potential community respiration of surface sediment was higher during the hypoxic period than during the normoxic period. In order to understand the BCC dynamics of the bay in more detail, I consider how diversity, richness and community
64
structure of sediment bacterial population shift in response to availability of DO in the water overlying the sediment.
I used a community fingerprinting method, automated ribosomal intergenic spacer analysis (ARISA), to assess bacterial community dynamics (Fisher and Triplett 1999).
Based on the combination of ARISA and Sanger sequencing of the bacterial 16SrRNA genes franked with ITS, I delineated the shift in the sediment bacterial community and demonstrated changes in their diversity in response to DO availability in the bottom water of Omura Bay. Based on these results, I discuss the significance of the findings in the context of the recent trend in coastal deoxygenation around the world.
Materials and methods Study site and sampling
Surface sediment samples and environmental parameters in the water column were collected in a center region of Omura Bay (St. 21, 32°55.390ʹN, 129°51.350ʹE, Fig. 1).
During summer months from 2011–2013 (Table 5), three replicate sediment cores were collected at each sampling from the center location at 20 ± 1 m depth with an acrylic pipe (31 cm long with 26 mm inner diameter) by scuba diving carefully avoiding the area with visible bioturbation. The DO, temperature, salinity and chlorophyll a of vertical profiles
65
in the sampling sites were obtained using a Conductivity Temperature Depth (CTD) profiler (AAQ, JFE-Advantec Co, Kobe, Japan).
All core samples were kept at in situ temperature and carefully brought to the laboratory within 3 h after sampling, during which time samples were handled carefully to avoid direct exposure to sunlight and other physical disturbances. Upon return to the laboratory, sediment was extruded from cores down to either 14 or 10 mm depth from the top. The top 0–7 mm (in 2011) or 0–5 mm (in 2012 and 2013) layer of sediment was pooled and regarded as the uppermost sediment layer, whereas the sediment within the depth range of 7–14 mm (in 2011) or 5–10 mm (in 2012 and 2013) was regarded as the subsurface layer. The pooled sediment layers from three cores were stored -20°C until DNA extraction and total organic carbon (TOC) measurement. Bacterial counting and TOC measurement were conducted as described in chapter III.
It was unfortunately not possible to keep the depth range of the uppermost surface and subsurface layers constant throughout the study, mainly because different devices were used to dissect the cores (Wada et al. 2012). Nevertheless, the uppermost sediment layers (0–7 mm) in 2011 mostly overlapped with those in the other 2 years (0–5 mm), and the number of cores in 2011 contributed only 30% of the total number (6 out of 20) of cores examined. Therefore, it would seem reasonable to assume that combining the 2011
66
data with those of 2012 and 2013 would not be significantly different from combining data among sediment layers of identical depth (0–5 mm). For the same reason, I also combined data from the uppermost and the subsurface layers of 2011 with those of 2012 and 2013 for analyzing impacts of DO availability on the microbial community in the surface sediment. However, when comparing the community structure and diversity between the two sediment layers, I combined and analyzed the data from 2012 and 2013 only to avoid possible complications arising from inclusion of 2011 samples.
Bacterial community analysis using ARISA
The DNA from each sediment layer was extracted using an ISOIL DNA extraction kit (NIPPON GENE, Tokyo, Japan), according to the manufacturer’s instructions. The extracted DNA was subjected to PCR amplification of the bacterial ribosomal ITS regions (rITS) using the universal primer sets (Cardinale et al. 2004): ITSF (5’-GTCGTAACAAGGTAGCCGTA-3’) with a 5’-end labeled with a 6-FAM, and ITSR-eub (5’-GCCAAGGCATCCACC-3’). The PCR was carried out in 12.5 μL reaction mixture using PCR buffer, 5U of KOD FX Neo polymerase (Toyobo, Osaka, Japan), 4 nmol of dNTP mix, 2 pmol of each primer and 10 ng of template DNA. The reaction was performed in a thermal cycler (DNA engine PTC-200, Bio-Rad, USA) using the following
67
program: 1 cycle at 94°C for 2 min, followed by 30 cycles of 98°C for 10 s, 59°C for 30 s and 68°C for 60 s. Products were then run on an ABI PRIZM 3730 automated sequencer operating as a fragment analyzer with 1200LIZ marker (Applied Biosystems, Foster City, CA, USA). The electropherograms were then analyzed using the Peak Scanner software (v1.0; Applied Biosystems) where each peak displayed in an electropherogram represented one operational taxonomic unit (OTU) and the peak height represented the amount of the OTU present. I used the cutoffs of the fragments described by Popa et al.
(2009) with slight modifications as follows: we eliminated (1) all fragments smaller than 82 bp, which was the sum between the size of the primers and the remaining parts of the 16S rRNA and 23S rRNA genes; (2) all fragments larger than 1000 bp; and (3) all fragments with fluorescent intensity less than 0.5% from the total fluorescence.
Subsequently, dynamic binning was performed for ARISA fragments (82–1000 bp) as described by Ruan et al. (2006). Maximum bin widths for fragments 82–700 and 701–
1000 bp were 3 and 5 bp, respectively. The relative abundance of the respective OTUs to the community was estimated as peak height divided by the cumulative peak height of a given sample.
68
Linking ARISA fragments to phylogenetic affiliations
In order to sequence some distinctive OTUs of ARISA fragments, I constructed clone libraries for the 16S-ITS region. Bacterial 16S-ITS regions were amplified from mixed genomic samples using PCR with the universal bacterial primer 1055f (5’-ATGGCTGTCGTCAGCT-3’) (Lane,1991) and ITSR-eub (5’-GCCAAGGCATCCACC-3’) (bacterial specific ITS) with the same PCR programs as described above. After purification (QIAquick PCR purification kit; Qiagen, Germany), the amplified genes were cloned into a plasmid vector pTAC-1 of a TA PCR cloning kit (BioDynamics Laboratory, Tokyo, Japan) according to the manufacturer’s instruction and transformed into competent Escherichia coli strain DH5α (BioDynamics Laboratory Tokyo, Japan).
The 16Sr-ITS gene sequences were determined with two primers, M13 BDFw
(5’-CAGGGTTTTCCCAGTCACGAC-3’) and M13 BDRev Primer (5’-
CGGATAACAATTTCACACAGG -3’), using a BigDye Terminator V3.1 Cycle Sequencing Kit (Applied Biosystems) and an ABI3730 PRISM Genetic Analyzer.
Sequence data were assembled into contiguous (contig) sequence using SeqMan (DNASTAR Inc., Madison, WI, USA). Phylogenetic affiliation of the clone was inferred by BLAST search on NCBI websites. The ITS region of each clone was reamplified with ITSF and ITSR-eub, electrophoresed and its length was estimated in the same way as
69
described above. If the length was different by more than 20 bp from that estimated by manual counting, I excluded the corresponding clone from the downstream analysis.
Nucleotide sequence accession numbers
All the sequence data obtained in this study have been deposited in the in DDBJ (accession nos. LC145643–LC145649, LC145651, LC145652, LC145654, LC145656, LC145657, LC145662, LC145664, LC145665, LC145669–LC145672, LC145675–
LC145682, LC145684–LC145688, LC145690, and LC145694–LC145696).
Statistical analysis and diversity
As previously stated, I combined the 2011 data with the 2012 and 2013 data when comparing the community structure and diversity across the bottom-water DO conditions.
However, when comparing between the uppermost and subsurface sediment layers, I combined and analyzed the 2012 and 2013 data only in order to avoid possible complications arising from the inclusion of 2011 samples.
To assess similarity between the bacterial communities, pairwise similarity matrices were calculated from relative OTU abundance data using the Bray–Curtis similarity and visualized using a nonparametric multidimensional scaling plot (nMDS) Analysis of
70
similarity (ANOSIM) was used to evaluate patterns of bacterial community similarity among groups of samples (Clarke and Warwick 2001). Similarity percentage analysis (SIMPER) was used to determine the OTUs that contributed most to similarity within, and dissimilarity among bottom-water DO groups (Clarke and Warwick 2001). Pearson rank correlation among the best combination of geochemical variables and Bray–Curtis pairwise resemblance (BIONENV) were calculated to access the relationship between community composition and z-transformed environmental variables (Clarke and Warwick 2001). Shannon (H’) and Simpson diversity (1-λ) indices were calculated to access both richness and evenness. All of these data analyses were performed with PRIMER 6 (PRIMER-E Ltd, Plymouth, UK) (Clarke and Gorley 2006). I also used a backward stepwise regression to select the most predictive variables for cells abundance and diversity indices of each sediment layer, which were performed using the statistical program R (R Core Team 2015). Environmental data were log10-transformed before performing the regression analysis. The Akaike information criterion (AIC), which balances the fit of a model against the number of parameters, was used to select the best fit model. In addition, a generalized linear model and quadratic regression models were compared to examine the relationship between environmental parameters, bacterial abundance and richness. Outliers in the data set were detected using Grubbs’ outlier test
71
at the 5% significance level before regression analysis.
Illumina MiSeq 16S rRNA Gene Amplicon Sequencing
Six genomic samples were taken from the uppermost sediment layer from June through December in 2011, which were amplified for Illumina MiSeq 16S rRNA Gene Amplicon Sequencing using a forward and reverse fusion primers. The primers specific for V3-V4 region of 16S rRNA gene which also incorporates the Illumina overhang adaptor 341F (5'- ACACTCTTTCCCTACACGACGCTCTTCCGATCT-CCTACGGGN GGCWGCAG-3') and 805R (5'- GTGACTGGAGTTCAGACGTGTGCTCTTCCGATC T-GACTACHVGGGTATCTAATCC-3') were used. The PCR was carried out in 25 μL reaction mixture using PCR buffer, 5 μl forward primer (1 μM), 5 μl reverse primer (1 μM), 12.5 μl 2X Kapa HiFi HotStart ReadyMix (Kapa Biosystems, Boston, MA, USA) and 10 ng of the extracted DNA. The program for 27 cycles of PCR had an initial step of 95°C for 3 min, followed by 95°C for 30 s, 55°C for 30 s, 72°C for 30 s, and finally 72°C for 5 min. The correct size of PCR products was verified by agarose gel electrophoresis, and the PCR products were further purified (Agencourt AMPure XP; Beckman Coulter, Indianapolis, USA). A second PCR reaction was performed on the purified PCR products (2.5μl) to index each of the samples. Two indexing primers (Illumina Nextera XT
72
indexing primers, Illumina, Sweden) were used per sample. The second PCR was carried out in 25 μL reaction mixture using PCR buffer, 5 μl forward index primer (1 μM), 5 μl reverse index primer (1 μM), 12.5 μl 22X Kapa HiFi HotStart ReadyMix (Kapa Biosystems, Boston, MA, USA), and 2.5μl of the first PCR products. PCRs were completed as described above, but only 8 amplification cycles were completed. The second PCR products were purified again (Agencourt AMPure XP; Beckman Coulter, Indianapolis, USA). After purification, the same volume and concentration of PCR amplicons were mixed every eight samples and sent to Bioengineering Lab. Co., Ltd. for the massively parallel sequence using an Illumina Miseq.
Processing and quality control of reads was performed using the R package dada2, version 1.6 (Callahan et al. 2016). Finally, amplicon sequence variant (ASVs), which are higher-resolution version of the OTU (see Callahan et al., 2017 for more detail) , were classified from the kingdom to the genus level using the Silva reference 16S rRNA gene database, version 132 resulting in the construction of an ASV table with read counts of all ASVs in all samples. All associated raw data is available at the DDBJ Sequence Read Archive under accession numbers DRA007315 (DRX140792-DRX140797).
73 Results
Environmental variables and bacterial abundance in overlying water and sediment
Figure 6 shows changes in dissolved oxygen (DO), temperature and salinity in the bottom water overlying the sediment (average value from 0 to 50 cm above sediment) at the center of Omura Bay during the field surveys in the years 2011–2013. In this study, the DO concentration range was divided into three categories according to the criteria used by Wright et al. (2012) with slight modifications: oxic (>90 μM O2), dysoxic (20–
90 μM O2) and suboxic or anoxic (<20 μM O2). Water temperature at the bottom ranged from 14 to 27.3°C and peaked in September each year. Salinity at the bottom ranged from 28.5 to 33.2, with a decreasing trend during summer months (June–September).
The TOC content ranged between 22.0-39.5 mg g−1 (dry weight sediment) in the uppermost layer and 33.6-39.7 mg g−1 in the subsurface layer of the sediment. The TOC content in the subsurface was significantly higher than that in the uppermost layer (Welch two-sample t-test; P < 0.05). The TOC content in the uppermost layer remained relatively high from June to late August, started to decline from late August to September, and reached minimum by the end of September each year. In contrast, the TOC content in the subsurface layer did not show such a trend (Table 5). However, there was a significant
74
correlation between TOC content in two layers (r = 0.70, P < 0.01). Relationship between TOC and DO in each sediment layer is shown in Figure 9. The TOC content of the surface sediment tended to increase under suboxic conditions.
Bacterial abundance in the uppermost layer ranged from 4.3 × 109 to 1.7 × 1011 cells g−1 and that in the subsurface from 1.6 × 109 to 9.3 × 1010 cells g−1, respectively. Bacterial abundance also covaried between the two layers (r = 0.98, P < 0.01). However, there was no clear change in bacterial abundance in response to DO availability. Multiple linear regression analysis was further used to demonstrate relationships between bacterial abundance and environmental variables. This analysis revealed that DO, temperature, salinity and TOC in combination explained more than half of the variability of bacterial abundance in the sediment layers (Table 6). Although any one of those factors alone poorly explained the variability of bacterial abundance across the sediment layers, TOC exhibited a weak linear relationship with the bacterial abundance (adjusted R2 = 0.19, P
= 0.05 for uppermost sediment, adjusted R2 = 0.17, P = 0.06 for subsurface sediment; see Supplementary Fig. A1).
Diversity of the bacterial community and its correlation with DO and
75 other environmental variables
A total of 164 unique OTUs were found in the surface sediment of the center of Omura Bay over the three consecutive years (2011–2013). The uppermost sediment layer harbored less OTUs (50.0 ± 8.1, average ± standard error) than the subsurface layer (58.1
± 8.4) across the studied DO conditions (Welch two-sample t-test; P < 0.01) (Table 5).
Comparison of the OTUs in each sediment layer between different DO availabilities revealed that 32 to 37% of the OTUs were shared between oxic and lower DO (dysoxic and suboxic) conditions, whereas a greater portion of the OTUs (39% for uppermost, and 47% for subsurface layer) were shared between dysoxic and suboxic conditions (Fig. 10).
In contrast, comparison between the two sediment layers showed the percentage of shared OTUs was highest under oxic conditions (on average 49.4%) and decreased towards less oxic (on average 38.7–40.5%) conditions (Fig. 10). However, no significant differences were detected among oxygen groups (Kruskal–Wallis, P = 0.45).
Multiple linear regression analysis between diversity indices (OTU number, H’ and 1-λ) and environmental variables revealed that DO alone would explain 19% of variation of OTU number and H’ in uppermost layer (Table 7). Though salinity and TOC were also identified as additional explanatory parameters for 1-λ, their relationship was insignificant (Table 7). It was further demonstrated that DO alone can explain partly, but significantly
76
the variation (Fig. 11), and that a quadratic regression model resulted in a unimodal pattern with greater R2 value than that achieved by a linear regression model. In all cases, the diversity indices peaked at around 11 µM O2. In contrast, other factors poorly explained the variability in diversity indices in the uppermost layer (see supplementary Fig. A2). However, combinations of DO and temperature could explain 27–38% of variation of H’ and 1-λ in the subsurface layer (Table 7). A linear and quadratic regression model analysis between diversity indices in the subsurface layer and DO were conducted but no significant relation was detected (Fig. 11).
Shifts in bacterial community of the surface sediment in response to DO availability at the center of Omura Bay were further depicted by nMDS based on Bray–Curtis similarity (Fig. 12). Separation of the bacterial community structure in the uppermost layer was significant when suboxic condition was compared with more oxygenated conditions (ANOSIM, R = 0.34, P < 0.05 for oxic vs suboxic, R = 0.324, P < 0.05 for dysoxic vs suboxic, Table 8), while no significant separation among oxygen groups was detected in the subsurface layer (Table 8a). Mean similarities of the microbial community between samples under dysoxic and suboxic conditions were greater than the mean similarity observed among the samples under oxic conditions (Fig. 13). Shifts in the bacterial community structure between the two sediment layers also became significant
77
only under suboxic conditions (ANOSIM, R = 0.327, P < 0.01, Table 8b). In support of this, similarity between bacterial communities of the two sediment layers tended to decrease towards less oxic conditions; the similarity was 60.2% on average under oxic conditions, which decreased to 57.3% and 55.4% under dysoxic and suboxic conditions, respectively (Fig. 13c).
I used BIOENV to examine the extent to which DO availability and other variables can explain the shifts in the bacterial community structure. I found that DO correlated best with the bacterial community structure of the uppermost sediment layers (ρ = 0.322, P < 0.01), whereas TOC alone was the best predictor for the subsurface bacterial community shifts (ρ = 0.327, P < 0.01) (Table 9).
Changes in relative abundance of the individual OTUs in response to DO availability
Figure 14 shows the rank abundance plots of the first 10 OTUs making up to 36–45%
cumulative contribution to the total abundance in upper and subsurface sediment layers under different DO conditions. An ITS fragment of 643 bp (“OTU-643”) was predominant in both sediment layers under all studied DO conditions. While more than half of the OTUs were consistently found within the rank positions irrespective of the DO
78
availability, their ranking order varied noticeably. For instance, “OTU-515” was ranked within the first four ranks under oxic conditions but was excluded under dysoxic and suboxic conditions, whereas “OTU-654” was ranked fifth under oxic condition and second or third under lower DO conditions. This observation prompted us to examine how each OTU changed its relative abundance in response to DO availability. I found that 87 out of the 164 unique OTUs became more abundant under higher DO (the OTU was presumably aerobic), of which 6% (5 out of 87) showed statistically significant positive coefficients (P < 0.05). However, 77 out of the 164 OTUs became abundant in lower DO condition and therefore had negative coefficients (the OTU was presumably anaerobic), of which 21% (16 out of 77) were significant (P < 0.05) (see supplementary Table A1).
Table 6 summarizes OTUs and the corresponding correlation coefficients in each sediment layer contributing up to 50% cumulative similarity under different DO conditions. For both layers, more than half (> 55%) of the OTUs under oxic conditions were positively correlated with DO, while less than half (27–33%) of OTUs under dysoxic and suboxic conditions showed positive correlation coefficients. I also noticed that nearly half of the OTUs (72 out of 164) showed contrasting correlation with DO between two sediment layers. For instance, “OTU-643” showed a positive correlation (r
79
= 0.4) in the upper sediment, whereas that in the lower sediment showed a negative correlation (r = −0.29) (Table A1).
I obtained a total of 56 clones for linking ARISA fragments to phylogenetic affiliations and assigned them into 26 distinctive phylotypes (OTUs) (see Supplementary Table A2). Sequencing the distinct 16SrDNA clones flanked with corresponding ITS fragments identified the predominant “OTU-643” as a close relative of Paraferrimonas sedimenticola, which belongs to Alteromonadales in Gammaproteobacteria. However,
due to our limited cloning effort, the sequence data were not used to describe shifts in phylogenetic compositions of the sediment bacteria, but were instead used to test whether the phylogenetic affiliation of OTUs might be reflected in the correlation between individual OTUs and DO availability. A total of 5 ARISA fragments (OTU-444, 510, 657, 662 and 686) were demonstrated to have multiple phylogenetic affiliations. This can happen in ARISA analysis as this technique is solely based on differences in the ITS fragment length. Accordingly, OTU-510 and OTU-657 each representing two different classes were not used for the following analysis. About 75% of OTUs (9 out of 12 clones) that were affiliated with Deltaproteobacteria (including members of Desulfobacterales and Desulfuromonadales) correlated negatively with DO. All of the OTUs (100%, 8 out of 8 clones) affiliated with Gammaproteobacteria (including members of
80
Alteromonadales, Aeromonadales, Oceanospirillales and Methylococcales) correlated positively with DO.
Relative abundance of the individual ASVs
Figure 15 shows relative abundance of bacterial groups in the uppermost sediment layer from June through December in 2011. The barplot shows the relative abundance of bacterial Class level (Fig. 15A) and the ten of most abundant bacterial Family (Fig. 15B).
Members of Gammaproteobacteria, Deltaproteobacteia and Bacteroidia predominated across the surface sediment. However, Gammaproteobacteria was less dominant in August (suboxic condition), while Deltaproteobacteria became dominant under decreasing DO conditions. Desulfobacteraceae was the most abundant bacterial family within all sediment samples (mean abundance; 10%), followed by Flavobacteriaceae (mean abundance; 6.7%) and Woeseiaceae (mean abundance; 6.4%).
Discussion
The primary goal of the present chapter was to gain quantitative insights into how seasonal bottom hypoxia would affect the bacterial community richness, diversity and composition in the surface sediment of Omura Bay, a typical enclosed bay in Japan. In a
81
previous study, Wada et al. (2012) revealed that the sediment microbial community respiration as measured by the reduction of tetrazolium salt, and the bacterial community structure as revealed by ARISA at the center of Omura Bay under persistent summer hypoxia were distinct from those in the reference site (the south fringe of the bay). The results seem to have reflected the interactions of the microbes with different DO regimes between the two sites in summer 2009. However, the spatio-temporal extent of hypoxia development and maintenance could change significantly within and between years (Nogami et al. 2000; Suzaki et al. 2013). It is therefore vital to conduct multi-year observations of hypoxia in the bay to infer seasonal and/or inter-annual responses of bacterial community to bottom hypoxia.
In the present study, I examined seasonal hypoxia in Omura Bay from 2011 to 2013, during which time I encountered a wide range of different patterns of hypoxia development and maintenance. There was a massive development of hypoxic water in 2012 in which I found transient but complete anoxia. No such development occurred in the other two years (2011 and 2013). Therefore, the intra- and inter-annual differences in hypoxia provide us a unique opportunity to examine impacts of a wide variety of DO availability on the microbial community in the surface sediment.
82
Effect of bottom-water hypoxia on organic matter preservation and bacterial abundance
It is well known that deposited organic matter is degraded by benthic fauna and microorganisms much more efficiently under aerobic conditions than under anaerobic conditions (Middelburg and Levin 2009; Jessen et al. 2017). Therefore, a substantial amount of organic matter can be preserved in soft-bottom sediment of the enclosed coastal area where the rate of DO consumption would often exceed that of DO supply. Hypoxic conditions would further intensify the preservation of organic matter (Middelburg and Levin 2009; Jessen et al. 2017). In the present study, I found sediment TOC was greater in the subsurface layer than in the uppermost layer of the Omura Bay center, and its amount tended to increase toward less oxic conditions (Fig. 9, Table 5). This is consistent with a previous report (Wada et al. 2012) and findings in other marine sediments (Mahmoudi et al. 2015; Jessen et al. 2017). Although the process leading to the accumulation of organic matter in subsurface sediment layer remains an open question, it may be attributable in part to bacterial chemoautotrophy. In fact, Lipsewers et al. (2017) found that the rate of bacterial carbon fixation linked with the oxidation of reduced sulfur compounds was slightly greater in the subsurface layer (5–10 mm) than in the uppermost layer (0–5 mm).
83
I found bacterial abundance within the surface sediment of the center of Omura Bay ranged from 1.6 × 109 to 1.7 × 1011 cells g−1 (Table 5). This is also consistent with a previous study in Omura Bay (Wada et al. 2012). The maximum abundance was observed in summer of 2012 when the bottom-water DO became anoxic (0 µM O2) and a dense bacterial mat containing filamentous sulfur oxidizers developed over the sediment surface during that period. Except for this episodic event, bacterial abundance of the sediment fell within the range of other reports (typically 108 to 109 cells g−1) (Parkes et al. 2000;
Schippers et al. 2012). However, similar to what was found in other marine sediment (Jessen et al. 2017), there was no clear correlation between bottom-water DO and bacterial abundance (see Supplementary Fig. A1). Nevertheless, I found a positive correlation between TOC and bacterial abundance (see Supplementary Fig. A1), which was consistent with other reports that showed a tight coupling between bacterial abundance and TOC (Mahmoudi et al. 2015). This suggests that sediment bacterial abundance in the center of Omura Bay is under bottom-up control.
Effect of bottom-water hypoxia on bacterial community structure and diversity
There are a number of different approaches to study microbial community in the
84
environment. I used a community fingerprinting method, ARISA, which has been widely used to assess bacterial community dynamics (Brown et al. 2005; Hewson, Jacobson-Meyers and Fuhrman 2007; Bertics and Ziebis 2009; Böer et al. 2009; Shade et al. 2011;
Wada et al. 2016). Although ARISA was more sensitive, less time consuming, and relatively low cost compared with other fingerprinting methods (Danovaro et al. 2006;
Cherif et al. 2008; Saro et al. 2014), it has become less frequently used as the advent of next-generation sequencing (NGS) techniques (van Dorst et al. 2014). Nevertheless, a number of recent studies employing both ARISA and NGS techniques to access microbial community demonstrated congruent results of both analyses regarding beta-diversity (Jami et al. 2014; van Dorst et al. 2014; Gobet et al. 2014; Salazar et al. 2016; Tytgat et al. 2016). In addition, alpha-diversity indices, such as Shannon’s index and Simpson evenness calculated from ARISA were highly correlated with those obtained from NGS data (Gobet et al. 2014).
The present results clearly demonstrated that the bacterial community richness, diversity and community structure in the surface sediment would change in response to the availability of bottom DO in Omura Bay. One of the most striking findings of our study is that there was a significant unimodal relationship between DO and sediment bacterial richness (number of OTUs) as well as diversity (H’ and 1-λ) in the surface
85
sediment (Fig. 11). This is consistent with the recent finding that the planktonic bacterial diversity peaked in upper portions of the OMZ, at the interface between well-oxygenated and low oxygen waters, resulting in a unimodal relationship between the richness and DO with peak at 20 μM O2 (Beman and Carolan 2013). Consequently, Beman and Carolan (2013) suggested the edge of the OMZ represents a transition zone where aerobes, anaerobes and microaerophiles may all co-exist. Similarly, the peak in diversity at an intermediate DO (11 μM O2) found in the present study may reflect a transitional condition under which aerobic and anaerobic bacteria can co-exist. As previous studies have suggested, hypoxic conditions overlying the sediment surface may allow physiologically and phylogenetically diverse types of bacteria to co-exist in surface sediment possibly due to the presence of a multitude of electron donor–acceptor couples (Zaikova et al. 2010; Beman and Carolan 2013). In fact, denitrification, anaerobic ammonium oxidation (annamox) and dissimilatory nitrate reduction to ammonium (DNRA) are anaerobic processes that can function under low DO conditions for up to 20–
90 μM (Lorenzen et al. 1998; Gao et al. 2010), 13–20 μM (Kalvelage et al. 2011; Jensen et al. 2011) and 110 μM (Jäntti and Hietanen 2012), respectively. In contrast, abundance and diversity of active aerobic bacteria would simply decrease with decreasing oxygen concentrations.
86
In addition, shifts in the sediment bacterial community associated with DO availability were clearly demonstrated (Fig. 12). The greater similarity in hypoxia than in normoxia regardless of the sediment layer is likely to reflect constraints that select for bacterial population adapting better to lower DO conditions (Fig. 13a, 13b). This is supported by the greater contributions of OTUs that were more abundant under lower DO conditions in the hypoxia community (Table 10).
Comparison between the two sediment layers revealed that the bacterial community similarity tended to decline towards low DO conditions. This may suggest changes in the extent of physical constraint that would homogenize the community across the two sediment layers. In a stable, enclosed coastal environment, macrobenthic as well as meiobenthic fauna are believed to play major roles in mixing sediment particles under normoxic conditions, thereby enhancing solute transport and formation of micro-niches suitable for the aerobic microbial community (Levin 2003). However, most of the bioturbating benthic animals are vulnerable to oxygen deficiency and thus reduce their abundance and irrigating activity during hypoxia. In fact, mass mortality of macrobenthos (polychaetes, bivalves and crustaceans, including mantis shrimp, Squilla oratoria De Haan) in Omura Bay, has been often reported during the summer (e.g., Mori et al. 1973).
Abundance of meiofauna of the bay was also found to decline under hypoxic conditions
87
(Nguyen et al. 2018). It is therefore highly likely that undisturbed and oxygen-depleted sediment would provide a relatively stable anoxic habitat for anaerobic prokaryotes, contributing to the increased diversification between two layers (Fig. 13c).
More than a third of OTUs in each sediment layer overlapped across all oxygen levels (Fig. 10). For example, OTU-643, affiliated with Gammaproteobacteria, (Alteromonadales) dominated across all oxygen levels (Fig. 14). This may be a reflection of resilience of the sediment bacterial population to disturbances caused by changing levels of bottom-water DO. However, each OTU’s relative abundance either positively or negatively responded to bottom-water oxygen concentration (see Supplementary Table A1). The phylogenetic affiliation of the OTUs found in this study was similar to that reported in other coastal sediments (Kawahara et al. 2009; Zinger et al. 2011) with a predominance of Gammaproteobacteria and Deltaproteobacteria (see Supplementary Table A2). It was demonstrated that relative abundance of the OTUs affiliated with Gammaproteobacteria correlated positively with DO, while that of the OTUs affiliated
with Deltaproteobacteria was inversely correlated with DO (see Supplementary Table A2). This suggests that OTUs affiliated with Gammaproteobacteria present in the surface sediment of Omura Bay were mainly composed of a group of bacteria that benefited most under oxic conditions, whereas most Deltaproteobacterial OTUs were composed of those
88
adapted better to hypoxic conditions. These results were supported by Miseq analysis of uppermost sediment samples in 2011, which shows tendency of Gammaproteobacteria to increase toward less oxic condition and opposite trend of Deltaproteobacteria (Fig. 15).
Our data are again consistent with the results of Mahmoudi et al. (2015). Mahmoudi et al.
(2015) reported that surface sediments under normoxia were dominated by Gammaproteobacteria, while those under hypoxia were dominated by Deltaproteobacteria, particularly sulfate-reducing bacteria. It has been widely reported
that Gammaproteobacteria strongly contribute to bacterial population in many marine sediment surfaces (e.g., Urakawa et al. 1999; Inagaki et al. 2003; Feng et al. 2009; Gobet et al. 2011; Orcutt et al. 2011), whereas Deltaproteobacteria often dominate under hypoxic conditions (Devereux et al. 2015; Jessen et al. 2017).
Taken together, it can be generalized that the relative abundance of OTUs affiliated with Gammaproteobacteria and Deltaproteobacteria in surface sediment would be sensitively and inversely related to DO availability of the overlying water.
Oxygen penetration depth (OPD) was not precisely determined in the study site. It is however possible to infer OPD under oxic conditions based on pictures of undisturbed sediment cores that were retrieved from the study site at the time of sampling at the oxic condition (May 2017) and dysoxic condition (August 2017) (see Supplementary Fig. A3).
89
During normoxic period, protruding tubes of polychaetes were abundant on the sediment surface and their burrows were well developed within 0-5 mm depth from the surface.
Under this condition, it is very much likely that OPD reached at least 5 mm depth. If the bioturbation by burrow building fauna is mostly confined within the uppermost sediment layer, and steep gradient of oxygen concentration is established below this depth, it may contribute either to oxidize organic matter in the uppermost layer or to preserve it in the subsurface layer of the sediment. This notion may also help explain the present results of BIOENV test of the relationship between bacterial community structure and environmental factors. Not only DO availability in the overlying water but the burrowing activities of macrofauna may have impact on the bacterial community structure in the uppermost layer. On the other hand, organic matter preserved in the subsurface sediment layer due to poor activities of infauna is likely to account for greater impacts on the bacterial community in the subsurface layer sediment.
90
Fig. 9 Box-and-whisker plots showing the variation in total organic carbon (TOC) obtained for the sediment samples analyzed. (a) and (b) show the TOC content under the same oxygen conditions for uppermost and subsurface layer, respectively. The line inside the box indicates median values. Lines extending from the boxes represent minimum and maximum values. Closed circles represent means.
91
Fig. 10 Percentages of shared OTUs between oxygen conditions from 2011–2013 are displayed as Venn diagrams for uppermost (a) and subsurface layer (b). Percentages of shared OTUs between sediment layers from 2012–2013 are displayed as Box-and-whisker plots (c). Numbers in Venn diagrams are percentages of shared OTUs. SD value is given in parenthesis. The line inside the box indicates median values. Lines extending from the boxes represent minimum and maximum values. Closed circles inside the box represent means.
92
Fig. 11 Scatter plots of diversity indices versus bottom water dissolved oxygen (DO) concentration. (a, d) represent the scatter plots of OTUs versus DO for uppermost (a) and subsurface layer (d). (b, e) represent the scatter plots of Shannon diversity indices versus DO for uppermost (b) and subsurface layer (e). (c, f) represent the scatter plots of OTUs versus DO for uppermost (c) and subsurface layer (f). Solid lines represent the best fit of the quadratic model to the data. Dashed lines represent the best fit linear model to the data. Filled circles represent samples from the uppermost sediment layer and open circles represent samples from the subsurface sediment layer. The DO data were log10 -transformed before analysis. Where infinity was formed by this transformation, the corresponding values were eliminated.