Climate change induces shifts in coastal Baltic Sea surface water microorganism stress and photosynthesis gene expression

The world’s oceans are challenged by climate change linked warming with typically highly populated coastal areas being particularly susceptible to these effects. Many studies of climate change on the marine environment use large, short-term temperature manipulations that neglect factors such as long-term adaptation and seasonal cycles. In this study, a Baltic Sea ‘heated’ bay influenced by thermal discharge since the 1970s from a nuclear reactor (in relation to an unaffected nearby ‘control’ bay) was used to investigate how elevated temperature impacts surface water microbial communities and activities. 16S rRNA gene amplicon based microbial diversity and population structure showed no difference in alpha diversity in surface water microbial communities, while the beta diversity showed a dissimilarity between the bays. Amplicon sequencing variant relative abundances between the bays showed statistically higher values for, e.g., Ilumatobacteraceae and Burkholderiaceae in the heated and control bays, respectively. RNA transcript-derived activities followed a similar pattern in alpha and beta diversity with no effect on Shannon’s H diversity but a significant difference in the beta diversity between the bays. The RNA data further showed more elevated transcript counts assigned to stress related genes in the heated bay that included heat shock protein genes dnaKJ, the co-chaperonin groS, and the nucleotide exchange factor heat shock protein grpE. The RNA data also showed elevated oxidative phosphorylation transcripts in the heated (e.g., atpHG) compared to control (e.g., atpAEFB) bay. Furthermore, genes related to photosynthesis had generally higher transcript numbers in the control bay, such as photosystem I (psaAC) and II genes (psbABCEH). These increased stress gene responses in the heated bay will likely have additional cascading effects on marine carbon cycling and ecosystem services.

Introduction

Climate change related warming is a severe environmental threat affecting all environments (IPCC, 2023) with a recent study suggesting a global warming related temperature rise above pre-industrial levels of 1.7 ± 0.1°C in 2020 (Mcculloch et al., 2024). In addition, the global ocean surface temperature has increased by almost 0.8°C since 1880 (Huang et al., 2017) and between 1971 and 2011 the temperature increased by ~0.11°C per decade (Rhein et al., 2013). Effects on aquatic ecosystems due to an increased water temperature include, e.g., deviations in algal bloom timings, changes in population growth due to increased metabolic rates, invasion of “exotic” species, and increased rate of disease outbreaks (Hoegh-Guldberg and Bruno, 2010). In addition to an increase in temperature, warmer waters contain less soluble oxygen that has resulted in a decrease in global ocean oxygen availability in the last five decades (Schmidtko et al., 2017). This can exacerbate the incidence and extent of oxygen deficient bottom waters, so called ‘dead zones’ (Diaz and Rosenberg, 2008). Environmental changes derived from climate change are therefore important for ecosystem services and functions. However, realistic studies of long-term climate changes on aquatic environments are rare, with many studies being based on short term and/or large changes in temperature, that both lack the influence of annual temperature variations.

Marine phytoplankton are estimated to account for ~45% of all photosynthetic incorporation of CO₂ into organic matter on the planet (Field et al., 1998), of which a large portion is eventually recycled by bacteria in the water column (Moran et al., 2016). An increase in temperature is typically related to a higher rate of photosynthesis in cyanobacteria (Coles and Jones, 2001). However, this pattern is less obvious in for example, marine Picocyanobacteria that are ubiquitous small photosynthetic bacteria for which an increase in temperature and/or CO₂ has shown varying results on photosynthesis (i.e., enhanced or stable) depending on the studied genera (Fu et al., 2007).

To investigate the effect of increased water temperature on aquatic communities, recent studies have used a study site located in the south of Sweden close to the city of Oskarshamn near a nuclear power plant (first in production in 1972). The power plant uses Baltic Sea water for cooling inside the reactor chamber (without coming into contact with the radioactivity) and this thermal water discharge has warmed the surface water in a temperature affected bay for more than 50 years by a yearly average of 5°C (Seidel et al., 2022b). This coastal site was used as a case study for predicted climate change effects such as an increased temperature (Pan et al., 2013). In addition, recent publications at this study site compared to a nearby unaffected control bay reveal differences of microbial communities in surface (Seidel et al., 2022b) and deeper sediment depths (Seidel et al., 2023b). The heated bay warming also leads to compression of sediment geochemical layers such that, e.g., the sulfate/methane transition zone occurs closer to the sediment surface along with shifts in seasonal bottom water microbial communities altering, e.g., the cyanobacterial bloom pattern (Seidel et al., 2022a). Furthermore, a laboratory-based thermal gradient study demonstrated the ability of the control and heated bay benthic community’s response to potential future marine heatwaves (Seidel et al., 2023a). The data show a broader thermal tolerance in the heated bay coupled with elevated RNA transcripts related to stress. However, that study lacks an in situ surface water investigation of the microbial community responses in the water column.

In this study, microorganisms were captured from surface waters from the temperature affected heated bay and the unaffected nearby control bay and their nucleic acids extracted and sequenced. The aim of this study was to investigate differences in the community structure of microbial communities and their RNA transcript-based activities in surface waters of the two bays. The results will give additional information regarding coastal Baltic Sea microbial communities in relation to potential future climate change warming related environmental changes.

Materials and methods

Field sampling

Two naturally fluctuating coastal systems were used that were chosen based upon their proximity to each other. The first was an enclosed bay affected by the warmed discharge water from the Oskarshamn nuclear power plant, Sweden (heated bay; 57.421750 N, 16.666967E) and the second a nearby bay unaffected by the discharge water from the power plant (control bay; 57.433517 N, 16.683700E; Figure 1). The two sites were sampled for surface water on the 16th and 17th April 2018. Temperature, dissolved oxygen, salinity, and pH were measured in the field with sensors (WTW Multiline). From each site, a total of approximately 29 L surface water (~10 cm below surface) was sampled using a 5 L acid-washed container (~1 M HCl, rinsed with deionized water and in situ water) and transferred into acid-washed plus surface water rinsed 20 L and 10 L containers. The water in the 20 L container was used for DNA extraction, while the 10 L container was used for RNA extraction and was mixed with 1 L fix solution to stop microbial activity (5% water-saturated phenol in absolute ethanol) (Feike et al., 2012). The sampled water was transferred to the laboratory (~1.5 h at 11°C–16.5°C) where it was stored overnight in climate rooms set to in situ temperatures (heated bay: 13.6°C, control bay: 10.4°C).

Figure 1

Figure 1. Map of the sampling site. Shown are the heated bay and unconnected control bay with the sampling sites marked with a star.

