Using 12S rRNA and mitochondrial cytochrome c oxidase subunit I (COI) gene metabarcoding, we examined diet composition, prey resources, and saxitoxin (STX) and parasite exposure in Steller sea lions (SSL) in the Gulf of Alaska in July and October 2022–2023 using 251 scat samples and 45 seawater samples. SSL scat and water were sampled from a broad geographic area in the the Gulf of Alaska, from Southeast Alaska through Kodiak Island. The COI gene region was amplified with the Leray-XT primers and the 12S gene region was amplified with the MiFish-U primers. Data available here include: 1. SSL diet and parasite and saxitoxin exposure information from scat samples, including metadata (for each scat sample: collection location and date, host sex determined via molecular methods, saxitoxin level), unfiltered results of taxonomic groups detected in each scat sample for each gene region, and filtered results of prey species detected which were used in our published study in Frontiers in Marine Science. And 2. Taxonomic diversity of the nearshore waters of the Gulf Alaska for initial examination of geographic variation in the nearshore prey field for SSL in July (informed by the 12S gene region metabarcoding) for comparison to prey composition determined from SSL scat.  Information about other taxonomic groups of interest in water samples included geographic variation in saxitoxin-producing Alexandrium spp., possible from COI gene metabarcoding which sampled well lower trophic levels.  Also included are 3. Information concerning reference databases used for bioinformatics for both gene regions and 4. raw sequence data for scat and water samples for each gene region (fastq files).

From manuscript:

2.1 Sample collection and preparation

2.1.1 SSL scat samples

Both in the field and laboratory, a rigorous quality control protocol was followed to prevent contamination of samples from exogenous DNA (King et al., 2008; USFWS, 2022).  In July and October 2022 and 2023, 251 scats were collected at 8 sites in 3 regions in the Gulf of Alaska: Cook Inlet, Kodiak Island, and PWS (Figure 1, Table 1).  Entire fresh scats (regardless of size or characteristics) were collected using clean, disposable gloves and spoon for each sample, placed in Ziploc^® bags and stored at -20°C within 6 hr of collection (Thomas et al., 2022).

Sample preparation followed protocols detailed for harbor seal diet metabarcoding studies (Thomas, 2015; Thomas et al., 2022).  To prepare samples for laboratory analysis, each frozen sample was partially thawed on a sterilized workspace and a ~4mL subsample was transferred to a 5mL cryovial using sterilize instruments and immediately stored at -80°C for STX analysis.  The remaining sample was submerged in 95% molecular-grade ethanol in 500mL or 1L Histoplex^® jars lined with 1gal-sized paint strainers made of 200μm nylon mesh and then, after fully thawed, manually homogenized with a clean wooden tongue depressor by gentle stirring for ~3–5min per sample.  The paint strainer with hard-part and non-scat items (e.g., rocks or debris if present) was removed from the jar, placed in a Ziploc^® bag, and archived at -20°C for future hard-part analysis.  The scat matrix ethanol mixture was stored at -20°C for 24 hr to allow the matrix to settle out of solution.  A clean, disposable pipette was used to transfer two subsamples of the matrix: 0.5mL to a 2mL cryovial (for DNA extraction) and 1mL to a 5mL cryovial (for archiving), which were topped off with 1.5mL and 3mL of fresh 95% molecular-grade ethanol (e.g., 1:3 ratio of matrix to ethanol), respectively, and stored at -80°C.  This procedure was used to ensure sample homogenization to maximize detection of diet items and reduce sample variance, as food items are not distributed equally within scat samples (Deagle et al., 2005).

2.1.2 Seawater samples to assess SSL prey field

