Sampling location and sediment coring
Samples were collected during IODP Exp. 382 ‘Iceberg Alley and Subantarctic Ice and Ocean Dynamics’ on-board RV Joides Resolution between 20 March and 20 May 2019. Specifically, we collected samples at Site U1534 (Falkland Plateau, 606 m water depth), U1536 (Dove Basin, Scotia Sea, 3220 m water depth), and Site U1538 (Pirie Basin, Scotia Sea, 3130 m water depth) (Fig. 1). Site U1534 is located at the Subantarctic Front on a contourite drift at the northern limit of the Scotia Sea. This setting is ideal to study the poorly understood role of Antarctic Intermediate Water (AAIC) and its impact on the Atlantic Meridional Overturning Circulation (AMOC) along the so-called ‘cold water route’ that connects to the Pacific Ocean through the Drake Passage, as opposed to the ‘warm water route’ that connects to the Indian Ocean via the Agulhas Current42. Sites U1536 and U1538 are located in the southern and central Scotia Sea, respectively, and were drilled to study the Neogene flux of icebergs through ‘Iceberg Alley’, the main pathway along which icebergs calved from the margin of the AIS travel as they move equatorward into the warmer waters of the Antarctic Circumpolar Current (ACC)23. sedaDNA samples collected at Site U1534 were from Hole C, at Site U1536 from Hole B, and at Site U1538 from Holes C and D (Table 1), and in the following we refer to site names only. IODP Expedition proposals undergo a rigorous environmental protection and safety review, which is approved by the IODP’s Environmental Protection and Safety Panel (EPSP) and/or the Safety Panel. The same procedure was applied to IODP Exp. 382 and approval was provided by the EPSP. Sediment samples for sedaDNA analyses were imported to Australia under Import Permit number 0002658554 provided by the Australian Government Department for Agriculture and Water Resources (date of issue: 19 September 2018), and were stored and extracted at a quarantine approved facility (AA Site No. S1253, Australian Centre for Ancient DNA). No ethical approval was required for this study.
Sample age determination
Age control for Site U1534 is based on tuning of benthic foraminifera δ18O to the LR04 stack43. Wherever present specimens of Uvigerina bifurcata were picked from samples at 10 cm intervals. During warmer periods when U. bifurcata was not present, Melonis affinis and/or Hoeglundina elegans were analysed. Sedimentation rates over the intervals sampled for sedaDNA typically range between 6 and 30 cm/kyr, with rates exceeding 100 cm/kyr during the Last Glacial Maximum ~20,000 years ago (20 ka). For our deepest sample, U1534C-10H-6_115cm (90.95 mbsf), we only have biostratigraphically assigned ages available (shipboard data), which date this sample as early Pleistocene (~2.5–0.7 million years ago, Ma44).
Low-resolution age control for both Sites U1536 and U1538 was established using shipboard magneto- and biostratigraphy21,23. Average sedimentation rates are ~10 cm/kyr for Site U1536, with elevated values (up to 20 cm/kyr) in the upper ~80 mbsf (the last ~400 ka). Site U1538 average sedimentation rates are twice as high, averaging ~20 cm/kyr. Especially in the upper ~430 mbsf (the last 1.8 Ma), rates are up to 40 cm/kyr. Higher resolution age models are based on dust climate couplings, correlating sedimentary dust proxy records such as magnetic susceptibility and sedimentary Ca and Fe records to ice-core dust proxy records over the last 800 ka45 and to a benthic isotopic stack26 before that. These age models were established for Site U1537 (adjacent to Site U1536) and provide orbital to millennial scale resolution. For this study we correlated sedimentary cycles of Sites U1536 and U1538 to U1537 to achieve similar resolution and to be able to determine if a sample originates from a glacial or interglacial period (Table 1).
Sampling of sedaDNA
A detailed description of sedaDNA sampling methods can be found in ref. 24. In brief, we used advanced piston coring (APC) to acquire sediment cores, which recovers the least disturbed sediments46,47,48 and is thus the preferred technique for sedaDNA sampling. All samples were taken on the ship’s ‘catwalk’, where, once the core was on deck, the core liners were wiped clean twice (3% sodium hypochlorite, ‘bleach’) at each cutting point. Core cutting tools were sterilised before each cut (3% bleach and 80% ethanol) of the core in 1 m sections. The outer ~3 mm of surface material were removed from the bottom of each core section to be sampled, using sterilised scrapers (~4 cm wide; bleach and ethanol treated). A cylindrical sample was taken from the core centre using a sterile (autoclaved) 10 mL cut-tip syringe, providing ~5 cm3 of sediment material. The syringe was placed in a sterile plastic bag (Whirl-Pak) and immediately frozen at −80 °C. The mudline (sediment/seawater interface) was transferred from the core liner into a sterile bucket (3% bleach treated), and 10 mL sample was retained in a sterile 15 mL centrifuge tube (Falcon) and frozen at −80 °C. Samples were collected at various depth intervals depending on the site to span the Holocene up to ~1 million years (Table 1). This lower depth/age limit was determined by switching coring system from APC to the extended core barrel (XCB) system.
To test for potential airborne contamination, at least one air control was taken during the sedaDNA sampling process per site. For this, an empty syringe was held for a few seconds in the sampling area and then transferred into a sterile plastic bag and frozen at −80 °C. The air controls were processed, sequenced and analysed alongside the sediment samples.
Contamination control using perfluoromethyldecalin tracers
As part of the APC process, drill fluid (basically, seawater) is pumped into the borehole to trigger the hydraulic coring system, therefore, the potential for contamination exists due to drill fluid making contact with the core liner. To assess the latter, we added the non-toxic chemical tracer perfluoromethyldecalin (PFMD) to the drill fluid at a rate of ~0.55 mL min−1 for cores collected at Sites U1534 and U153649. As we found that PFMD concentrations were very low at these sites (Results section), the infusion rate was doubled prior to sedaDNA sampling at Site U1538 to ensure low PFMD concentrations represent low contamination and not delivery failure of PFMD to the core. At each sedaDNA sampling depth, one PFMD sample was taken from the periphery of the core (prior to scraping, to test whether drill fluid reached the core pipe), and one next to the sedaDNA sample in the centre of the core (after scraping, to minimise differences to the sedaDNA sample, and testing if drill fluid had reached the core centre). We transferred ~3 cm3 of sediment using a disposable, autoclaved 5 mL cut-tip syringe into a 20 mL headspace vial with metal caps and Teflon seals. We also collected a sample of the tracer-infused drill fluid at each site, by transferring ~10 mL of the fluid collected at the injection pipe on the rig floor via a sterile plastic bottle into a 15 mL centrifuge tube (inside a sterile plastic bag) and freezing it at −80 °C. These drill fluid controls were processed and analysed in the same way as the sedaDNA samples including sequencing. Samples were analysed using gas chromatography (GC-µECD; Hewlett-Packard 6890).
A detailed description of the PFMD GC measurements is provided in ref. 24. Briefly, PFMD measurements were undertaken in batches per site for U1534, U1536 and U1538. This included the analyses of PFMD samples collected at two additional holes at these sites, U1534D and U1536C, from which we also collected sedaDNA samples but that are not part of this study. PFMD is categorised as the stereoisomers of PFMD (C11F20), which add up to 87-88% (and with the remaining 12% being additional perfluoro compounds unable to be separated by the manufacturer). We exclusively refer to the first and measurable PFMD category, calibrating for the 88% in bottle concentrations during concentration calculations. Each GC analysis run included the measurement of duplicate blanks and duplicate PFMD standards. Due to a large sample number, PFMD at Site U1538 was measured in three separate runs, with the first and last run including triplicate blank and triplicate PFMD standards (duplicates in the second run), and the last run also containing a drill-fluid sample. To blank-correct PFMD concentrations, we subtracted the average PFMD concentration of all blanks per run from PFMD measurements in that run. To determine the detection limit of PFMD, we used three times the standard deviation of the average blank PFMD values per run; due to all blank values for the U1538 runs being 0, we used three times the standard deviation of the lowest PFMD standard for this site in this calculation. This provided us with a PFMD detection limit of 0.2338 ng mL−1. Any PFMD measurements of samples below this limit were rejected.
sedaDNA extractions and metagenomic library preparations
A total of 80 sedaDNA extracts and metagenomic shotgun libraries (Table 1) were prepared following8,10. For the sedaDNA extractions, we randomised our samples and controls and extracted sedaDNA in batches of 16 extracts/libraries at a time, with each batch including at least one air control and one extraction blank control (EBC), and the last batch including mudline and PFMD samples to avoid contamination of the sedaDNA samples. In brief, we used 20 µL sedaDNA extracts in a repair reaction (using T4 DNA polymerase, New England Biolabs, USA; 15 min, 25 °C), then purified the sedaDNA (MinElute Reaction Cleanup Kit, Qiagen, Germany), ligated adaptors (T4 DNA ligase, Fermentas, USA, where truncated Illumina-adaptor sequences containing two unique 7 base-pair (bp) barcodes were attached to the double-stranded DNA; 60 min, 22 °C), purified the sedaDNA again (MinElute Reaction Cleanup Kit, Qiagen), and then added a fill-in reaction with adaptor sequences (Bst DNA polymerase, New England Biolabs, USA; 30 min, 37 °C, with polymerase deactivation for 10 min, 80 °C). We amplified the barcoded libraries using IS7/IS8 primers50 (8 replicates per sample, where each replicate was a 25 µL reaction containing 3 µL DNA template; using 22 cycles), purified (AxyPrep magnetic beads, Axygen Biosciences, USA; 1:1.8 library:beads) and quantified them (Qubit dsDNA HS Assay, Invitrogen, Molecular Probes, USA). We amplified the libraries (8 replicates per sample, 13 amplification cycles) using IS4 and GAII Indexing Primers50, purified (AxyPrep magnetic beads, at a ratio of 1:1.1 library:beads), quantified and quality-checked using Qubit (dsDNA HS Assay, Invitrogen, USA) and TapeStation (Agilent Technologies, USA). We combined the libraries into an equimolar pool (volume of 68 µL in total), diluted this pool with nuclease-free H2O to 100 µL, and performed a ‘reverse’ AxyPrep clean-up to retain only the small DNA fragments typical for ancient DNA (≤ 500 bp; initial library:beads ratio of 1:0.6, followed by 1:1.1, and double-eluted in 30 µL nuclease-free H2O8,51). We added one more AxyPrep clean-up to remove primer-dimer (library:beads ratio of 1:1.05) and checked sedaDNA quantity and quality via TapeStation and qPCR (QuantStudio, Applied Biosystems, USA). The libraries sequenced at the Garvan Institute for Medical Research, Sydney, Australia (Illumina NovaSeq 2 × 100 bp).
sedaDNA data processing
The sequencing data was processed and filtered as described in detail in refs. 8, 10. Briefly, data filtering involved the removal of sequences <25 bp (AdapterRemoval v. 2.1.7-foss-2016a, Schubert et al., 2016), removal of low-complexity (Komplexity software52,–threshold 0.55) and duplicate reads (‘dedupe’ tool in BBMap v.37.36.), and quality control after each step (FastQC v.0.11.4, Babraham Bioinformatics; MultiQC v.1.0.dev053).
Using MALT (version 0.4.0; semiglobal alignment54), we ran our data against different curated reference databases that are commonly used to identify marine eukaryotes, including the SILVA small (16S/18S, SSU) and large subunit (23S/28S, LSU) ribosomal RNA database (version NR 132, https://www.arb-silva.de/), as well as a combined SSU + LSU database. For the latter, we merged the initial SSU and LSU databases (using the ‘cat’ command) into one database prior to the alignments, with the intention to maximise taxonomic resolution and number of reads (as has been shown before by using both SILVA SSU and LSU7, see Supplementary Information and Supplementary Data 2). To further maximise the number of reads (to enable robust downstream and sedaDNA damage analysis, see below), we avoided normalising (rarefying) our SSU/LSU data and worked with relative abundances. For an estimate of total abundance of phytoplankton, we also ran our data against a recently developed database for the single-copy photosynthetic gene psbO, which is present in both prokaryotes and eukaryotes, mainly in one copy per genome55. The latter was initially performed with non-subsampled (non-rarefied) data to be able to determine an adequate subsampling depth56 (representative of the diversity in the data) for subsequent quantitative analyses (excluding potential artifacts due to differences in library sizes). We determined 1.1 Mio reads as an adequate subsampling depth (lowest number of reads in libraries that contained psbO, see Supplementary Information), and subsampled all samples to this depth for quantitative psbO analyses. The resulting.blastn files (SSU + LSU, psbO) were converted to .rma6 files using the Blast2RMA tool in MEGAN (version 6_18_9), with the default settings except for a minimum support percent of zero (‘off’) and a minimum percent identity of 95%.
Subtractive filtering (i.e., subtracting reads for species identified in EBCs, air, and PFMD controls from samples) was conducted with MEGAN CE57 (version 6.21.12). The filtered read counts per taxon (eukaryotes only and phylum-level for SSU and LSU datasets, as well as species level for diatoms to investigate this class in more detail, see Results; and species level (all species) for the psbO dataset) and site were exported for downstream analysis. All taxa determined in the controls post-SSU + LSU alignments are listed in Supplementary Data 6; no contaminants were detected in the psbO dataset. We used the default settings for analysing the taxonomic data in MEGAN CE, i.e., following the NCBI taxonomy and using a LCA (lowest common ancestor) approach, where every read is assigned to a taxon and if a read has significant matches to two different taxa a and b, where a is an ancestor of b, then only the more specific match to b will be used57. Where a taxon was determined multiple times with different reference-groups (e.g., ‘Stramenopiles’, ‘Stramenopiles incertae sedis’ ‘Stramenopiles environmental sample’, Stramenopiles unclassified’) then they were grouped (e.g., ‘Stramenopiles’). For visualisation purposes of the large SSU-LSU phylum-level dataset, we kept all taxa that were abundant with at least 1% on average across all samples as separate taxa (combining the two fungal groups Ascomycota and Basidiomycota), while taxa less abundant than that were grouped (‘Eukaryota.rare’; Results section, Fig. 3). To facilitate visualisation of the SSU-LSU derived diatom dataset, we grouped broader taxa (i.e., ‘Bacillariophyta’, ‘Bacillariophyceae’, ‘Bacillariophycidae’, ‘Bacillariaceae’ as ‘Bacillariophyta’, and ‘Coscinodiscophyceae”, Coscinodiscophycidae’, ‘Coscinodiscales’, ‘Coscinodiscaceae’ as ‘Coscinodiscophyceae’ (Results section, Fig. 5).
sedaDNA damage analysis
As a means of authenticating our sedaDNA, we tested sedaDNA damage using the ‘MALTExtract’ and ‘Postprocessing’ tools of the HOPS v0.33-2 pipeline25. We used the same configurations as in Armbrecht et al. (2021), using the taxalist ‘Eukaryota’ (i.e., only specifying the term ‘Eukaryota’, which captures all eukaryotes) for the SSU + LSU dataset, and a taxalist containing all taxa identified via the non-subsampled psbO data, see Supplementary Information) for the psbO dataset. We used the non-subsampled psbO data because the subsampled data provided too few reads (<5025, Results section) for damage analysis. MaltExtract outputs a read summary, with reads are categorised as ancient (showing damage) or default (passing stringent filtering criteria but not showing damage) for the taxalist-specified taxa25. Based on the latter (ancient vs. default reads), the proportion of sedaDNA damage per taxon was determined for each dataset (SSU + LSU and non-subsampled psbO). A summary of reads classified as ancient and default are provided with the Supplementary Data 3 and 4.
Statistics and correlations
To ensure that our new combined SSU + LSU database generates comparable taxonomic composition to SSU and LSU alone, we performed regression and correlation analyses on the average relative abundance per taxon (phylum-level) acquired when using the SSU, the LSU and the combined SSU + LSU database (using the software PAST Software v.4.0358; and with results provided in Supplementary Information Fig. 1, Supplementary Information Table 1).
Lastly, we investigated potential relationships between sedaDNA damage (determined by HOPS post SSU-LSU alignment), relative abundances of eukaryotes post SSU + LSU alignment (phylum-level) and downcore temperature and porewater geochemical parameters, which were measured as part of the shipboard sampling during IODP Exp. 382. Downhole formation temperatures were calculated based on the temperature gradient obtained at each Site with the Advanced Piston Corer Temperature Tool (APCT)24. Key geochemical parameters for the assessment of organic matter degradation rates (ammonium, alkalinity, phosphate, sulfate, mMol L−1), and total silicon (µM), were routinely determined for each site for the characterisation regarding the intensity of redox zonation in each core24. To investigate relationships between taxonomic composition and cold and warm climate phases, we added benthic δ18O data from26 corresponding to the ages assigned to our samples to this correlation analysis. Mudline samples, as well as the two samples 23162 (U1438 57.55 mbsf) and 23165 (U1538 62.85 mbsf), were removed from the following analyses as no eukaryotes were determined in these samples. Pearson correlation analysis was performed in PAST v.4.0358.
Reporting summary
Further information on research design is available in the Nature Research Reporting Summary linked to this article.
Source: Ecology - nature.com