Cell capture used two types of 142 mm polycarbonate in-line filter holders (Geotech, United States) with different pore-sizes to capture ultra-small bacteria (0.1 μm Omnipore, hydrophilic polytetrafluoroethylene, Merck and 0.2 μm hydrophilic polycarbonate; Isopore, Merck). The filters were sterilized by ethanol and subsequently autoclaved. Cell capture for DNA extraction was via filtering 2 L heated bay water (n = 3 for each filter type) and 1–1.4 L water from control bay (n = 3 for each filter type) with the differing volumes due to eventual clogging of the filters. In the same manner, filters for RNA extraction were collected after filtering 1 L of water mixed with fix solution (n = 3 for each site and filter type; with no clogging). Filters were transferred into sterile 10 mL cryotubes (VWR Collection, VWR) and stored at −80°C until nucleic acids extraction.

Nucleic acids extraction and sequencing

DNA was extracted from the frozen filters with the DNeasy PowerWater Kit (QIAGEN) while RNA was extracted and DNase treated with the RNeasy PowerWater Kit (QIAGEN) according to the manufacturer’s instructions. RNA was synthesized into cDNA and subsequently purified using the Ovation RNA-SeqSystemV2 (NuGEN) and MinElute Reaction Cleanup kit (QIAGEN), respectively. The quality of extracted nucleic acids was measured on a NanoDrop 2000 spectrophotometer while concentrations were measured on a Qubit 2.0 fluorometer.

16S rRNA gene amplicon sequencing was conducted by amplifying the partial 16S rRNA gene from the extracted DNA as described previously (Broman et al., 2017a,b). In brief, primers 341F and 805R (Herlemann et al., 2011) were used with PCR programs according to Hugerth et al. (2014). 16S rRNA gene amplicon library preparation was conducted according to Lindh et al. (2015) and sequenced on a Illumina MiSeq flow cell with a 2 × 300 bp setup. All community DNA and RNA (i.e., cDNA on two HiSeq 2,500 lanes with a 2 × 125 bp setup) samples were sequenced at the Science for Life Laboratory, Stockholm, Sweden. One of the 0.2 μm cDNA filters from the control bay was not successfully sequenced by the facility and therefore, this sample group only has two replicates. A complete list of sample names, sequence facility IDs, number of reads obtained before and after merging paired-ends etc., is available in Supplementary Table S2.

Bioinformatics and statistics

The raw 16S rRNA gene amplicon sequencing data were trimmed, denoised, merged, and chimeras removed by using the nf-core-pipeline for amplicon sequencing called ampliseq (v.2.5.0) within Nextflow (v.23.04.1) using the default settings (Straub et al., 2020). In brief, the forward sequences were trimmed at 269 bp, while the reverse were trimmed at 209 bp and the double removal of the primer sequences were set to true to make sure all primers were removed. The Ampliseq pipeline incorporates DADA2 (v.1.22.0) (Callahan et al., 2016) for denoising and quality control of the sequences followed by generation of amplicon sequence variants (ASVs). Taxonomic classification of the ASVs used SBDI-GTDB (Sativa curated 16S GTDB database, R06-RS202-1; FigShare doi: 10.17044/scilifelab.14869077.v3). Unknown cyanobacterial sequences were blasted against the ‘Phytoref’ database (2017-04-04) to reveal potential eukaryotic diatoms. Further data analysis was conducted using R (v.4.3.1) and RStudio (v.2023.06.2) (Team, 2019). Rarefaction curves (Supplementary Figure S1) were calculated using the ‘vegan’ (v.2.6–4) package using the step size set to 20 and suggested that most microbial diversity was covered in the samples from both bays (Oksanen et al., 2019). Principal component analysis (PCA) based on Bray–Curtis dissimilarities was calculated using relative abundances (normalized counts) and the ‘vegdist’ (distance = bray, default settings) function in vegan plus the ‘prcomp’ function (default settings) in stats (v.4.3.1), while a perMANOVA (999 permutations) comparing the communities between the bays was conducted using the ‘adonis2’ (vegan) function. Alpha diversity analysis (Shannon’s H) was performed using the normalization method of scaling with ranked subsampling within the SRS package using the default settings (v.0.2.3) (Beule and Karlovsky, 2020). A non-parametric Kruskal–Wallis (‘Kruskal.test’ function, default settings) test to check for significant differences between the bays in Shannon’s H diversity used stats (v.4.3.1). Differential abundance analysis of the ASV counts between the two bays was performed using the ALDEx2 package (‘aldex’ function, conditions = bay, default settings, v.1.32.0) (Gloor, 2015). After the analysis, only significant ASVs (based on Benjamin-Hochberg corrected p-values) with an effect size greater than 1/−1 were retained.

The metadenovo (v.2.5.1dev, default settings) nf-core-pipeline for metatranscriptomic data within Nextflow (v.23.04.3) was used for the raw RNA sequences. Briefly, trimming and adapter removal was done by using cutadapt (v.3.4) and trimgalore (v.0.6.7), while sequences were filtered using BBduk (v.39.01). De novo assembly was performed with Megahit (v.1.2.9) while the metadenovo orf-caller option Prokka (which uses prodigal) (Hyatt et al., 2010) within the Prokka software tool (v.1.14.6) was used as open read frame caller (Seemann, 2014). For functional annotation eggnog within Eggnog-mapper (v.2.1.9) (Cantalapiedra et al., 2021) was chosen while EUKulele (v.2.0.3) (Krinos et al., 2021) has been used to annotate taxonomical classification. For prokaryote annotation, GTDB was used within EUkulele while a second annotation with PhyloDB (default) was used to specify eukaryotic transcripts. After this, the eukaryotic transcripts were filtered out. The final mRNA transcript output was analyzed using R within RStudio and expressed as transcripts per million (TPM). Alpha diversity (Shannon’s H) was calculated using TPM values (García-Nieto et al., 2022) and compared to each other using the non-parametric Kruskal–Wallis test. Differential expression analysis was conducted using the DESeq2 (v.1.40.2) package (Love et al., 2014). Before the analysis comparing the two bays, low abundant counts (at least 10 counts) were filtered out that were not detected within the smallest group sample size (here bay; n = 6). The ‘DESeq’ function (settings; test = “LRT,” reduced = ~1, sfType = “poscounts,” minReplicatesForReplace = Inf, fitType = “local”) was used as well as the ‘results’ function (settings; contrast = “heated,” “control”). A PCA was done by log-transforming (‘rlog’ function within the DESeq2 package) the transcripts to account for the compositional nature of relative abundance data (Aitchison, 1986) as well as the ‘prcomp’ function within the stats package. All Plots were generated using the ggplot2 (v.3.5.0) package (Wickham, 2016). A full overview of the used code can be found online, see the data availability statement for further information.

Results

Field data