From 24 June–23 July 2023, seawater samples were collected via a 1.7L Niskin bottle (Model 1010, General Oceanics, Inc., Miami, FL, USA) at a depth of 10m during Gulf-wide SSL population surveys.  Our goal for this sampling was to acquire initial data as a pilot study to determine if a fast, efficient seawater sampling protocol would provide useful data on nearshore fish assemblages as an efficient addition to our regular, annual SSL surveys.  A standard depth of 10m was chosen for the following reasons:  (1) Dives of adult female SSL from the western population were relatively shallow with median or average dive depths of 21–28m during summer and winter (Merrick and Loughlin, 1997; Lander et al., 2020); up to 54% of dives of adult females were between 4–10m and >90% were <50m (Merrick and Loughlin, 1997).  (2) This depth allowed us to sample below the water’s surface but was shallow enough to allow quick, hand-deployment of the Niskin bottle; hand deploying to deeper depths, especially in strong currents, could be prohibitive. (3) We desired information on the nearshore environment and did not assume foraging areas or foraging depths of SSL; a rigorous, comprehensive study evaluating the utility of eDNA to study vertebrate diversity found equivalent information across 10–80m sampling depths suggesting vertical mixing may provide information on fish diversity from throughout the upper water column, even from surface water collections (Closek et al., 2019).

 Three 1.7L replicates per site were collected at 15 sites, at least 0.5km distant from SSL aggregations, from SEAK to Kodiak Island (Figure 1).  Seawater was filtered directly from the Niskin bottle in the field using a Smith-Root eDNA Citizen Scientist Sampler and 5µm self-preserving polyethersulfone filters (Smith-Root, Inc., Vancouver, WA, USA).  Five negative controls of distilled water were processed in triplicate at weekly intervals throughout the sampling period to test for contamination of the field sampling system.  Filters were stored in a cool, dry, and dark location until subsequent DNA extraction within 3 months post sample collection.

2.2 DNA metabarcoding

2.2.1 DNA extraction

All laboratory work for scat and seawater samples was performed at the Washington Department of Fish and Wildlife’s Molecular Genetics Laboratory (WDFW-MGL) in AirClean 600 Workstations equipped with HEPA filtered air and UV light irradiation.  Work surfaces and equipment were sterilized with 10% bleach and exposed to UV-C light for one hour to neutralize exogenous DNA.  Seawater filter DNA was extracted using DNeasy^TM Blood & Tissue Kits (Qiagen, Hilden, Germany) following the methods of Pilliod et al. (2013) with one extraction negative control processed per batch of 47 filters.  SSL scat DNA was extracted using QIAamp^TM Fast DNA Stool Mini Kits (Qiagen) following a customized protocol for pinniped scat (Deagle et al., 2005) with one extraction negative control processed per batch of 23 scats.

Library preparation -  We assessed prey diversity using two genes, mitochondrial cytochrome c oxidase subunit I (COI) and 12S rRNA MiFish (12S) genes, to provide better biodiversity coverage.  The 12S rRNA MiFish-U primer was optimal for bony fishes (Actinopterygii) and informative of cartilaginous fishes (Chondrichthyes; Miya et al., 2015), whereas the COI primer amplifies all metazoans, including cephalopods (Arfin et al., 2023) which are potentially important SSL prey (Pitcher, 1981).  For pinniped scat diet studies, the COI gene (although a different fragment than used by in our study) has been used in particular to distinguish salmon species (Lewis, 2022; Thomas et al., 2022; Trzcinski et al., 2024).  The COI and 12S genes were amplified in separate reactions.  The COI gene region was amplified with the Leray-XT primers (Wangensteen et al., 2018) and the 12S gene region was amplified with the MiFish-U primers (Miya et al., 2015).  The Leray-XT primers amplified a 313bp gene fragment and included the forward primer mlCOIintF-XT 5′-GGWACWRGWTGRACWITITAYCCYCC-3′ (Wangensteen et al., 2018), modified from the original mlCOIintF primer developed by Leray et al. (2013), and the reverse primer jgHCO2198 5′-TAIACYTCIGGRTGICCRAARAAYCA-3′ (Geller et al., 2013).  The MiFish-U primers amplified a 170bp gene fragment and included the forward primer 5′-GTCGGTAAAACTCGTGCCAGC-3′ and the reverse primer 5′-CATAGTGGGGTATCTAATCCCAGTTTG-3′ (Miya et al., 2015).  

