{"paper_id":"0a8f8d3e-1a55-47f0-b3d5-dc4a61165162","body_text":"1 \nMulti-millennial genetic resilience of Baltic diatom populations disturbed in 1 \nthe past centuries 2 \nAuthors: Alexandra Schmidt1,2*, Sarah Bolius3, Anna Chagas1, Juliane Romahn4,5,6, Jérôme Kaiser3, Helge 3 \nW. Arz3,  Miklós Bálint4,5,6, Anke Kremp3, Laura S. Epp1* 4 \n1. Environmental Genomics, Department of Biology, University of Konstanz, Constance, Germany 5 \n2. international Max Planck Research School Quantitative Behaviour Ecology & Evolution, 6 \nConstance, Germany 7 \n3. Leibniz Institute for Baltic Sea Research Warnemünde, Rostock, Germany 8 \n4. Senckenberg Biodiversity and Climate Research Centre, Frankfurt am Main, Germany 9 \n5. Loewe Center for Translational Biodiversity Genomics (LOEWE -TBG), Frankfurt am Main, 10 \nGermany 11 \n6. Institute for Insect Biotechnology, Justus Liebig University, Giessen, Germany 12 \n 13 \n* Correspondence 14 \nAlexandra Schmidt, University of Konstanz, Mainaustr. 252, 78464 Constance, Germany 15 \nEmail: Alexandra.schmidt@uni-konstanz.de  16 \nLaura S. Epp, University of Konstanz, Mainaust. 252, 78464 Constance, Germany 17 \nEmail: Laura.epp@uni-konstanz.de  18 \n 19 \nAbstract 20 \nLittle is known about the genetic diversity and stability of natural populations over millennial time 21 \nscales, although the current biodiversity crisis calls for heightened understanding. Marine 22 \nphytoplankton, the primary producers forming the basis of food webs in the oceans, play a pivotal role 23 \nin maintaining marine ecosystems health and serve as indicators of environmental change. This study 24 \nexamines the genetic diversity and shifts in allelic composition in the diatom species Skeletonema 25 \nmarinoi over ~ 8000 years  in the Baltic Sea by analyzing chloroplast and mitochondrial genomes. 26 \nAncient environmental DNA (aeDNA) from sediment cores demonstrates stability and resilience of 27 \ngenetic composition and diversity of this species across millennia in the context of major climate 28 \nevents. Accelerated change in allelic composition is observed from historical periods onwards, 29 \ncoinciding with times of intensifying human activity, like the Roman Empire, the Viking Age, and the 30 \nHanseatic Age, suggesting that anthropogenic stressors have profoundly impacted this species for the 31 \nlast two millennia. The data indicate a very high natural stability and resilience of the genomic 32 \ncomposition of the species and underscore the importance of uncovering genomic disruptions caused 33 \nby human impact on organisms, even those not directly exploited, to better predict and manage future 34 \nbiodiversity. 35 \n 36 \nIntroduction 37 \nHuman activities have had a profound impact on natural environments, resulting in the endangerment 38 \nof a multitude of ecosystems, including marine environments 1. Phytoplankton forms the basis of 39 \nmarine food webs, and, as crucial components of marine ecosystems, changes in this group are central 40 \nto understanding ecosystem shifts2.  Among the phytoplankton, diatoms, which play a significant role 41 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n2 \nin global biogeochemical cycles, are sensitive bioindicators 3, as the different species have specific 42 \nrequirements. Skeletonema marinoi , a prominent marine diatom species 4, is influenced by 43 \ntemperature, migration, and human activities5. For this species, changes in bloom timing6 and optimal 44 \ngrowth temperature7 have already been observed as a result of recent climate change. Investigation 45 \nof marine phytoplankton response to current global changes are mostly based on taxonomic diversity, 46 \nbut recent studies show that diatoms also quickly adapt to environmental changes through acclimation 47 \nand genetic adaptation on a population genomic level 8. In fact, intraspecific genomic variation is an 48 \nimportant factor in the response and resilience of diatoms to environmental perturbations 9. Thus, 49 \ninvestigating the genomic composition of diatom populations and its changes is crucial for 50 \nunderstanding their adaptability to environmental changes. 51 \nThe Baltic Sea, despite its relatively brief history, has undergone significant transformations, making it 52 \nsuitable for studies addressing organism responses to environmental changes. For the last ~ 10,000 53 \nyears, these transformations include a transition between fresh and brackish water around 8,500 -54 \n8,000 years ago 10 and major climate changes since then 11. Due to its shallow enclosed nature and a 55 \nsteep salinity gradient, modern biodiversity is relatively low 12. Current pressures such as pollution, 56 \neutrophication and global warming present significant challenges to the biota13. 57 \nSediments can serve as ecological archives when undisturbed, offering insights into past 58 \nenvironmental changes and biodiversity shifts 14, and provide long time series of phytoplankton 59 \ndynamics. Phytoplankton archives in sediments include biomarkers, microfossils, (living) resting stages, 60 \nand sedimentary ancient DNA (sedaDNA). These remains contain information on biodiversity and can 61 \nreveal changes in adaptive traits 7,15. Here, we develop a time series of population -level responses of 62 \nthe diatom S. marinoi in the Baltic Sea, by employing a targeted approach to enrich chloroplast and 63 \nmitochondrial genomes using hybridization enrichment. This allows us to analyze genetic diversity and 64 \nshifts in allelic composition across complete organellar genomes. We investigate 1) the degree of 65 \nstability, response and resilience  of S. marinoi to natural environmental changes in the Baltic Sea over 66 \nthe last ca. 8000 years, and 2) identify changes in the genetic composition of populations in recent 67 \ncenturies of heightened anthropogenic activities. 68 \n 69 \nResults 70 \nWe analyzed sediment samples of two cores located in the Baltic Sea, the Eastern Gotland Basin (EGB) 71 \nand the Gulf of Finland (GOF) (Fig. 1A), covering a time span up to approximately the Ancylus Lake 72 \n(extrapolated age: 8148 cal yr BP) (BP; with present = 1950 Common Era) (Fig. 1B). We investigated 73 \nthe population structure of S. marinoi over this period at the two locations by focusing on Single 74 \nNucleotide Polymorphisms (SNPs) on the organelle genomes (Fig. 1C). This was achieved by target 75 \nenrichment of the chloroplast and the mitochondrion.   76 \nThis study provides an analysis of the genetic dynamics of S. marinoi populations in the EGB and GOF 77 \nover millennia, revealing insights into the effects of both natural environmental changes and 78 \nanthropogenic activities. The use of RNA baits for hybridization enrichment resulted in a higher degree 79 \nof resolution, enabling the identification of a greater number of SNPs in both mitochondrial and 80 \nchloroplast genomes. The analysis revealed that climate events exert an influence on genetic 81 \nvariations, with specific climate events aligning with temporal variations in the population’s genetic 82 \nmakeup, but that the genetic composition reverted back to a state that was mostly stable across 83 \nmillennia in between these events. Longer lasting patterns of genetic change over time were observed 84 \nat both sites in the past two millennia, coinciding with periods of intensified anthropogenic activities. 85 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n3 \nThese findings imply high resilience of S. marinoi populations across time and climate events across 86 \nmillennia, along with more recent changes. 87 \n 88 \nFig. 1. A) Location of the coring sites and corresponding cores in the central Baltic Sea. Gulf of Finland (GOF; core 89 \nEMB262/12-3GC) and Eastern Gotland Basin (EGB; core EMB262/6 -30GC). B) Sample distribution (in cm) in the sediment 90 \ncores and corresponding ages. C)  Number of SNPs called per 100 base pairs across the genome for each sample. The 91 \nmitochondrion and chloroplast genomes of S. marinoi are shown. The gaps indicate the repeat masked regions. 92 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n4 \nGeneral data assessment 93 \nDNA vs. RNA baits 94 \nIn order to study past S. marinoi population dynamics in the EGB and GOF, we employed DNA and RNA 95 \nbaits for hybridization enrichment of full chloroplast and mitochondrial genomes. The use of RNA baits 96 \nresulted in a higher degree of resolution. Following trimming, DNA and RNA baits yielded 34.9 million 97 \nand 51.8 million sequences, respectively. Of these, 34.8 million (DNA) and 35.3 million (RNA) reads 98 \nwere successfully mapped to reference organelles (Fig. 2A and Supplementary Table S1). In the case 99 \nof RNA baits, 6.41 million reads were mapped to the mitochondrion, while 28.9 million were mapped 100 \nto the chloroplast. For DNA baits, the numbers were 6.27 million and 28.6 million, respectively. The 101 \nRNA baits identified a greater number of SNPs in both the mitochondrial (301 vs. 92) and chloroplast 102 \n(1716 vs. 403) genomes. Consequently, the analysis focuses on the results obtained from RNA baits. 103 \nThe results of the DNA bait analysis can be found in the supplementary section 2.1. 104 \nThe average coverage achieved for the experimental samples was generally high, with an average 105 \nmitochondrial coverage of 78.77 and an average chloroplast coverage of 575.5. This coverage 106 \ndemonstrates the efficiency of the hybridisation capture approach and ensures reliability of the 107 \ngenomic data collected. Furthermore, the Extraction Blanks and Library Blanks produced low to non -108 \nmeaningful coverage (see Supplementary Table S2). All information about the samples, corresponding 109 \nage, location and other metadata is available in Supplementary Table S3. 110 \nAncient DNA – damage pattern 111 \nWe analyzed damage patterns in ancient DNA (aDNA) samples, focusing on C-to-T substitutions. These 112 \nsubstitutions are key aDNA authentication markers due to their pervasiveness in aDNA and their 113 \nfunction in indicating cytosine deamination, a prevalent form of DNA damage. No significant 114 \ncorrelation (R = 0.18, p = 0.315) was found between mapping coverage and substitution frequency, 115 \nsuggesting that damage pattern does not affect the mapping coverage. The degree of damage was 116 \nfound to significantly increase with the age of the samples (R = 0.65, p = 0.0003, Fig. 2B). A greater 117 \nnumber of sequences were mapped to the chloroplast due to its larger genome size. This allows a 118 \ngreater number of reads to be used to detect damage patterns, increasing the overall value of the 119 \ndata. 120 \nGenetic diversity and spatial differentiation 121 \nA comparison of the genetic differences between S. marinoi  retrieved from the EGB and the GOF 122 \nrespectively, revealed that while the two locations exhibited some differences, they also shared a 123 \nconsiderable degree of similarity. This was evidenced by the FST value, a measure used in population 124 \ngenetics to quantify the genetic differentiation between populations, which ranged between 0.05 and 125 \n0.1 (see Supplementary Material, Fig. S6). Notably, the genomic data from EGB shows a higher level of 126 \ngenetic diversity than the one from GOF. This is presented by different measures of distinct variant 127 \ndiversity, with significant differences observed (p-values: 0.002, 0.005, 0.006, see Fig. 2 C-E).  128 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n5 \n 129 \nFig. 2:  Data assessment: A) Sequencing Yield and Single Nucleotide Polymorphisms (SNPs): Post-trimming sequencing yield 130 \nfor DNA and RNA baits, and the number of SNPs detected in mitochondrial and chloroplast genomes for each bait set. B) 131 \nC-to-T Substitution: Analysis of C -to-T substitutions over time for each organelle from all data (GOF & EGB). C -E) 132 \nComparative analysis of normalized nucleotide diversity and diversity indices. C) Relative Nucleotide Diversity (PI), D) 133 \nShannon, and E) Simpson between Eastern Gotland Basin (EGB) and Gulf of Finland  (GOF). Each bar graph shows data for 134 \nboth sites with significant p-values. 135 \n 136 \nEffects of environmental change 137 \nOur study on S. marinoi  organelles at two locations in the Baltic Sea reveals that across several 138 \nmillennia the genetic composition of the population remained stable, punctuated by differences in 139 \nspecific periods, putatively coinciding with environmental changes. A Principal Component Analysis 140 \n(PCA), conducted to investigate differences in allelic composition between sediment layers, showed 141 \nthat the majority of samples are situated within a primary cluster, but a number of smaller clusters are 142 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n6 \nfound that correspond with periods marked by specific environmental events (Fig. 3A). Together, PC1 143 \nand PC2 describe 42.6% of the data variance. 144 \nThe retrieved specific clusters (Fig. 3A) correspond to specific climate periods, including the Holocene 145 \nThermal Maximum (HTM, 10,000-6,000 cal yr BP), Little Ice Age/Late Antique Little Ice Age (LIA, 550 -146 \n100 cal yr BP/LALIA, 1,414 -1,290 cal yr BP), Medieval Climate Anomaly (MCA), and Modern Warm 147 \nPeriod (MWP). We verified that existing variation in allelic composition is not a result of the C -to-T-148 \nsubstitution rate by performing a PERMANOVA, which showed that the C -to-T-substitution rate 149 \naccounted for approximately 6.7% of the variation in the PCA space (R² = 0.067, p = 0.19), though this 150 \nwas not statistically significant (see Supplement, Fig. S5). Additionally, we performed a PERMANOVA 151 \nto assess the impact of Total Organic Carbon (TOC) on the variation in the PCA space, which showed 152 \nthat TOC accounted for approximately 4.5% of the variation (R² = 0.045, p = 0.328), also not statistically 153 \nsignificant. 154 \nIn the EGB, the period under consideration extends to the Ancylus Lake (extrapolated age: 8,148 cal yr 155 \nBP), and responses are observed during the HTM, LALIA, MCA, and LIA periods. The GOF record, which 156 \ncovers the last ca. 4,300 years, displays changes during the LIA and MWP periods. During periods 157 \nwithout profound climate events, we observe a consistency in allelic composition, which points 158 \ntoresilience of the algal populations over time. However, within warm periods (e.g. HTM and MCA), 159 \nthere is a notable degree of variability, particularly for EGB. This pattern of long -term stability and 160 \nresilience, interrupted by distinct variation, is visualized in Fig. 3B, depicting changes along PC1 in both 161 \ncores. This shows temporal variations that align with specific climate events, suggesting that the 162 \npopulation’s genetic makeup is undergoing changes in response to these events. Additionally, 163 \nnucleotide diversity increased significantly during the transition from the Littorina Sea to the Modern 164 \nBaltic Sea, around 4,000 cal yr BP,  and then successively decreased again (Supplementary Material, 165 \nFig S7). 166 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n7 \n 167 \nFig. 3: Principal component analyses of the allelic composition of Skeletonema marinoi  organelles, demonstrating the 168 \ngenetic diversity and population structure. A) PC1 and PC2 are categorized by climate events, total organic carbon (TOC) 169 \nand the age of each sample is given for both sites. B) PC1 versus time for both sites. The red dashed line marks the transition 170 \nfrom the Littorina Sea to the Modern Baltic Sea. Climate events include the Holocene Thermal Maximum (HTM), Late 171 \nAntique Little Ice Age (LALIA), Medieval Climate Anomaly (MCA), Little Ice Age (LIA), Modern Warm Period (MWP) and NA 172 \n= No event 173 \n 174 \nAllele turnover, a measure of the rate at which new alleles replace old ones, provides insights into 175 \ngenetic variation over time. Figure 4 displays patterns of genetic change over time. A Generalized 176 \nAdditive Model (GAM) was used to analyze allele turnover as a function of time. 177 \nAt EGB, an increase in allele turnover is observed around 1,500 -1,000 cal yr BP, indicating a potential 178 \nfor heightened genetic change from the Medieval period onwards. At GOF, turnover remains 179 \nconsistent until around 94 to -34 cal yr BP. While industrialization and anthropogenic impacts may be 180 \nassociated with recent genetic shifts, establishing a causal link between these changes and earlier 181 \nhistorical periods requires caution. 182 \n 183 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n8 \n 184 \nFig. 4: Temporal Analysis of Allele Turnover in Skeletonema marinoi. This figure presents the Generalized Additive Model 185 \n(GAM) fitted to the allele turnover over time (BP). The data includes both Eastern Gotland Basin (EGB) and Gulf of Finland 186 \n(GOF) sites. Key environmental events, such as the Holocene Thermal Maximum (HTM), Late Antique Little Ice Age (LALIA), 187 \nMedieval Climate Anomaly, Little Ice Age (LIA), and Modern Warm Period (MWP), are noted. Additionally, anthropogenic 188 \nperiods and associated shipping activity estimated from literature sources are shown. 189 \n 190 \nDiscussion 191 \nThe millennial-long stability, spanning from the Ancylus Lake (approx. 8,148 calibrated years BP), to 192 \nLittorina Sea and Modern Baltic Sea to 1,500 –1,000 cal yr BP in the Eastern Gotland Basin (EGB), and 193 \nfrom around 4300 to 94 –34 cal yr BP in the Gulf of Finland (GOF), aligns with findings from other 194 \nphytoplankton studies, albeit with a different temporal perspective. The findings of these studies 195 \nindicate that long -term genetic stability is frequently maintained in marine populations due to their 196 \nsubstantial population size and high dispersal capacity, even in the face of significant environmental 197 \nchanges 16. These factors contribute to a buffering effect against genetic drift and local extinctions, 198 \nallowing populations to maintain genetic diversity over long periods of time. The ability of S. marinoi 199 \nto cope with changing environmental conditions through reversible shifts in genetic composition as 200 \nshown here for the Baltic Sea further supports this stability. This resilience is a common trait among 201 \nphytoplankton, allowing them to respond to environmental changes without losing overall genetic 202 \ndiversity17. The observed stability is expected given the ecological and evolutionary characteristics of 203 \nphytoplankton populations. Large population sizes, high reproductive rates, and the potential for gene 204 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n9 \nflow between populations contribute to the maintenance of genetic diversity and stability over long 205 \nperiods of time 18. Other findings, such as studies of diatoms and other phytoplankton species have 206 \nreported patterns of genetic stability and adaptability of phytoplankton to environmental change 19. 207 \nThese findings imply that the observed stability in S. marinoi populations is not only expected, but also 208 \nindicative of the broader ecological and evolutionary dynamics of marine phytoplankton. Our data 209 \ndemonstrates that this stability can be upheld across several millennia. 210 \nHowever, this stability is not without punctuated interruptions. Our results show that climatic events 211 \nhave an influence on the allelic composition of S. marinoi  populations. For example, shifts in the 212 \ngenetic composition of the EGB population were observed during the Holocene Thermal Maximum, 213 \nthe Late Antique Little Ice Age, and the Medieval Climate Anomaly. Similarly, the GOF population 214 \nshowed changes during the Little Ice Age and the Modern Warm Period. These shifts highlight the 215 \nability of S. marinoi to cope with different climatic events. The discrepancy between EGB and GOF may 216 \nbe due to the slightly higher genetic diversity in the EGB (Fig. 2C -E), which may facilitate adaptation. 217 \nThe differences between the two sites can also be attributed to their distinct geographical 218 \ncharacteristics. The EGB, located in the Baltic Proper, may experience an influx of migrants from other 219 \npopulations, increasing genetic diversity and resilience. In contrast, the genetic composition of the GOF 220 \npopulation could be influenced by its unique environment, characterized by freshwater influx and 221 \nlower salinity levels20, as well as more extensive ice cover in the northern regions of the Baltic Sea21. In 222 \nparticular, samples from the Little Ice Age and the Modern Warm Period show more pronounced 223 \nresponses in the GOF. The allelic composition remains relatively constant in the absence of climatic 224 \nevents, indicating a stable genetic structure under stable conditions and once a shift has manifested. 225 \nAfter this long period of punctuated genetic stability, recent centuries reveal a novel pattern of rapid 226 \nallelic turnover, and this coincides with increased human activity. This rapid turnover is particularly 227 \nevident in the GOF during periods such as the pre-Roman Iron Age, the Viking Age, the Hanseatic Age, 228 \nand the Industrial Revolution. These results highlight the possible impact of human activities on the 229 \ngenetic dynamics of S. marinoi populations. This increased activity is manifested by increased shipping 230 \nactivity in the Baltic Sea 22,23. Ballast water exchange 24 in the last two centuries could have been a 231 \nspecific agent of increased changes in haplotype composition, causing rapid translocation of 232 \npopulations.  This rapid allelic turnover contrasts with the more gradual and reversible shifts observed 233 \nin response to natural climatic events, highlighting the influence of anthropogenic factors on the 234 \ngenetic composition of these populations. Due to its coastal position, the GOF is more exposed to these 235 \nchanges, and perhaps the phytoplankton is also more directly affected. 236 \nThe rapid allelic turnover observed in recent centuries, particularly in response to human activities, 237 \nunderscores the influence that anthropogenic factors may have on the genetic dynamics of 238 \npopulations. Our study on the population dynamics of S. marinoi populations provides valuable insights 239 \ninto the long-term resilience and adaptability of marine phytoplankton. The observed genetic stability 240 \nover several millennia, punctuated by reversible shifts in response to climatic events, underscores the 241 \ninherent resilience of these populations. However, the recent rapid allelic turnover highlights the 242 \npotential impact of human activities on population dynamics and structure, and thus presents an 243 \nopportunity for further investigation.  244 \nAt the same time it is important to acknowledge the limitations of our study. The limited number of 245 \nsamples included may affect the overall strength and generalizability of our findings. In addition, the 246 \nlack of precise tests limits our ability to draw definitive conclusions about long-term genetic trends and 247 \nthe full extent of human impact. To date, we also lack a sufficient understanding of the distribution 248 \nand taphonomy of ancient eDNA, but in sediments, DNA of small organisms, such as phytoplankton, 249 \nhas been shown to give a representative signal of a water body 25. The data of the two cores display a 250 \nhigh level of congruency, and despite the limitations, this indicates a link between environmental 251 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n10 \nchange, human activity, and genetic dynamics. Future studies with larger sample sizes and e.g. more 252 \nlocations are needed to further explore the complex interactions between environmental change and 253 \ngenetic responses. 254 \nThe results of our study provide a solid base for future research on the effect of current global change 255 \non the genetic composition of natural populations across millennia. By understanding the genetic 256 \nstability and adaptability of  S. marinoi  populations, we can better understand ecological and 257 \nevolutionary dynamics of marine phytoplankton. This knowledge is crucial for predicting how these 258 \npopulations may respond to current and future environmental changes, including climate change and 259 \ndirect anthropogenic pressures, and thus offer a valuable foundation for future conservation efforts. 260 \nOur study provides novel insights into the genetic resilience and adaptability of S. marinoi populations. 261 \nIt demonstrates that these populations maintain stability over millennia and respond dynamically to 262 \nboth natural climate fluctuations and human-induced changes. In sum, it highlights the importance of 263 \nconsidering both natural and anthropogenic drivers and their relative impacts when assessing and 264 \npredicting genetic resilience of marine ecosystems. 265 \n 266 \n 267 \nMaterial & Methods 268 \nStudy area 269 \nThe Baltic Sea is a relatively young and shallow brackish water system currently facing a significant risk 270 \nof increasing hypoxia, a condition characterized by low oxygen levels and bottom water anoxia. This 271 \nphenomenon is primarily attributable to the reduction in dissolved oxygen and nutrient discharge 272 \nresulting from warmer water temperatures, which in turn promote the formation of algal blooms 26. 273 \nThese blooms decompose and consume oxygen, leading to the development of hypoxia. In recent 274 \nhistory, the Baltic Sea has experienced a notable increase in eutrophication, a process driven by the 275 \nexcessive input of nutrients such as nitrogen and phosphorus from agricultural runoff, wastewater 276 \ndischarge, and industrial activities 13. The enrichment of nutrients has resulted in the proliferation of 277 \nphytoplankton and algal blooms, which, upon decomposition, have further depleted oxygen levels in 278 \nthe water. The Baltic's unique position as a land-enclosed sea with high human and industrial activity, 279 \ncoupled with its low biodiversity, makes it particularly vulnerable to these changes. 280 \nThe  Baltic Sea's history is well known. It has undergone substantial environmental changes in the past. 281 \nAfter the end of the last glaciation approximately 10,000 years ago, it went through several fresh- and 282 \nsaltwater stages, starting as a large meltwater lake, then transitioning into various states of salinity 283 \ndue to glacial retreat and land rebound27. These stages included the Baltic Ice Lake, the Yoldia Sea, the 284 \nAncylus Lake, and the Littorina Sea 28. The Baltic Sea’s temperature and salinity have fluctuated over 285 \ntime, with notable increases during the Holocene Thermal Maximum and the Medieval Climate 286 \nAnomaly, and a recent increase since 185027,29,30. 287 \nWithin the Baltic Sea this study covers a time series of two distinct locations: The Eastern Gotland 288 \nBasin, located in the Baltic Proper with a profound water depth of max. 249 meters, and the 289 \ncomparatively shallower Gulf of Finland, reaching a depth of 81 meters. At both locations, we collected 290 \nsediment cores to inform on the effects of environmental changes on the population genetic 291 \ncomposition of S. marinoi. 292 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n11 \nSampling  293 \nTwo sediment cores were taken (Fig. 1A) in April 2021, during expedition EMB262. 1) Eastern Gotland 294 \nBasin (EGB, 57°17.004'N, 020°07.244'E, 241 m water depth, core EMB262/6-30GC); 2) Gulf of Finland 295 \n(GOF, 59°34.443’N, 023°36.461’E, 81 m water depth, core EMB262/12-3GC). The cores (Fig. 1B, ca. 500 296 \ncm) were taken using a gravity corer (GC). The sampling for the EGB began at a depth of 32 cm. The 297 \ncores were sampled on board of the research vessel using sterile syringes according to Epp et al.31 and 298 \nimmediately frozen for storage. 299 \nCore dating was performed as described in Schmidt et al. 32. This was done by correlating organic 300 \ncarbon records from post -Littorina transgression sediments in the Baltic Sea. The chronology of the 301 \ncores is determined by the relative S content and the Br/K ratio. The Br/K ratio reflects changes in the 302 \nbulk organic carbon content of the sediments. The data are visually matched with XRF and organic 303 \ncarbon data from dated sediment cores from the same or nearby locations for different historical 304 \nperiods. For the sediment core from the Gulf of Finland, an independent Bayesian age model was 305 \nconstructed based on radiocarbon dates of bulk organic matter in the sediments. This allows the 306 \nsediments to be accurately dated to around 2300 BC. The same core samples were used as in Schmidt 307 \net al.32 and in this study. 308 \n 309 \nLaboratory work 310 \nDNA extraction 311 \nThe PowerSoil Pro Kit from Qiagen (Hilden, Germany) was used to extract DNA from both sediment 312 \ncores (n=26, see supplementary Table SX) covering various depths (ages). For the extraction 0.5 g of 313 \nsediment was used per sample, resulting in 26 samples and 2 extraction blanks. The extraction process 314 \nwas carried out according to the manufacturer's protocol with some modifications including an 315 \novernight incubation at 56 °C with the addition of 20 μL proteinase K (20 mg/ml) to enhance sample 316 \nlysis. The washing steps were performed using a Qiagen Vacuum Pump (Hilden, Germany). 317 \nSubsequently, samples were centrifuged for 3 minutes. After centrifugation, the elution process was 318 \nperformed in two steps. Each part involved adding 75 μL of elution buffer to the sample, allowing it to 319 \nincubate for 5 minutes, and then collecting the eluate. Both eluates were collected in the same tube. 320 \nLibrary preparation 321 \nA single -stranded library approach according to Gansauge et al. 33 was performed to maximize the 322 \nretrieval of genetic information from these degraded sedaDNA samples. This method, optimized for 323 \nfragmented DNA, has been shown to significantly increase the sequence yield and preserves unique 324 \nmolecules, thereby providing a more comprehensive genomic analysis. The library preparation was 325 \nconducted with 40 ng input DNA.  All 26 samples and 2 extraction blanks were processed in the library 326 \npreparation. Additionally, 2 library blanks were created.  327 \nAfter the indexing PCR, each library was cleaned using the MinElute PCR purification Kit (50,250; 328 \nQiagen) resulting in a final amount of 20 μL per sample (For more detailed info, see Supplement Section 329 \n1.1.). The samples were then analyzed on a 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, 330 \nUSA) using the Agilent Bioanalyzer High Sensitivity DNA Analysis Kit (Agilent Technologies, Santa Clara, 331 \nCA, USA). Due to adapters still being present in the samples, we performed an additional purification 332 \nusing the HighPrep TM bead cleanup (MagBio Genomics Inc., Gaithersburg, MD, USA). 1.6x of 333 \nHighPrepTM PCR reagent was used. The DNA concentration was measured using ds -DNA HS Assay Kit 334 \nand the Qubit® 4 fluorometer (Invitrogen Thermo Fisher, Waltham, MA, USA). The libraries were 335 \ncombined into pools of four, each with the same DNA amount. The extraction and library blanks were 336 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n12 \nincorporated in the same volume as the sample with the lowest concentration. The libraries were then 337 \nfurther used for hybridization enrichment of chloroplast and mitochondrial genomes.  338 \nHybridization enrichment 339 \nIn this study, a hybridization enrichment approach was employed. Hybridization enrichment allows for 340 \ntargeted isolation of specific sequences, improving sequencing efficiency by reducing the presence of 341 \nnon-target DNA. Thus, providing a better representation of the studied species, enabling a more 342 \ndetailed and precise analysis of its genetic composition. 343 \nBait design 344 \nThe design of the bait set was carried out for two target reference organelle genomes: the 345 \nmitochondrion and the chloroplast (GenBank PRJNA493755). Together, these genomes consist of 346 \n170,788 nucleotides (nt), with 127,202 nt from the Chloroplast and 43,586 nt from the mitochondrion, 347 \nand an average GC content of 30.9%. To improve the accuracy of the baits, the contigs were 348 \nsoftmasked for simple and low-complexity repeats. 349 \nSubsequently, baits were designed with a length of 80 nt and a tiling coverage of 4x. This means that, 350 \non average, each nucleotide in the target sequence is covered by four different baits, thereby ensuring 351 \nredundancy and enhancing the probability of successful hybridization. This process led to the 352 \ngeneration of 10,000 probes. The baits were then filtered based on softmasking for simple and low 353 \ncomplexity repeats. Only those with 35% or less softmasking retained. A total of 7,985 baits passed 354 \nthis filtering step. 355 \nEnrichment 356 \nThis bait design was then used for two different bait sets: 1) RNA baits (BioCat, Heidelberg, Germany); 357 \n2) DNA baits (Integrated DNA Technologies, IDT). The libraries were enriched with both RNA and DNA 358 \nbaits, resulting in two sets of libraries: 1) RNA baits enriched; 2) DNA baits enriched. General steps are 359 \nexplained in the following. More detailed information on the used protocols can be found in the 360 \nsupplementary data (section 1.2.).  361 \nThe enrichment with RNA baits was conducted according to the myBaits v.5.02 manual. A two -round 362 \nenrichment strategy was performed. The hybridization mix was prepared and transferred to a rotation 363 \noven once the 24-hour incubation at a hybridization temperature of 63°C started. In the following bead 364 \ncleanup, the libraries were resuspended and amplified. Each 27.4 μL PCR reaction included the 365 \nfollowing nine components: 1) 3.45 μL of DEPC treated H 2O; 2) 2,5 μL 10X HiFi PCR Buffer (Invitrogen 366 \nThermo Fisher, Waltham, MA, USA); 3) 0.25 dNTPs (25mM); 4) 1 μL BSA (20 mg/ml); 5) 1 μL MgSO4 (50 367 \nmM); 6) 0.2 μL Platinum ™ Taq DNA -Polymerase High Fidelity (5 U/ μL) (Invitrogen Thermo Fisher, 368 \nWaltham, MA, USA); 7) 2 μL IS5_bridge_P5 (10 μM); 8) 2 μL IS5_bridge_P5 (10 µM), 9) 15 μL library 369 \npool. The amplification was run using the following settings: 1 minute initiation at 94°C, 14 cycles of 370 \n15 seconds denaturation at 94°C, 20 seconds annealing at 60°C and 1 minute extension at 68°C, 371 \nfollowed by 2 minutes of final elongation at 68°C and afterward the sample was stored at a 372 \ntemperature of 20°C. Afterwards, the samples were purified with beads and eluted in 10 μL. The 373 \nsecond round of enrichment was performed in a similar way. The PCR was conducted in 8 cycles and 374 \nthe pools were eluted in 17 μL after the bead purification (see more detailed in supplement section 375 \n1.2.2.). 376 \nEnrichment with DNA baits was also performed in two rounds following the IDT (Integrated DNA 377 \nTechnologies, Löwen, Belgium) xGen TM hybridization capture of DNA libraries protocol (option: 378 \nAMPure XP Bead DNA concentration protocol), with the following modification. The hybridization 379 \ntemperature was set to 63°C. Once this temperature was reached on the thermal cycler, the tubes 380 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n13 \nwere moved to a rotation oven and incubated for 16 hours (see more detailed in supplement section 381 \n1.2.1.). 382 \nAll resulting pools were quantified using Qubit® ds -DNA HS Assay Kit and Agilent Bioanalyzer HS Kit. 383 \nAfterwards, the pools enriched with RNA baits were collectively combined into a single pool, while 384 \nthose enriched with DNA baits were combined into another separate pool, each in equal DNA amounts. 385 \nThis resulted in two distinct library pools. Both of these pools were sent for Illumina sequencing 386 \n(paired-end, 2 x 150 bp, MiSeq-v3) at Fasteris SA (Geneva, Switzerland). 387 \n 388 \nBioinformatics 389 \nThe sequences were trimmed, mapped against references, and converted to finally calculate mapping 390 \ncoverage and for variant calling. First, sequence file names were modified to fit the sample identifier 391 \nusing python v3.10.13. Then raw sequence reads were processed using Autotrim v0.6.1, a tool that 392 \nautomates the quality control and trimming of high -throughput sequence data 34. Autotrim employs 393 \nTrimmomatic v0.3935, to trim sequences based on quality and remove adapters, FastQC v0.12.1 36 to 394 \nassess the quality of the raw sequence data, and MultiQC v1.14 37 to aggregate the results into a single 395 \nreport. Following the initial processing, the trimmed reads were mapped to the reference chloroplast 396 \nand mitochondrial genomes (GenBank: PRJNA493755) using BWA -MEM v0.7.1738. BWA is a software 397 \npackage for mapping low-divergent sequences against a reference genome, which is particularly useful 398 \nfor precisely aligning sequencing reads.  399 \nThe mapped reads were then converted to bam format using Samtools v1.1539, a suite of programs for 400 \ninteracting with high-throughput sequencing data. In addition, the per-sample mapping coverage was 401 \ncalculated. After careful consideration, we decided against deduplicating our sequence data as 402 \ndeduplication resulted in a significant loss of data. This is probably due to the nature of our short and 403 \nfragmented sequences. In addition, each variation, even those present in duplicates, may hold 404 \nbiologically significant information. Therefore, removing duplicates may have eliminated critical 405 \ngenetic diversity, compromising the integrity of our analysis.  406 \nVCFtools v0.1.1640 and BCFtools v1.1939 were used for calling and filtering variants from the alignment. 407 \nVariant calling is the process of identifying differences, such as single-nucleotide polymorphisms (SNPs) 408 \nand insertions/deletions (indels), between the sequenced sample and the reference genome. Variants 409 \nwere filtered to include only those with a quality score (Q) of 30 or higher. The identified variants were 410 \nthen phased into haplotypes using Beagle v5.441. Haplotyping, or phasing, is the process of determining 411 \nthe distribution of alleles of multiple variants along each chromosome. In this context, it can be said it 412 \naccounts for unique combinations of alleles at variant sites (distinct variations) across the chloroplast 413 \nand mitochondrial genome. 414 \nThe pattern of ancient DNA damage was assessed using mapDamage2 v2.2.2 42. This tool quantifies 415 \npatterns of DNA damage in ancient samples, which can provide insights into the preservation and 416 \nauthenticity of ancient DNA. 417 \n 418 \nAnalyses 419 \nThe bioinformatic output was converted into a table format and for further analyses with R v4.3.1: The 420 \nmapping coverage (i.e., the number of reads mapped against reference genome), allele frequencies, 421 \ndistinct variants, nucleotide diversity and C -to-T-substitution rate were combined with the sample 422 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n14 \nmetadata (Supplement, Tab. S3) using the R package “tidyverse” 43. The data was visualized using the 423 \npackages “tidyverse´”, “ggpubr” v0.6.044 and “tidypaleo” v0.1.345.  424 \nIn order to facilitate a comparison of the genetic diversity between GOF and EGB,  nucleotide diversity 425 \n(resulting from Samtools) and distinct variants (resulting from Beagle) were normalized based on the 426 \nmapping coverage for mitochondria and plastid. Second, Pearson's correlation was employed to verify 427 \nthe consistency of nucleotide diversity between the mitochondrion and chloroplast. For comparison, 428 \nthe dataset was subset to the overlapping time period of both locations. Subsequently, a t -test was 429 \napplied to compare and test the significance of the nucleotide diversity (mitochondrial and chloroplast 430 \ndata combined) between the two sites. In addition, each distinct genetic variant was assigned a unique 431 \nidentifier for easy reference. These identifiers were utilized to calculate two diversity indices: the 432 \nShannon index and the Simpson index. These indices were compared between the EGB and GOF sites 433 \nusing a t-test,  434 \nTo further identify spatial differences between GOF and EGB, we analyzed allele frequencies (resulting 435 \nfrom  VCFtools & BCFtools) and calculated the Fixation Index (FST) using corresponding samples from 436 \nboth locations. The FST is a measure of population differentiation attributable to genetic structure. 437 \nThe FST value ranges from 0 to 1, with 0 indicating complete interbreeding (i.e., the two populations 438 \nare freely intermixing) and 1 indicating that all genetic variation is explained by the population 439 \nstructure (i.e., the two populations do not share any genetic diversity). We implemented a custom R 440 \nfunction to calculate the FST over time. This process involved sorting the data by age, subsetting 441 \nsamples within a specified time window, and computing FST based on allele frequency variance. The 442 \nallelic composition, derived from parsing VCFtools output, was analyzed using Principal Component 443 \nAnalysis (PCA) to examine the genetic structure.. We also performed Pearson’s correlation and 444 \nPERMANOVA to assess the impact of C-to-T substitution rates and Total Organic Carbon (TOC) on the 445 \nvariation in the PCA space, specifically the first two principal components (PCs). These analyses helped 446 \ndetermine the extent to which these factors influenced genetic composition. For PCA and 447 \nPERMANOVA we also used the functions provided by the package “vegan” v2.6.446. 448 \nThe C -to-T-substitutions at the first position from the mapDamage2 output were examined to 449 \nauthenticate the ancient DNA (aDNA) and analyze the damage patterns. A Pearson's correlation test 450 \nwas conducted to assess the influence of mapping coverage on these substitutions. This was done to 451 \nevaluate if coverage has an impact on the observed damage. The Baltic Sea Mn/Ti ratio data were also 452 \nintegrated into the analysis. 453 \nTo assess the genetic diversity and adaptability of the population, we calculated the rate of allele 454 \nturnover. This measure, derived from the ratio of significant allele changes (greater than 1%) to total 455 \ngenerations, reflects the estimated rate of genetic change over time based on the present data. Given 456 \nthe low FST values, we inferred that the two sites, EGB and GOF, belong to a single population. The 457 \ncalculation covered all of the data, from the oldest sample from EGB to the youngest sample from GOF, 458 \nproviding insight into the evolving genetic makeup of the population. 459 \nTo identify phases of change in the allele turnover across the entire dataset, a generalized additive 460 \nmodel (GAM) was fitted using the “gam” function from the “mgcv” package47. This model was applied 461 \nto the complete dataset, including periods where we did not have direct samples, as the allele turnover 462 \nwas estimated based on the available data. We calculated the mean allele turnover for the present 463 \ndata points, which were then included in the GAM plot. Allele turnover was used as a response variable 464 \nto fit a gam model as a function of “age”. 465 \nCorresponding metadata to this analysis can be found in supplementary material (Table S3). Data on 466 \nshipping activity in the Baltic Sea were obtained from three main sources: Kontny 48, Mägi 49, and 467 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n15 \nRytkönen et al. 50. Kontny's work provided insights into maritime contacts during the Roman and 468 \nMigration Periods (1st -7th centuries AD), while Mägi's research focused on the role of the Eastern 469 \nBaltic in Viking Age communication across the Baltic Sea. Rytkönen et al.'s research on the statistical 470 \nanalysis of Baltic Sea shipping revealed a significant increase in maritime traffic during the late 20th 471 \nand early 21st centuries. These sources were used to establish a timeline of shipping intensity from 472 \n1200 BP to the present. 473 \n 474 \nAcknowledgements 475 \nThis study was funded by the K314/2020 grant of the Collaborative Excellence Programme of the 476 \nLeibniz Association. The authors acknowledge support by the High Performance and Cloud Computing 477 \nGroup at the Zentrum für Datenverarbeitung of the University of Tübingen, the state of Baden -478 \nWürttemberg through bwHPC and the German Research Foundation (DFG) through grant no INST 479 \n37/935-1 FUGG. AS acknowledges funding from the International Max Planck Research School for 480 \nQuantitative Behaviour, Ecology and Evolution. MB is supported by the Hessian Ministry of Higher 481 \nEducation, Research and the Arts through the LOEWE Centre for Translational Biodiversity Genomics.  482 \nFor help with preparing reagents for library preparation and enrichment experiments we thank Clara 483 \nLinde Schlimbach and Anna Weis. For setting up the ancient DNA clean laboratory in Frankfurt where 484 \nall DNA extractions were done by AS and JR, we thank Leonie Schardt. The authors declare no conflicts 485 \nof interest. 486 \nAI statement 487 \nIn the preparation of this manuscript, we employed artificial intelligence (AI) tools to enhance the 488 \nquality of our work and to optimize our research process. 489 \nR Script Optimization: We utilized an AI -based code optimization to refine our R scripts. The tool 490 \nprovided suggestions for code enhancement, identified potential bugs, and recommended more 491 \nefficient coding practices. This not only improved the performance of our scripts but also ensured the 492 \nreproducibility and reliability of our results. 493 \nWe believe that the use of AI in our research process has improved the quality of our work. However, 494 \nwe stress that the final decisions on manuscript content and the interpretation of the results were 495 \nmade by the authors. 496 \nData availability 497 \nThe datasets generated and analyzed during this study are publicly available. Data and scripts can be 498 \naccessed through the GitHub repository at  https://github.com/Alex-132/seda_enrichment. 499 \nAdditionally, sequence data have been deposited in the European Nucleotide Archive (ENA) under 500 \nthe accession number XXXXX (will be uploaded there). 501 \n 502 \nReferences 503 \n1. Goulletquer, P., Gros, P., Boeuf, G. & Weber, J. The Impacts of Human Activities on Marine 504 \nBiodiversity. in Biodiversity in the Marine Environment 15–20 (Springer Netherlands, 505 \nDordrecht, 2014). doi:10.1007/978-94-017-8566-2_2. 506 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n16 \n2. Panja, P., Kar, T. & Jana, D. K. Impacts of global warming on phytoplankton–zooplankton 507 \ndynamics: a modelling study. Environ Dev Sustain (2024) doi:10.1007/s10668-023-04430-3. 508 \n3. Benoiston, A.-S. et al. The evolution of diatoms and their biogeochemical functions. Phil. 509 \nTrans. R. Soc. B 372, 20160397 (2017). 510 \n4. Johansson, O. N. et al. Skeletonema marinoi as a new genetic model for marine chain-511 \nforming diatoms. Sci Rep 9, 5391 (2019). 512 \n5. Nakweya, G. Southern Ocean phytoplankton changes affect global climate regulation. Nat 513 \nAfrica d44148-023-00242–9 (2023) doi:10.1038/d44148-023-00242-9. 514 \n6. Hjerne, O., Hajdu, S., Larsson, U., Downing, A. S. & Winder, M. Climate Driven Changes in 515 \nTiming, Composition and Magnitude of the Baltic Sea Phytoplankton Spring Bloom. 516 \nFrontiers in Marine Science 6, (2019). 517 \n7. Hattich, G. S. I. et al. Temperature optima of a natural diatom population increases as 518 \nglobal warming proceeds. Nature Climate Change 1–8 (2024). 519 \n8. Rynearson, T. A., Bishop, I. W. & Collins, S. The Population Genetics and Evolutionary 520 \nPotential of Diatoms. in The Molecular Life of Diatoms (eds. Falciatore, A. & Mock, T.) 29–57 521 \n(Springer International Publishing, Cham, 2022). doi:10.1007/978-3-030-92499-7_2. 522 \n9. Godhe, A. & Rynearson, T. The role of intraspecific variation in the ecological and 523 \nevolutionary success of diatoms in changing environments. Phil. Trans. R. Soc. B 372, 524 \n20160399 (2017). 525 \n10. Björck, S. The late Quaternary development of the Baltic Sea basin. in Assessment of 526 \nclimate change for the Baltic Sea Basin 398–407 (Springer, 2008). 527 \n11. Zillén, L., Conley, D. J., Andrén, T., Andrén, E. & Björck, S. Past occurrences of hypoxia in the 528 \nBaltic Sea and the role of climate variability, environmental change and human impact. 529 \nEarth-Science Reviews 91, 77–92 (2008). 530 \n12. Ojaveer, H. et al. Status of Biodiversity in the Baltic Sea. PLoS ONE 5, e12467 (2010). 531 \n13. Reusch, T. B. H. et al. The Baltic Sea as a time machine for the future coastal ocean. Sci. 532 \nAdv. 4, eaar8195 (2018). 533 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n17 \n14. Smol, J. P. The power of the past: using sediments to track the effects of multiple stressors 534 \non lake ecosystems. Freshwater Biology 55, 43–59 (2010). 535 \n15. Ellegaard, M. et al. Dead or alive: sediment DNA archives as tools for tracking aquatic 536 \nevolution and adaptation. Commun Biol 3, 1–11 (2020). 537 \n16. Rengefors, K., Kremp, A., Reusch, T. B. H. & Wood, A. M. Genetic diversity and evolution in 538 \neukaryotic phytoplankton: revelations from population genetic studies. J. Plankton Res. 539 \nplankt;fbw098v1 (2017) doi:10.1093/plankt/fbw098. 540 \n17. Laso-Jadart, R., O’Malley, M., Sykulski, A. M., Ambroise, C. & Madoui, M.-A. Holistic view of 541 \nthe seascape dynamics and environment impact on macro-scale genetic connectivity of 542 \nmarine plankton populations. BMC Ecol Evo 23, 46 (2023). 543 \n18. Brand, L. E. Review of Genetic Variation in Marine Phytoplankton Species and the Ecological 544 \nImplications. Biological Oceanography (1980) doi:10.1080/01965581.1988.1074954. 545 \n19. Salmaso, N. & Tolotti, M. Phytoplankton and anthropogenic changes in pelagic 546 \nenvironments. Hydrobiologia 848, 251–284 (2021). 547 \n20. Myrberg, K. & Soomere, T. The Gulf of Finland, Its Hydrography and Circulation Dynamics. in 548 \nPreventive Methods for Coastal Protection: Towards the Use of Ocean Dynamics for 549 \nPollution Control (eds. Soomere, T. & Quak, E.) 181–222 (Springer International Publishing, 550 \nHeidelberg, 2013). doi:10.1007/978-3-319-00440-2_6. 551 \n21. Leppäranta, M. History and Future of Snow and Sea Ice in the Baltic Sea. in Oxford Research 552 \nEncyclopedia of Climate Science (2023). doi:10.1093/acrefore/9780190228620.013.891. 553 \n22. Salimi, P. A. et al. A review of the diversity and impact of invasive non-native species in 554 \ntropical marine ecosystems. Mar Biodivers Rec 14, 11 (2021). 555 \n23. Wang, Q., Cheng, F., Xue, J., Xiao, N. & Wu, H. Bacterial community composition and 556 \ndiversity in the ballast water of container ships arriving at Yangshan Port, Shanghai, China. 557 \nMar Pollut Bull 160, 111640 (2020). 558 \n24. Andrés, J. et al. Environment and shipping drive environmental DNA beta‐diversity among 559 \ncommercial ports. Molecular Ecology 32, 6696–6709 (2023). 560 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n18 \n25. Wang, Y., Wessels, M., Pedersen, M. W. & Epp, L. S. Spatial distribution of sedimentary DNA 561 \nis taxon-specific and linked to local occurrence at intra-lake scale. Communications Earth 562 \n& Environment 4, 172 (2023). 563 \n26. Carstensen, J. et al. Hypoxia in the Baltic Sea: Biogeochemical Cycles, Benthic Fauna, and 564 \nManagement. AMBIO 43, 26–36 (2014). 565 \n27. Snoeijs-Leijonmalm, P. & Andrén, E. Why is the Baltic Sea so special to live in? in Biological 566 \nOceanography of the Baltic Sea (eds. Snoeijs-Leijonmalm, P., Schubert, H. & Radziejewska, 567 \nT.) 23–84 (Springer Netherlands, Dordrecht, 2017). doi:10.1007/978-94-007-0668-2_2. 568 \n28. Andrén, T. et al. The Development of the Baltic Sea Basin During the Last 130 ka. in The 569 \nBaltic Sea Basin (eds. Harff, J., Björck, S. & Hoth, P.) 75–97 (Springer Berlin Heidelberg, 570 \nBerlin, Heidelberg, 2011). doi:10.1007/978-3-642-17220-5_4. 571 \n29. Kabel, K. et al. Impact of climate change on the Baltic Sea ecosystem over the past 1,000 572 \nyears. Nature Climate Change 2, 871–874 (2012). 573 \n30. Seppä, H., Bjune, A. E., Telford, R. J., Birks, H. J. B. & Veski, S. Last nine-thousand years of 574 \ntemperature variability in Northern Europe. Clim. Past 5, 523–535 (2009). 575 \n31. Epp, L. S., Zimmermann, H. H. & Stoof-Leichsenring, K. R. Sampling and Extraction of 576 \nAncient DNA from Sediments. in Ancient DNA (eds. Shapiro, B. et al.) vol. 1963 31–44 577 \n(Springer New York, New York, NY, 2019). 578 \n32. Schmidt, A. et al. Decoding the Baltic Sea’s past and present: A simple molecular index for 579 \necosystem assessment. Ecological Indicators 166, 112494 (2024). 580 \n33. Gansauge, M.-T., Aximu-Petri, A., Nagel, S. & Meyer, M. Manual and automated preparation 581 \nof single-stranded DNA libraries for the sequencing of DNA from ancient biological remains 582 \nand other sources of highly degraded DNA. Nat Protoc 15, 2279–2300 (2020). 583 \n34. Waldvogel, A.-M. et al. The genomic footprint of climate adaptation in Chironomus riparius. 584 \nMolecular Ecology 27, 1439–1456 (2018). 585 \n35. Bolger, A. M., Lohse, M. & Usadel, B. Trimmomatic: a flexible trimmer for Illumina sequence 586 \ndata. Bioinformatics 30, 2114–2120 (2014). 587 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n19 \n36. Andrews, S. FastQC: A quality control analysis tool for high throughput sequencing data. 588 \n(2021). 589 \n37. Ewels, P., Magnusson, M., Lundin, S. & Käller, M. MultiQC: summarize analysis results for 590 \nmultiple tools and samples in a single report. Bioinformatics 32, 3047–3048 (2016). 591 \n38. Li, H. & Durbin, R. Fast and accurate short read alignment with Burrows–Wheeler transform. 592 \nBioinformatics 25, 1754–1760 (2009). 593 \n39. Danecek, P. et al. Twelve years of SAMtools and BCFtools. GigaScience 10, giab008 (2021). 594 \n40. Danecek, P. et al. The variant call format and VCFtools. Bioinformatics 27, 2156–2158 595 \n(2011). 596 \n41. Browning, B. L., Tian, X., Zhou, Y. & Browning, S. R. Fast two-stage phasing of large-scale 597 \nsequence data. The American Journal of Human Genetics 108, 1880–1890 (2021). 598 \n42. Jónsson, H., Ginolhac, A., Schubert, M., Johnson, P. L. F. & Orlando, L. mapDamage2.0: fast 599 \napproximate Bayesian estimates of ancient DNA damage parameters. Bioinformatics 29, 600 \n1682–1684 (2013). 601 \n43. Wickham, H. et al. Welcome to the Tidyverse. JOSS 4, 1686 (2019). 602 \n44. Kassambara, A. ggpubr: ‘ggplot2’ Based Publication Ready Plots. (2023). 603 \n45. Dunnington, D. W., Libera, N., Kurek, J., Spooner, I. S. & Gagnon, G. A. tidypaleo : 604 \nVisualizing Paleoenvironmental Archives Using ggplot2. J. Stat. Soft. 101, (2022). 605 \n46. Oksanen, J. et al. vegan: Community Ecology Package. R package version 2.5-6. 2019. 606 \nCommunity Ecol 8, 732–740 (2020). 607 \n47. Wood, S. N. Fast Stable Restricted Maximum Likelihood and Marginal Likelihood Estimation 608 \nof Semiparametric Generalized Linear Models. Journal of the Royal Statistical Society Series 609 \nB: Statistical Methodology 73, 3–36 (2011). 610 \n48. Kontny, B. Maritime contacts across the Baltic Sea during the Roman and Migration Periods 611 \n(1st-7th centuries AD) in the light of archaeological sources from the central-European 612 \nperspective. (2023). 613 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint \n\n20 \n49. Mägi, M. In Austrvegr: The Role of the Eastern Baltic in Viking Age Communication across the 614 \nBaltic Sea. Brill; Leiden. vol. 84 (2018). 615 \n50. Rytkönen, J., Siitonen, L., Riipi, T., Sassi, J. & Sukselainen, J. Statistical Analyses of the Baltic 616 \nMaritime Traffic. (VTT Technical Research Centre of Finland, 2002). 617 \n 618 \n.CC-BY 4.0 International licensemade available under a \n(which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in perpetuity. It is \nThe copyright holder for this preprintthis version posted March 13, 2025. ; https://doi.org/10.1101/2025.03.10.642313doi: bioRxiv preprint","source_license":"CC-BY-4.0","license_restricted":false}