Differences in surface and bottom water temperatures and salinities between the bays for multiple sampling times over several years are provided in Supplementary Table S1 [data previously published in Seidel et al. (2022a)]. The surface water temperatures fluctuated less in the heated bay and the differences between the bays were more pronounced in winter compared to summer. In addition, the bottom water salinities were similar in both bays over several years while the surface water differences were more pronounced. At the time of sampling for this study, the heated bay sampling site had a water depth of 3.6 m, similar temperature (~13.6°C) and salinity (~6.5 PSU) between bottom and surface waters, plus 12.8 mg/L oxygen and pH 8.40 in the bottom waters. The control bay site was unaffected by the thermal discharge and therefore, had colder surface and bottom water temperatures (10.4°C and 8.1°C, respectively). This site was located at a depth of 2.8 m, it had a lower salinity in the surface than in the bottom water (3.8 PSU and 6.2 PSU, respectively), an oxygen concentration of 15.7 mg/L, and the pH was 7.93 (Table 1).

Table 1

Table 1. Overview of in situ environmental parameters of surface waters in the heated and control bays (both n = 1).

Nucleic acid sequencing data

The Illumina MiSeq 16S rRNA gene sequencing yielded on average 165,783 ± 42,500 reads per sample (one standard deviation shown) with an average of 107,687 ± 44,187 reads after quality control and merging (Supplementary Table S2). Illumina HiSeq 2,500 RNA-seq yielded on average 54 ± 22 million read-pairs per sample. On average 1.09 ± 0.55 million reads per sample were obtained and 14,835 ± 2,472 open read frames were identified (Supplementary Table S2). No statistical differences on 16S rRNA gene taxonomy assignments or RNA transcript counts for annotated gene functions were detected for the different filter pore-sizes (i.e., 0.1 μm and 0.2 μm) in either the 16S rRNA gene amplicon sequencing and metatranscriptome dataset (Supplementary Table S3). Therefore, the samples were treated as replicates.

Differences in community compositions between the two bays

Distinct differences were detected between the microbial communities of the heated and control bay’s surface waters with the Shannon’s H index in both bays varying across the replicates but showing a slightly higher diversity in the control compared to the heated bay that fell short of the traditional statistical significance level (heated 4.35 ± 0.13, control 4.56 ± 0.18, non-parametric Kruskal–Wallis-Test between bays, H₁ = 3.69, p = 0.054; Figure 2). Multivariate PCA analysis using Bray–Curtis dissimilarity on the ASVs showed a statistically significant difference between the two bays (perMANOVA, F_1,9 = 45.4, p = 0.003; Figure 2). A clear separation of the two bays microbial communities was observed on the first axis (95.3% explanation) while a small variation between replicates within the bays was shown on the second axis (2.9%; Figure 2).

Figure 2

Figure 2. 16S rRNA gene amplicon sequencing alpha diversity (mean ± std., n = 6; A), a PCA based on Bray–Curtis dissimilarities (B), and stacked-bar plot on family level with relative abundance per sample >1% (C). Families with a relative abundance <1% make up the remainder of the abundance. Heated bay is shown in orange while control bay is shown in blue.

After filtering phyla with low relative abundances (<0.5% per sample) Actinobacteriota, Bacteroidota, Cyanobacteria, and Proteobacteria were found to be dominant phyla in both bays (Supplementary Figure S2). Comparing the phyla occurrence within each bay showed a higher relative abundance of Proteobacteria in the control bay compared to the heated bay (average relative abundance per bay, heated 17.7 ± 3.6%, control 35.11 ± 3.4%) as well as higher relative abundance of Cyanobacteria (heated 0.7 ± NA %, control 2.5 ± 0.5%). Cyanobacteria was only identified in one replicate of the heated bay samples with a relative abundance >0.5%. In contrast, the Bacteroidota had higher relative abundances in the heated bay compared to the control (heated 31.8 ± 5.1%, control 12.2 ± 2.4%). Differences in microbial communities between the two bays were more distinct on lower taxonomical levels (Figure 2). The top family members with relative abundances >1% included the ecologically diverse Burkholderiaceae (heated bay 4.9 ± 1.2%, control 24.2 ± 3.9%); Microbacteriaceae that are identified in aquatic milieu and in association with, e.g., plants (4.8 ± 3.1, 8.6 ± 3.2); Nanopelagicaeae that can be ultra-small and with streamlined genomes such that they are dependent on a symbiotic partner [1.3 (n = 1), 3.7 ± 0.9]; the predominantly aerobic Flavobacteriaceae found in diverse habitats including marine environments where some strains can harvest light energy (12.7 ± 3.4, 1.8 ± 5.5); and the marine Rhodobacteraceae that are involved in carbon and sulfur cycling (8.3 ± 3.0, 4.5 ± 1.1). In addition, Crocinitomicaceae in the Flavobacteriales order that were originally isolated from polar regions (4.7 ± 1.9) and Cyclobacteriaceae that have a range of properties and habitats including temperature and salinity ranges (6.8 ± 1.5) were only present in the heated bay while the heterocytous cyanobacteria Nostocaceae family (2.6 ± 0.5) was only identified in the in the control bay.

Filtering significant ASV differential abundance (p < 0.05) between the bays with at least an effect size >1/−1 showed clear differences on the family level (Figure 3). For example, ASVs with significantly higher abundance in the heated bay included the Ilumatobacteraceae that were originally isolated from desert soil; the chemoorganotrophic Salibacteraceae that require sea salts for growth; along with the candidate families UBA1997 (Rickettsiales bacterium that are obligate intracellular bacteria), UBA3031 (Burkholderiales), UBA5976 (Nanopelagicales), and UBA955 from the Bacteroidetes that comprise up to 30% of marine coastal bacterioplankton. In contrast, ASVs within for example Burkholderiaceae, the photoheterotrophic Halieaceae family of coastal marine Gammaproteobacteria, the methylotrophic Methylomonadaceae that are also implicated in anaerobic oxidation of methane in interaction with other bacteria, Nostocaceae, the sulfur oxidizing Sulfurmonadaceae, and chemolithoautotrophic sulfur oxidizing Sulfurovaceae were significantly higher in abundant in the control bay. Finally, some families had ASVs that were significantly abundant in either one or the other bay, e.g., the methanol (but not methane) oxidizing Methylophilaceae, Rhodobacteraceae, and the UBA953 that is an Opitutae bacteria found in many environments including seawater and marine sediment.

Figure 3

Figure 3. Differential abundance analysis between the two bays on amplicon sequence variant (ASV) level. The figure shows the sum of significant differential abundant ASVs on the level of family. Orange shades indicate higher effect sizes (standardized mean differences between groups) in the heated bay while the blue shades show higher effect sizes in the control bay.

Differences in RNA transcript levels between bays

Differences in RNA transcripts were observed between the two bays, showing a similar pattern as seen before within the community ASV compositions. Shannon’s H diversity was calculated using TPMs as the base for normalization to investigate potential differences between the diversity of transcripts between the bays. The results showed no difference in the heated bay (5.96 ± 1.1) compared to the control bay (6.7 ± 0.3) that had a slightly lower variation [Kruskal–Wallis-Test, insignificant H(1) = 0, p = 1; Figure 4]. Furthermore, log-transformed transcript counts visualized in a PCA (Figure 4) showed clustering of the samples on the first axis according to bay explaining 66.0% of the variation, while on the second axis a clustering among replicates occurred (explained 6.5% variation).