PCRs were performed in 30μL volumes using the Multiplex PCR Kit (Qiagen).  COI reactions contained 15μL (1X) of multiplex master mix, 0.8μM of each primer and 2μL of template DNA with thermal cycling conditions as follows: 95°C for 10min, 35 cycles of 95°C for 60s, 50°C for 60s, and 72°C for 60s, and a final extension of 72°C for 5min. 12S reactions contained 15μL (1X) of multiplex master mix, 2μM of each primer and 1μL of template DNA with thermal cycling conditions as follows: 95°C for 10min, 14 touchdown cycles of 94°C for 30s, 69.5-50°C for 30s, and 72°C for 90s, 25 cycles of 94°C for 30s, 50°C for 30s, and 72°C for 45s, and a final extension of 72°C for 10min. A PCR negative and positive control were included on each 96-well PCR plate.  The negative control consisted of sterile molecular grade water in lieu of template DNA and kangaroo DNA was used as a positive amplification control.

PCR products were size selected using Mag-Bind^® TotalPure NGS (Omega Bio-tek, Inc., Norcross, GA, USA) beads.  A ratio of beads to product of 0.8X was used for the COI amplicon and 1.2X was used for the 12S amplicon.  Sample amplicons were then indexed with Nextera DNA unique dual (UD) indexes (IDT^® for Illumina^®), normalized using the SequalPrep^™ Normalization Plate Kit (Invitrogen) and each 96-well plate was subsequently pooled.  Plate libraries were bead cleaned at a 0.8X bead ratio, quantified with a Qubit^® fluorometer (Life Technologies), and pooled by amplicon.  Amplicon libraries were quantified and normalized to 4nM prior to loading at a 1:1 12S to COI ratio on the NextSeq^™ 1000 platform with 25% PhiX Control v3 Library (Illumina, Inc., San Diego, CA, USA).  Seawater filters were sequenced using the NextSeq^™ 1000 P1 (300 cycles) Reagent Kit for single end reads and scats were sequenced using the NextSeq^™ 1000 P1 (600 cycles) Reagent Kit for paired end reads (Illumina, Inc.).

2.2.2 Bioinformatics

