Marine ecosystem shifts with deglacial sea-ice loss inferred from ancient DNA shotgun sequencing

Anthropogenic global warming significantly threatens high-latitudinal ocean ecosystems, primarily due to sea-ice loss. This loss causes rising ocean temperatures, increased light transmission, stronger water column stratification, and reduced nutrient availability, ultimately impacting ecosystem services like food supply and the biological carbon pump. The Bering Sea exemplifies this, experiencing a continuous decline in sea-ice duration, projected to worsen in the coming decades. Extreme events, like the record-low sea-ice extent in 2018, highlight the immediate ecosystem consequences of sea-ice loss, including cascading effects through the food web due to reduced energy transfer between trophic levels. While short-term responses are documented, long-term rearrangements in pelagic and benthic communities remain unclear. Understanding past ecosystem changes during the transition from seasonal sea-ice to ice-free conditions is crucial for predicting future changes. Marine sediments provide a valuable archive of climate history, with proxies like biomarkers and microfossils offering insights into past sea-ice coverage. However, these traditional methods have limitations in analyzing non-fossilizing organisms or those with poorly understood fossil records. Sedimentary ancient DNA (sedaDNA) metagenomic shotgun sequencing offers a powerful alternative, enabling broader taxonomic analysis and circumventing PCR biases. This approach, coupled with co-occurrence network analysis (using Spearman networks and Gaussian copula graphical models (GCGMs)), allows for a comprehensive investigation of ecosystem shifts in response to sea-ice loss. The study focuses on sediment core SO201-2-12KL from the western Kamchatka continental slope, leveraging multi-proxy reconstructions (diatoms and IP25 biomarker) to establish a paleoenvironmental context.

Existing literature highlights the vulnerability of high-latitude marine ecosystems to climate change and sea-ice loss. Studies on the Bering Sea demonstrate short-term ecological responses to varying sea-ice conditions, including changes in phytoplankton blooms, copepod abundance, and fish populations (Duffy-Anderson et al., 2017; Stabeno & Bell, 2019; Duffy-Anderson et al., 2019; Siddon et al., 2020). Long-term trends, however, remain poorly understood due to data limitations. Previous paleo-oceanographic research has employed biomarkers (like IP25) and microfossils (diatoms) to reconstruct past sea-ice coverage in the Bering Sea (Belt & Müller, 2013; Matul, 2017; Méheust et al., 2018; Méheust et al., 2016). These studies provide valuable context but lack the taxonomic breadth needed to fully understand ecosystem-level changes. Recent advancements in sedaDNA metabarcoding and metagenomic shotgun sequencing have expanded the potential for reconstructing past marine ecosystems (De Schepper et al., 2019; Zimmermann et al., 2020; Zimmermann et al., 2021; Armbrecht et al., 2021; Armbrecht et al., 2021; More et al., 2019). Co-occurrence network analyses, particularly GCGMs, offer a refined approach to disentangling environmental effects from intrinsic species interactions in metagenomic data (Popovic et al., 2019).

Sediment samples from core SO201-2-12KL (Western Kamchatka continental slope) were collected and processed for sedaDNA extraction using the DNeasy PowerMax Soil Kit. Single-stranded library preparation, optimized for highly degraded ancient DNA, was employed, followed by Illumina sequencing. Bioinformatic processing involved quality control, deduplication (using FastUniq and clumpify), trimming, and merging of overlapping reads (Fastp). Taxonomic classification was performed with Kraken2 using the NCBI nt database, focusing on families of phototrophic bacteria, photo- and heterotrophic protists, marine macrophytes, and Metazoa. Families present in at least three samples with at least ten counts were retained. Resampling (500 times) to a consistent number of counts per sample addressed bias from varying sequencing depths. Negative controls were included to assess contamination. Damage pattern analysis (HOPS) validated the ancient origin of the DNA. Co-occurrence networks were constructed using Spearman rank correlations (igraph) and GCGMs (ecoCopula), incorporating IP25 and SST data as covariates. Statistical analyses included Spearman rank correlation, Benjamini-Hochberg correction, and a Jaccard-index comparison between the Spearman and ecoCopula networks. Interpolated SST and IP25 values were used to assess correlations between families and environmental variables.

Shotgun sequencing yielded 918,186,452 paired-end reads, with ~70.76% passing quality checks. A high proportion (~98.5%) remained unclassified, reflecting the incompleteness of marine species databases. Analysis focused on family-level assignments, revealing 167 families. Pelagic taxonomic assignments were dominated by fish (66.8%), phototrophic protists (10.2%), and phototrophic bacteria (12.8%). Cyanobacteria, Chlorophyta, and diatoms were the most abundant phototrophic families. Zooplankton representation was low. Among metazoans, Salmonidae, Serranidae, and Gadidae were prevalent, along with marine mammals (4.5%). Co-occurrence network analysis (Spearman) identified two modules: a 'sea-ice module' (42 families positively correlated with IP25) and a 'sea-ice-free module' (120 families showing a trend or positive correlation with SST). The ecoCopula network, accounting for environmental effects and mediator taxa, revealed similar modular structure but with a less strict co-occurrence pattern. The sea-ice module included cold-adapted algae (Bacillariaceae, Stephano-discaceae, Thalassiosiraceae) and copepods, linked to Gadidae (Pacific cod, walleye pollock, Polar cod). The ice-free module was characterized by phototrophic bacteria and chlorophytes, with higher abundances of Clupeidae (Pacific herring) and Salmonidae. Benthic sedaDNA analysis revealed a lower read count compared to pelagic, potentially due to lower benthic productivity or database gaps. Bivalves, starfish, priapulid worms, corals, and sea urchins were detected. The pelagic:benthic ratio was higher during glacial periods and lower during warmer periods. Zosteraceae (seagrass) showed increased abundance in the Holocene, positively correlated with SST. Laminariaceae (kelp) was more abundant during glacial periods, co-occurring with Pacific cod.

The findings demonstrate significant ecosystem shifts associated with deglacial sea-ice loss in the Western Bering Sea. The transition from a sea-ice adapted ecosystem dominated by diatoms and cold-water species to an ice-free ecosystem characterized by cyanobacteria, herring, and salmon mirrors predictions for future changes in ice-free regions. The shift from larger to smaller phytoplankton may alter carbon export and benthic food supply, impacting ecosystem services. The northward expansion of salmon and herring, observed in both the sedaDNA data and modern observations, indicates a potential reconfiguration of fishing grounds. The lower abundance of cod and pollock under warmer, ice-free conditions further suggests a significant restructuring of the food web. The study highlights the potential of sedaDNA metagenomics to assess long-term ecosystem responses to climate change, providing insights that are not readily accessible through traditional methods. The results underscore the importance of considering taxonomic shifts in paleoproductivity estimations.

This study provides novel insights into long-term marine ecosystem responses to sea-ice loss using sedaDNA metagenomics. The observed shifts in community composition, particularly the transition from diatom-dominated to cyanobacteria-dominated primary production, have significant implications for carbon cycling and benthic food webs. The northward expansion of commercially important species like salmon and herring is also revealed. Future research should focus on expanding the sedaDNA database, improving temporal resolution, and investigating the functional consequences of these taxonomic shifts in more detail.

The study is limited by the incompleteness of marine organism databases, affecting taxonomic resolution. The temporal resolution, though relatively high for sedaDNA studies, could be improved with increased sampling density. While the use of both Spearman and ecoCopula networks helps to account for indirect interactions, the interpretation of correlations as direct ecological interactions should be cautious, considering the temporal averaging inherent in sedaDNA data. The reliance on interpolated SST data and the potential regional biases in this proxy need to be acknowledged.