Figure 4

Figure 4. RNA transcript alpha diversity (Shannon’S H Index, mean ± std., n = 5 to 6) between the two bays (A), a PCA based on log transformed transcripts (B), and the most abundant annotated families per sample with the remainder summarized in ‘Other’ (C). Heated bay is shown in orange while control bay is shown in blue.

Taxonomic annotation of the transcripts resulted in 1124 families with the top five within each sample plotted in Figure 4. While there was a large proportion of genes classified as “other,” distinct family groups were present in each bay. The heated bay showed more variation of different families within the higher abundant transcripts with, e.g., Lentimicrobiaceae that is a member of the phylum Bacteroidetes of which the genus Lentimicrobium is an anaerobic slow-growing bacterium, UBA8524 from the Bacteroidota, and the acidophilic and mostly thermophilic Thermoplasmata. In contrast, the control bay had higher TPM counts within Ketobacteraceae, UBA3072 (WOR-3), Nocardioidaceae, and Methylomirabilaceae. Similar results were also seen within the differential expression analysis, showing significant higher expression for transcripts related to those families mentioned above, such as Lentimicrobiaceae (heated bay) or Ketobacteraceae (control bay; Supplementary Figure S3). Furthermore, most cyanobacterial related RNA transcripts were higher in the control compared to the heated bay (Supplementary Figure S3).

Stress related RNA transcripts

Statistically different RNA transcript numbers showed more transcripts related to stress in the heated compared to the control bay (Figure 5). These transcripts were particularly assigned to the Crocinitomicaceae, Flavobacteriaceae, and Sphingobacteriaceae that are likely adapted to low temperatures and included the heat shock protein HSP70 genes dnaKJ [differences between the levels of expression for each gene between the two groups in logarithmic scale (Log₂FoldChange)] on family level 1.5, 1.3, and 1.1 for Crocinitomicaceae, Flavobacteriaceae, and Sphingobacteriaceae, respectively; the co-chaperonin groS that aids in stress-induced mutagenesis (1.6 and 2.4); and nucleotide exchange factor heat shock protein grpE (2.8 and 3.8). In addition, other statistically different heated bay stress related RNA transcripts included the heat shock-specific RNA polymerase gene groL in family D2472 (3.9) plus the heat shock protein HSP83 (2.1) and HSP90 (2.5) genes in the Myxococcales (unclassified). In contrast, the Burkholderiaceae showed significantly higher stress related gene transcripts in the control bay for atpB (−6.1), groL (−4.8) that is activated by the heat shock-specific RNA polymerase Eσ³², and groS (−7.4); the genome streamlined Nanopelagicaceae identified in freshwater for the chaperonin dnaK (−2.5); and the Microbacteriaceae for dnaK (−4.6).

Figure 5

Figure 5. Differential RNA transcript numbers related to stress according to KEGG Orthology and their annotated potential expressed organisms on family level. Blue shows higher RNA transcript numbers expression in the control bay while orange shows higher in the heated bay. The Log₂FoldChange shows the level of expression of a gene within the tested groups (heated & control) where positive values (orange shades) show higher expression in the heated bay and negative values (blue shades) indicate higher expression in the control bay.

Energy metabolism

The differences in RNA transcript numbers between the two bays were also shown within gene transcripts related to energy metabolism (Figure 6 and Supplementary Figure S3). Oxidative phosphorylation transcripts showed higher numbers in, e.g., atpHG coding for the δ and γ subunits, respectively (Log₂FoldChange for the whole bay 2.2 and 2.3) in the heated bay compared to atpAEFB coding for the α, c, b, and a subunits (−1.8, −3.2, −1.1, and -5.3) in the control bay. This split in ATPase genes between the bays was predominantly along the F₁ component (subunits γ and δ) in the heated bay attributed to the Flavobacteriaceae family compared to the membrane spanning F₀ sub-units (c, b, and a) in the control bay from the, e.g., Cyanobiaceae family. Furthermore, genes related to photosynthesis had generally higher transcript numbers in the control bay, such as the photosystem I genes psaAC (−3.1 and −2.5) along with the photosystem II genes psbABCEH (−5.3, −3.4, −1.9, −3.9, and -4.3) that were for example, attributed to the Elainellaceae, Microcystaceae, and Ketobacteriaceae families. The control bay showed higher transcript numbers of genes related to nitrogen metabolism such as glutamate synthase gltB (−1.2) and glutamine synthetase glnA (−2.1) both potentially involved in ammonia assimilation by the D2472 and Nanopelagicaceae families (Figure 6).

Figure 6

Figure 6. Differential expressed RNA transcript numbers related to energy metabolism according to KEGG Orthology (A) and the associated families (B). The blue color (negative values) shows higher RNA transcript numbers (higher Log₂FoldChange) in the control bay while orange (positive values) shows higher in the heated bay. The x-axis shows the genes transcripts (annotated genes), while the y-axis shows the different metabolism the gene transcripts are higher expressed in (A) and in which taxa they have been higher expressed and associated with (B).

A noticeable result was the distribution of RNA transcripts attributed to methanogens and methanotrophs in the surface waters of the two bays. For instance, the Log₂FoldChange values for methane monooxygenase/ammonia monooxygenase genes pmoC-amoC (−9.1) and pmoB-amoB (−3.8) that based on the transcripts were assigned to the methanotrophic Methylomonadacea (likely pmoBC) and the nitrifying (likely amoBC) Nitrosopumilaceae families in the control bay (Figure 6 and Supplementary Figure S3). These increased mRNA transcripts were reflected in higher relative abundances of average ASVs within the control bay of methanotrophic microorganisms such as Methylobacter (heated 0.04 ± NA control 0.9 ± 0.2% average relative abundance per bay) and Methyloprofundus (heated 0.1 ± NA control 0.1 ± 0.003). The opposite pattern was observed with transcripts annotated to genera of archaeal methanogenic organisms that had higher average TPM values in the heated versus control bay such as Methanobacteriales (20.8 heated bay), Methanobacteriaceae (10.1 heated bay), and Methanocullaceae (17.4 and 4.1).

Discussion

Long-term climate change related temperature increases lead to a global scale decrease in marine diversity (Colby et al., 2020; Chaudhary et al., 2021). This decrease in diversity is detrimental to community health as a high diversity (combination of different species) buffers against environmental changes (Elmqvist et al., 2003). For instance, an increase in temperature influences the microbial community composition and their energy- and nutrient cycling, which further impacts marine carbon cycling (Abirami et al., 2021; Seidel et al., 2023a). A decline in microbial diversity in the heated bay surface waters of this investigation is supported by previous studies of the heated bay benthic microbial communities (Seidel et al., 2022a).