Amplicon sequence data were analyzed separately (COI and 12S) using either stand-alone QIIME 2 (Bolyen et al., 2019) and DADA2 (Callahan et al., 2016) or with Tourmaline (https://github.com/aomlomics/tourmaline), a Snakemake pipeline that wraps QIIME 2 and DADA2, providing reproducible metabarcoding analysis.  Adapters and primers were trimmed from demultiplexed FASTQ reads using Cutadapt (Martin, 2011).  The program DADA2 was used to quality filter reads.  Reads were truncated to a common length (210bp for the COI amplicon and 160bp for the 12S amplicon) with a maximum number of expected errors=2, chimeras were removed (consensus method), and amplicon sequence variants (ASVs) were exported.  

To assign COI taxonomy, a custom reference database was generated that included existing sequence data for published Alaska SSL prey observed in the Resource Assessment and Conservation Engineering (RACE) Gulf of Alaska summer biennial bottom trawl surveys (https://www.fisheries.noaa.gov/foss/f?p=215%3A28) from 2015–2023, the MIDORI database (Machida et al., 2017), and all mitochondrial COI sequences in the National Center for Biotechnology Information (NCBI) nucleotide database.  To assign 12S taxonomy, a rCRUX (Curd et al., 2024) generated reference database was used that included all 12S sequences in the NCBI nucleotide database, an additional custom database comprised of all Actinopterygii (bony fishes) mitogenomes (Gold et al., 2023), and sequences available for published SSL prey and for fishes and cephalopods in the RACE surveys.  If voucher specimen sequence data were publicly available, they were included in the reference library.  If sequences were identical, assignments were produced at a higher taxonomic resolution.

Global taxonomic alignments between query and reference sequences were performed using the VSEARCH consensus taxonomy classifier (Rognes et al., 2016) and matches with ≥ 97% identity were retained.  Read counts for each operational taxonomic unit (OTU) per sample were tallied.

2.3 Sex identification using SSL scat samples

DNA extract from 250 scat samples was provided to the Alaska Department of Fish and Game's Gene Conservation Laboratory to determine SSL sex identification by adapting existing protocols (Lewis, 2022; Gard et al., 2024).  Briefly, sex-specific oligonucleotides developed for SSL were used to assess the presence of X- and Y-chromosomally linked genes (ZFX and SRY) in the scat DNA.  DNA was amplified on a GeneAmp™ PCR System 9700 with a Multiplex PCR Kit (Qiagen), 0.2μM of each primer, and 4μL of DNA using Gene Conservation Laboratory’s standard 14 cycle preamplification PCR protocol with thermal cycling conditions as follows: 95°C for 15min, 14 cycles of 95°C for 15s and 60°C for 4min.  The individual primers and probes from the two qPCR assays were then combined into a single 80X sex ID genotyping assay (16μM each probe, 72μM each primer).  Final concentrations of the probes and primers in the PCR reaction were 0.2µM and 0.9µM.  Genotypes were collected on an Applied Biosystems™ QuantStudio 12K Flex using their 2X TaqMan™ GTXpress™ Master Mix, 1X genotyping assay, and 2.5μL of 0.1X preamplification PCR product according to the PCR protocol of Schwarz et al. (2018).  Genotypes were scored at 42 cycles.

2.4 Saxitoxin (STX) in SSL scat samples

 In addition to scat samples from the current study, 8 samples from the Kenai Peninsula collected by the Alaska SeaLife Center (Jan–May 2017–2019) and 12 samples from SEAK (July 2016) were included.  STX in scat samples was measured by the Northwest Fisheries Science Center's Wildlife Algal-Toxin Research and Response Network (WARRN-West) laboratory, following the protocol of Lefebvre et al. (2022).  Briefly, the ~4mL raw scat samples provided to the lab were thawed slowly in a small cooler, stirred thoroughly, and ~1mL per sample was aliquoted into 14mL polypropylene screw-cap tubes (Falcon-BD).  For STX extraction, 50% methanol was added to samples in a 1:4 wt/v ratio (1 part sample, 3 parts solvent) and samples were vortexed.  Samples were homogenized for 60s using an Omni GLH 850 homogenizer and the homogenized sample was then centrifuged at 3,082xg (Jouan CR3i centrifuge) for 20min at 4°C.  The supernatant was added to a 0.22μm Durapore™ membrane filter (Millipore Ultrafree-MC centrifugal concentration device) and filtered in a desk-top microcentrifuge (AccuSpin Micro 17, Fisher Scientific) for 3min at 12,000rpm.  Sample extracts were stored at 4°C until analysis by enzyme-linked immunosorbent assay (ELISA).

STX in samples was quantified using a commercially-available Abraxis saxitoxin ELISA kit (PN 52255B: Gold Standard Diagnostics, Horsham, PA).  Standards solutions provided in these kits were loaded along with the samples on plates and a standard curve was derived for each plate.  Manufacturer instructions were followed with modifications for matrix effects for marine mammal feces determined by Hendrix et al. (2021) with sample dilution of 1:50 filtered sample extract:sample diluent solution provided with the ELISA kit.  The ELISA kit was designed to measure STX with some limited cross-reactivity to several other paralytic shellfish poison (PSP) toxins.  Consequently, all PSP concentrations are listed as STX equivalents and may underestimate the presence of other congeners.  

2.5 Statistics

2.5.1 Diet and prey data summaries

We included taxa for bony and cartilaginous fishes for the 12S data and for these fishes, cephalopods, and other potentially interesting taxa (e.g., parasitic worms) for the COI data.  For both metabarcoding datasets for scat samples, we performed the following minimum sequence threshold filtering steps prior to statistical analyses: (1) samples with <20 total prey assigned reads were removed (Trzcinski et al., 2024), (2) OTUs were assigned into prey groups at the species level or at the genus level for groups in which significant assigned reads were only to genus level (Deagle et al., 2019), (3) global (summed across all samples) total assigned reads per prey group were calculated and groups with <100 total assigned reads were removed, and (4) prey groups <0.01% of the per sample proportion of total prey assigned reads were removed to normalize detection rates among samples which may vary in read depth (Pornon et al., 2016; Richardson et al., 2019) to assess region*season variation in diet. For seawater samples, all assigned fish reads were used without filtering data, and assigned reads for the three replicates per site were summed to provide sufficient sample size.