Large differences were also observed between the heated and control bays for the 16S rRNA gene based microbial communities. For example, higher relative abundances of populations belonging to Proteobacteria and Cyanobacteria were detected in the control bay compared to the heated bay. These shifts in relative abundance between the bays was also shown in previous studies where these phyla showed the highest response to temperature changes (Seidel et al., 2023a). In contrast, Paerl and Huisman (2009) show that for example, Cyanobacteria increase in warming waters and one reason for this discrepancy between the studies could be due to the shift in blooming pattern between the bays with an earlier peak in Cyanobacterial blooms within the heated bay (Seidel et al., 2022a), leading to potentially higher relative abundances in the control bay at the time of the sampling (spring bloom April–May). However, this potential shift in bloom time was only based on observations of four time point measurements for bottom water and surface sediment analyses. Another reason for the difference could be the continuous mixing of the surface waters in the heated bay (due to the influx of warm water from the power plant), reducing the relative abundance of cyanobacterial communities in surface waters (Visser et al., 1996). A shift in spring bloom time could also potentially explain the higher relative abundances of Bacteroidetes in the heated bay, as they are known to be major a responder to phytoplankton blooms (Smith et al., 2017). A second difference in microbial community’s between the bays was higher relative abundance of Flavobacteria in the heated bay that supports the increase in bottom water Flavobacteriales (Seidel et al., 2022a). This was potentially due to the increased abundance of families within this phylum playing a role in organic matter degradation (Bissett et al., 2008; Broman et al., 2019) that is likely increased in the warmer waters.

Microbial community differences between the bays were also demonstrated in, e.g., increased RNA transcripts related to stress in the heated bay as compared to elevated energy related photosynthesis gene transcripts in the control bay. The increase in stress related RNA transcripts within the heated bay has also been observed in surface sediments (Seidel et al., 2022b) and benthic waters (Seidel et al., 2022a). Heat-shock proteins are synthesized within the cell when the ambient temperature is just a few degrees centigrade above the cells optimum (Richter et al., 2010). This supports that the microbial community may not have fully adapted to temperatures in the heated bay after over 50 years of acclimation or alternatively, the stable-state protein synthesis of heat-shock proteins is a part of their adaptation to the warmer surface water [reviewed in Tomanek (2010)]. The increased numbers of control bay photosynthesis related gene transcripts may have been related to the higher relative abundance of Cyanobacteria, which were potentially still in spring bloom state at the time of sampling compared to the heated bay.

Furthermore, a study by Mai et al. (2021) showed that higher temperature leads to a reduction of photosynthetic capacity in phytoplankton in tropical oceans. Another study conducted within the Baltic Sea showed that the cyanobacteria Rivularia growth rates decrease with increasing temperatures (Kibria, 2024) while other, e.g., eukaryotic green algae growth increases. This suggests that long-term warmer temperatures can affect Cyanobacteria activity negatively. The data in this study suggested that oxygenic photosynthesis genes had higher transcript numbers in the control bay across different taxa compared to the heated bay while anoxygenic photosynthesis genes were slightly increased in the heated bay. The decreased oxygenic photosynthesis could be a result of the decreased activity of organisms such as Cyanobacteria, as previously described by Bennett et al. (2020) where higher temperatures led to a decrease in oxygenic photosynthesis. The present study could also explain the higher anoxygenic photosynthesis RNA transcripts at warmer temperatures such as resulting from the increased transcripts attributed to Flavobacteria. It has to be mentioned that the mRNA transcript taxonomical annotations were based on the best bitscore, which could result in most likely erroneous assignments such as photosynthesis genes being identified in Candidatus Pelagibacter, which could have been likely a contamination of an assembled genome in the database or overlap due to homologous genes, e.g., CO₂-fixing SAR11 bacteria (Jing et al., 2022). Due to this, these cases of photosynthesis gene mis-annotation were not discussed.

Differences between the two bay’s community structure and RNA transcripts were further underlined by variations in abundances of surface water methanotrophs and methanogens. While there were higher relative abundances as well as transcripts of methanotrophic related microorganisms found in the control bay, higher amounts of methanogenic archaea were found in the surface waters of the heated bay. While the study was carried out in surface waters and both bays were oxic (Supplementary Table S2), these contrasting relative abundances and RNA transcript-based activities might be due to the lower oxygen solubility at raised temperatures favoring anaerobic methanogens (Buan, 2018) versus aerobic methanotrophs (Guerrero-Cruz et al., 2021) in the heated and control bays, respectively. However, it cannot be clarified if these observations were due to potential stratification due to lack of mixing in the control bay or resuspension of sediment due to mixing in the heated bay. Either way, the methanotrophs in the control bay bottom water could lead to methane oxidation while mixing of waters in the heated bay may have resuspended methanogens close to the sediment surface (Seidel et al., 2022b). Alternatively, while methanogenesis most commonly occurs in anoxic conditions, there are records of methane generation in surface waters under oxic conditions in a freshwater lake (Grossart et al., 2011) and in the Baltic Sea (Stawiarski et al., 2019). In general, these findings were in accordance with the discrepancies in microbial community structure between the bays and show that the two bays represented two ecologically contrasting systems.

Two different pore-size filters (0.1 and 0.2 μm) were used without any pre-filtering or a sequential filtration step to capture the surface water microbial communities in the heated and control bays. The lack of significant differences in microbial taxonomical results between the filter pore sizes was likely derived from the filter pores becoming clogged with biomass and particulate matter such that differences between free-living and particle associated microbes could not be discerned. However, the approach removed biases in prokaryotic community structure investigations resulting from sequential filtration of seawater and ensured the whole microbial communities were captured (Padilla et al., 2015).

Finally, it is noteworthy that the variation in alpha diversity within the replicate samples for each bay was high, which could be a result of spatial differences in physical characteristics. Working with natural habitats has many obstacles as even small spatial changes can result in gradients of, e.g., temperature and oxygen concentration creating contrasting environmental niches. Compared to the heated bay, the control bay had a lower salinity, and the control bay salinity was also lower in the surface versus the bottom waters. This indicated an inflow of freshwater in the control bay, as observed in previous studies (Seidel et al., 2022b), and is likely an effect of a freshwater plume from the nearby catchment area to the west. The difference at the time of sampling could be for example due to a predominantly transient or seasonal effect due to snow melt and/or riverine input. Salinity is known to influence bacterial diversity plus community composition (Herlemann et al., 2011) and increases the rate of phytoplankton photosynthesis (Qasim et al., 1972). While the two bays constitute a complex, natural system (that is one of the study’s strengths), the transient shift in salinity was rather small which may also affect microbial populations (Herlemann et al., 2011) but likely not at the observed degree which lead to the conclusion that warmer waters could lead to overall shifts in bacterial community composition over the year and their nutrient cycling. In summary, temperature and its accompanying factors together interplay with location specific environmental conditions to shape microbial communities and result in a complex interlaced relationship that needs to be further investigated to understand future climate change related effects.

Conclusion

The results of the study added information about how a long-term temperature increase could affect microbial communities and their energy-cycling in coastal surface waters that were in accordance with previous results in benthic and surface sediment microbial communities. However, conclusions must be tentative as the surface water is highly variable, and with different seasonal factors playing key roles for microbial community structure that cannot be supported by a single time-point sample. Hence, sampling over a longer time-period would give a better resolution. With this caveat in mind, the difference in microbial composition and higher stress in the heated bay will likely result in a more fragile ecosystem in the near future.

Data availability statement

The original contributions presented in the study are included in the article/Supplementary material, further inquiries can be directed to the corresponding author. The R script detailing the analysis can be downloaded at https://github.com/laseab/surf_comp_bay for more detailed information.

Author contributions

LS: Data curation, Formal analysis, Investigation, Validation, Writing – original draft, Writing – review & editing. EB: Investigation, Validation, Writing – review & editing, Writing – original draft. MS: Investigation, Writing – review & editing, Validation. KB: Investigation, Writing – review & editing, Validation. AF: Conceptualization, Funding acquisition, Validation, Writing – review & editing. SH: Conceptualization, Validation, Writing – review & editing. MK: Conceptualization, Validation, Writing – review & editing. MD: Conceptualization, Funding acquisition, Project administration, Supervision, Writing – original draft, Writing – review & editing.

Funding

The author(s) declare that financial support was received for the research, authorship, and/or publication of this article. MD and MK thank the Swedish Research Council for Sustainable Development, Formas (contracts FR-2020/0008 and 2022–01016); AF thanks the Swedish Research Council, Vetenskaprådet (contract 2020–03519); and EB thanks the King Carl XVI Gustaf 50th Anniversary Fund for Science, Technology and the Environment.

Acknowledgments

The authors wish to thank OKG AB (Oskarshamns Nuclear Energy) for permission to work at the heated bay and surrounding areas. Emelie Nilsson is acknowledged for computational assistance. Bioinformatics analyses were carried out utilizing the Uppsala Multidisciplinary Center for Advanced Computational Science (UPPMAX) at Uppsala University (projects NAISS 2023/22-893 & 2023/6-261). The computations were enabled by resources provided by the National Academic Infrastructure for Supercomputing in Sweden (NAISS) and the Swedish National Infrastructure for Computing (SNIC) at UPPMAX, Uppsala University partially funded by the Swedish Research Council through grant agreements no. 2022-06725 and no. 2018-05973. The authors acknowledge support from the National Genomics Infrastructure in Stockholm funded by Science for Life Laboratory, the Knut and Alice Wallenberg Foundation and the Swedish Research Council, and SNIC/Uppsala Multidisciplinary Center for Advanced Computational Science for assistance with massively parallel sequencing and access to the UPPMAX computational infrastructure.

Conflict of interest

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.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb.2024.1393538/full#supplementary-material

References

Abirami, B., Radhakrishnan, M., Kumaran, S., and Wilson, A. (2021). Impacts of global warming on marine microbial communities. Sci. Tot. Environ. 791:147905. doi: 10.1016/j.scitotenv.2021.147905

Crossref Full Text | Google Scholar

Aitchison, J. (1986). The statistical analysis of compositional data. London: Chapman & Hall.

Google Scholar

Bennett, A. C., Murugapiran, S. K., and Hamilton, T. L. (2020). Temperature impacts community structure and function of phototrophic Chloroflexi and Cyanobacteria in two alkaline hot springs in Yellowstone National Park. Environ. Microbiol. Rep. 12, 503–513. doi: 10.1111/1758-2229.12863

PubMed Abstract | Crossref Full Text | Google Scholar

Beule, L., and Karlovsky, P. (2020). Improved normalization of species count data in ecology by scaling with ranked subsampling (SRS): application to microbial communities. PeerJ 8:e9593. doi: 10.7717/peerj.9593

PubMed Abstract | Crossref Full Text | Google Scholar

Bissett, A., Bowman, J., and Burke, C. M. (2008). Flavobacterial response to organic pollution. Aquat. Microb. Ecol. 51, 31–43. doi: 10.3354/ame01174

Crossref Full Text | Google Scholar

Broman, E., Asmala, E., Carstensen, J., Pinhassi, J., and Dopson, M. (2019). Distinct coastal microbiome populations associated with autochthonous- and allochthonous-like dissolved organic matter. Front. Microbiol. 10:2579. doi: 10.3389/fmicb.2019.02579

PubMed Abstract | Crossref Full Text | Google Scholar

Broman, E., Sachpazidou, V., Pinhassi, J., and Dopson, M. (2017a). Oxygenation of hypoxic coastal Baltic Sea sediments impacts on chemistry, microbial community composition, and metabolism. Front. Microbiol. 8, 1–15. doi: 10.3389/fmicb.2017.02453

Crossref Full Text | Google Scholar

Broman, E., Sjöstedt, J., Pinhassi, J., and Dopson, M. (2017b). Shifts in coastal sediment oxygenation cause pronounced changes in microbial community composition and associated metabolism. Microbiome 5, 1–18. doi: 10.1186/s40168-017-0311-5

Crossref Full Text | Google Scholar

Buan, N. R. (2018). Methanogens: pushing the boundaries of biology. Emerg. Topics Life Sci. 2, 629–646. doi: 10.1042/ETLS20180031

Crossref Full Text | Google Scholar

Callahan, B. J., Mcmurdie, P. J., Rosen, M. J., Han, A. W., Johnson, A. J. A., and Holmes, S. P. (2016). DADA2: high-resolution sample inference from Illumina amplicon data. Nat. Methods 13, 581–583. doi: 10.1038/nmeth.3869

PubMed Abstract | Crossref Full Text | Google Scholar

Cantalapiedra, C. P., Hernández-Plaza, A., Letunic, I., Bork, P., and Huerta-Cepas, J. (2021). eggNOG-mapper v2: functional annotation, orthology assignments, and domain prediction at the metagenomic scale. Mol. Biol. Evol. 38, 5825–5829. doi: 10.1093/molbev/msab293

PubMed Abstract | Crossref Full Text | Google Scholar

Chaudhary, C., Richardson, A. J., Schoeman, D. S., and Costello, M. J. (2021). Global warming is causing a more pronounced dip in marine species richness around the equator. Proc. Natl. Acad. Sci. USA 118:e2015094118. doi: 10.1073/pnas.2015094118

PubMed Abstract | Crossref Full Text | Google Scholar

Colby, G. A., Ruuskanen, M. O., St Pierre, K. A., St Louis, V. L., Poulain, A. J., and Aris-Brosou, S. (2020). Warming climate is reducing the diversity of dominant microbes in the largest high Arctic Lake. Front. Microbiol. 11:561194. doi: 10.3389/fmicb.2020.561194

PubMed Abstract | Crossref Full Text | Google Scholar

Coles, J. F., and Jones, R. C. (2001). Effect of temperature on photosynthesis-light response and growth of four phytoplankton species isolated from a tidal freshwater river. J. Phycol. 36, 7–16. doi: 10.1046/j.1529-8817.2000.98219.x

Crossref Full Text | Google Scholar

Diaz, R. J., and Rosenberg, R. (2008). Spreading dead zones and consequences for marine ecosystems. Science 321, 926–929. doi: 10.1126/science.1156401

PubMed Abstract | Crossref Full Text | Google Scholar

Elmqvist, T., Folke, C., Nyström, M., Peterson, G., Bengtsson, J., Walker, B., et al. (2003). Response diversity, ecosystem change, and resilience. Front. Ecol. Environ. 1, 488–494. doi: 10.1890/1540-9295(2003)001[0488:RDECAR]2.0.CO;2

Crossref Full Text | Google Scholar

Feike, J., Jurgens, K., Hollibaugh, J. T., Kruger, S., Jost, G., and Labrenz, M. (2012). Measuring unbiased metatranscriptomics in suboxic waters of the central Baltic Sea using a new in situ fixation system. ISME J. 6, 461–470. doi: 10.1038/ismej.2011.94

PubMed Abstract | Crossref Full Text | Google Scholar

Field, C. B., Behrenfeld, M. J., Randerson, J. T., and Falkowski, P. (1998). Primary production of the biosphere: integrating terrestrial and oceanic components. Science 281, 237–240. doi: 10.1126/science.281.5374.237

PubMed Abstract | Crossref Full Text | Google Scholar

Fu, F.-X., Warner, M. E., Zhang, Y., Feng, Y., and Hutchins, D. A. (2007). Effects of increased temperature and CO₂ on photosynthesis, growth, and elemental ratios in marine Synechococcus and Prochlorococcus (Cyanobacteria)¹. J. Phycol. 43, 485–496. doi: 10.1111/j.1529-8817.2007.00355.x

Crossref Full Text | Google Scholar

García-Nieto, P. E., Wang, B., and Fraser, H. B. (2022). Transcriptome diversity is a systematic source of variation in RNA-sequencing data. PLoS Comput. Biol. 18:e1009939. doi: 10.1371/journal.pcbi.1009939

PubMed Abstract | Crossref Full Text | Google Scholar

Gloor, G. (2015). ALDEx2: ANOVA-like differential expression tool for compositional data. ALDEX Manual Modular 20, 1–11.

Google Scholar

Grossart, H.-P., Frindte, K., Dziallas, C., Eckert, W., and Tang, K. W. (2011). Microbial methane production in oxygenated water column of an oligotrophic lake. Proc. Natl. Acad. Sci. USA 108, 19657–19661. doi: 10.1073/pnas.1110716108

PubMed Abstract | Crossref Full Text | Google Scholar

Guerrero-Cruz, S., Vaksmaa, A., Horn, M. A., Niemann, H., Pijuan, M., and Ho, A. (2021). Methanotrophs: discoveries, environmental relevance, and a perspective on current and future applications. Front. Microbiol. 12:678057. doi: 10.3389/fmicb.2021.678057

PubMed Abstract | Crossref Full Text | Google Scholar

Herlemann, D. P., Labrenz, M., Jurgens, K., Bertilsson, S., Waniek, J. J., and Andersson, A. F. (2011). Transitions in bacterial communities along the 2000 km salinity gradient of the Baltic Sea. ISME J. 5, 1571–1579. doi: 10.1038/ismej.2011.41

PubMed Abstract | Crossref Full Text | Google Scholar

Hoegh-Guldberg, O., and Bruno, J. F. (2010). The impact of climate change on the world’s marine ecosystems. Science 328, 1523–1528. doi: 10.1126/science.1189930

Crossref Full Text | Google Scholar

Huang, B., Thorne, P. W., Banzon, V. F., Boyer, T., Chepurin, G., Lawrimore, J. H., et al. (2017). Extended reconstructed sea surface temperature, version 5 (ERSSTv5): upgrades, validations, and Intercomparisons. J. Clim. 30, 8179–8205. doi: 10.1175/JCLI-D-16-0836.1

Crossref Full Text | Google Scholar

Hugerth, L. W., Wefer, H. A., Lundin, S., Jakobsson, H. E., Lindberg, M., Rodin, S., et al. (2014). DegePrime, a program for degenerate primer design for broad-taxonomic-range PCR in microbial ecology studies. Appl. Environ. Microbiol. 80, 5116–5123. doi: 10.1128/AEM.01403-14

PubMed Abstract | Crossref Full Text | Google Scholar

Hyatt, D., Chen, G.-L., Locascio, P. F., Land, M. L., Larimer, F. W., and Hauser, L. J. (2010). Prodigal: prokaryotic gene recognition and translation initiation site identification. BMC Bioinformat. 11:119. doi: 10.1186/1471-2105-11-119

PubMed Abstract | Crossref Full Text | Google Scholar

IPCC (2023). AR6 synthesis report: climate change 2023.). IPCC, Geneva.

Google Scholar

Jing, X., Gong, Y., Xu, T., Davison, P. A., Macgregor-Chatwin, C., Hunter, C. N., et al. (2022). Revealing CO₂-fixing SAR11 Bacteria in the ocean by Raman-based single-cell metabolic profiling and genomics. BioDesign Res. 2022:9782712. doi: 10.34133/2022/9782712

PubMed Abstract | Crossref Full Text | Google Scholar

Kibria, G. (2024). Consequences of increased temperature on the abundances of epilithic cyanobacteria and macroalgae in the Baltic Sea. Discover Environ. 2:42. doi: 10.1007/s44274-024-00063-z

Crossref Full Text | Google Scholar

Krinos, A. I., Hu, S. K., Cohen, N. R., and Alexander, H. (2021). EUKulele: taxonomic annotation of the unsung eukaryotic microbes. J. Open Source Software 6:2817. doi: 10.21105/joss.02817

Crossref Full Text | Google Scholar

Lindh, M. V., Figueroa, D., Sjöstedt, J., Baltar, F., Lundin, D., Andersson, A., et al. (2015). Transplant experiments uncover Baltic Sea basin-specific responses in bacterioplankton community composition and metabolic activities. Front. Microbiol. 6:223. doi: 10.3389/fmicb.2015.00223

Crossref Full Text | Google Scholar

Love, M. I., Huber, W., and Anders, S. (2014). Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15:550. doi: 10.1186/s13059-014-0550-8

PubMed Abstract | Crossref Full Text | Google Scholar

Mai, G., Liu, J., Xia, X., Pang, X., Li, B., Yu, L., et al. (2021). Acutely rising temperature reduces photosynthetic capacity of phytoplankton assemblages in tropical oceans: a large-scale investigation. Front. Mar. Sci. 8:710697. doi: 10.3389/fmars.2021.710697

Crossref Full Text | Google Scholar

Mcculloch, M. T., Winter, A., Sherman, C. E., and Trotter, J. A. (2024). 300 years of sclerosponge thermometry shows global warming has exceeded 1.5 °C. Nat. Clim. Chang. 14, 171–177. doi: 10.1038/s41558-023-01919-7

Crossref Full Text | Google Scholar

Moran, M. A., Kujawinski, E. B., Stubbins, A., Fatland, R., Aluwihare, L. I., Buchan, A., et al. (2016). Deciphering Ocean carbon in a changing world. Proc. Natl. Acad. Sci. USA 113, 3143–3151. doi: 10.1073/pnas.1514645113

PubMed Abstract | Crossref Full Text | Google Scholar

Oksanen, J., Guillaume Blanchet, F., Friendly, M., Kindt, R., Legendre, P., Mcglinn, D., et al. (2019). Vegan: community ecology package. Available at: CRAN.R-project.org/package=vegan.

Google Scholar

Padilla, C. C., Ganesh, S., Gantt, S., Huhman, A., Parris, D. J., Sarode, N., et al. (2015). Standard filtration practices may significantly distort planktonic microbial diversity estimates. Front. Microbiol. 6:547. doi: 10.3389/fmicb.2015.00547

PubMed Abstract | Crossref Full Text | Google Scholar

Paerl, H. W., and Huisman, J. (2009). Climate change: a catalyst for global expansion of harmful cyanobacterial blooms. Environ. Microbiol. Rep. 1, 27–37. doi: 10.1111/j.1758-2229.2008.00004.x

PubMed Abstract | Crossref Full Text | Google Scholar

Pan, J., Marcoval, M. A., Bazzini, S. M., Vallina, M. V., and Marco, S. (2013). “Coastal marine biodiversity: challenges and threats” in Marine ecology in a changing world. eds. A. H. M. Arias and C. Menéndez (Boca Raton, FL: CRC Press), 43–67.

Google Scholar

Qasim, S., Bhattathiri, P., and Devassy, V. (1972). The influence of salinity on the rate of photosynthesis and abundance of some tropical phytoplankton. Mar. Biol. 12, 200–206. doi: 10.1007/BF00346767

Crossref Full Text | Google Scholar

Rhein, M., Rintoul, S. R., Aoki, S., Campos, E., Chambers, D., Feely, R. A., et al. (2013). “Observations: ocean” in Climate change 2013: the physical science basis. Contribution of working group I to the fifth assessment report of the intergovernmental panel on climate change. eds. T. F. Stocker, D. Qin, G.-K. Plattner, M. Tignor, S. K. Allen, and J. Boschung, et al. (Cambridge, United Kingdom and New York, NY, USA: Cambridge University Press), 255–315.

Google Scholar

Richter, K., Haslbeck, M., and Buchner, J. (2010). The heat shock response: life on the verge of death. Mol. Cell 40, 253–266. doi: 10.1016/j.molcel.2010.10.006

Crossref Full Text | Google Scholar

Schmidtko, S., Stramma, L., and Visbeck, M. (2017). Decline in global oceanic oxygen content during the past five decades. Nature 542, 335–339. doi: 10.1038/nature21399

PubMed Abstract | Crossref Full Text | Google Scholar

Seemann, T. (2014). Prokka: rapid prokaryotic genome annotation. Bioinformatics 30, 2068–2069. doi: 10.1093/bioinformatics/btu153

PubMed Abstract | Crossref Full Text | Google Scholar

Seidel, L., Broman, E., Nilsson, E., Ståhle, M., Ketzer, M., Pérez-Martínez, C., et al. (2023a). Climate change-related warming reduces thermal sensitivity and modifies metabolic activity of coastal benthic bacterial communities. ISME J. 17, 855–869. doi: 10.1038/s41396-023-01395-z

PubMed Abstract | Crossref Full Text | Google Scholar

Seidel, L., Broman, E., Ståhle, M., Nilsson, E., Turner, S., Hendrycks, W., et al. (2022a). Long-term warming of Baltic Sea coastal waters affects bacterial communities in bottom water and sediments differently. Front. Microbiol. 13:873281. doi: 10.3389/fmicb.2022.873281

PubMed Abstract | Crossref Full Text | Google Scholar

Seidel, L., Ketzer, M., Broman, E., Shahabi-Ghahfarokhi, S., Rahmati-Abkenar, M., Turner, S., et al. (2022b). Weakened resilience of benthic microbial communities in the face of climate change. ISME Commun. 2:21. doi: 10.1038/s43705-022-00104-9

PubMed Abstract | Crossref Full Text | Google Scholar

Seidel, L., Sachpazidou, V., Ketzer, M., Hylander, S., Forsman, A., and Dopson, M. (2023b). Long-term warming modulates diversity, vertical structuring of microbial communities, and sulfate reduction in coastal Baltic Sea sediments. Front. Microbiol. 14:1099445. doi: 10.3389/fmicb.2023.1099445

PubMed Abstract | Crossref Full Text | Google Scholar

Smith, M. W., Herfort, L., Fortunato, C. S., Crump, B. C., and Simon, H. M. (2017). Microbial players and processes involved in phytoplankton bloom utilization in the water column of a fast-flowing, river-dominated estuary. Microbiology 6:e00467. doi: 10.1002/mbo3.467

PubMed Abstract | Crossref Full Text | Google Scholar

Stawiarski, B., Otto, S., Thiel, V., Gräwe, U., Loick-Wilde, N., Wittenborn, A. K., et al. (2019). Controls on zooplankton methane production in the Central Baltic Sea. Biogeosciences 16, 1–16. doi: 10.5194/bg-16-1-2019

Crossref Full Text | Google Scholar

Straub, D., Blackwell, N., Langarica-Fuentes, A., Peltzer, A., Nahnsen, S., and Kleindienst, S. (2020). Interpretations of environmental microbial community studies are biased by the selected 16S rRNA (gene) amplicon sequencing pipeline. Front. Microbiol. 11:550420. doi: 10.3389/fmicb.2020.550420

PubMed Abstract | Crossref Full Text | Google Scholar

Team (2019). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. Available at: https://www.r-project.org/ (Accessed version 3.5.3).

Google Scholar

Tomanek, L. (2010). Variation in the heat shock response and its implication for predicting the effect of global climate change on species’ biogeographical distribution ranges and metabolic costs. J. Experimen. Biol. 213, 971–979. doi: 10.1242/jeb.038034

PubMed Abstract | Crossref Full Text | Google Scholar

Visser, P., Ibelings, B., Van Der Veer, B., Koedood, J., and Mur, R. (1996). Artificial mixing prevents nuisance blooms of the cyanobacterium Microcystis in Lake Nieuwe Meer, the Netherlands. Freshw. Biol. 36, 435–450. doi: 10.1046/j.1365-2427.1996.00093.x

Crossref Full Text | Google Scholar

Wickham, H. (2016). ggplot2: elegant graphics for data analysis. New York: Springer-Verlag.

Google Scholar