We performed initial targeted sampling in Salish Sea waters from November 2023 to February 2025, mostly from the PNNL-Sequim floating dock and pier, and less frequently at four nearby locations: Cline Spit, south of John Wayne Marina near an estuary outflow, open waters North of Port Townsend, and in Discovery Bay (Table 1; Figure 1). Sampling events, including filtering and collection (Figures 2, 3) occurred near-monthly at the PNNL-Sequim dock over the course of a year, with offsite sampling and multiple within-month sampling concentrated during peak summer months. A total of 30 primary sampling events were performed on 25 separate days at ebb, flood, and slack tides, with the majority occurring during flood tides. Water levels, temperature, salinity, and collection dates are shown in Figure 4.

Samples were collected using two pilot methods, either 1) time-discrete samples using a horizontal water sampler, Niskin bottle, or direct surface sampling using a sterilized carboy; or 2) time-extended in situ in-line filtrations at the sampling site. The in-situ apparatus consisted of a custom-designed, weighted metal frame attached to flexible tubing that we tested using both peristaltic and vacuum pumps at approximately 550–680 kPa and 130 kPa, respectively (Supplementary Table 1). We configured the tubing as an in-situ mechanism for two- and three-step serial filtrations using in-line filters that could be deployed from height and immediately pump raw seawater directly through two or three sequentially smaller filter sizes to capture target species, with the final filtrate collected into a sterile vessel. Alternatively, we used a direct water sampler (a Niskin bottle, horizontal water sampler, and/or a sterilized carboy) to collect discrete samples followed by filtrations in a laboratory. The in-line system collected samples over extended periods of time (up to 3 hours with especially dense water samples through a three-step filtration), compared to sampling durations of less than one minute using a direct water sampler.

Prior to sampling events, all sample collection vessels (carboys, Niskin bottles, and horizontal water samplers) were sterilized by autoclave or disinfected using a 10% bleach solution. We primed the collection vessel, tubing, and priming filters prior to every sample collection with water from the site, then discarded the priming water away from the collection point to avoid contamination. After collecting a pre-determined target volume, we immediately recorded metadata and sampling event context and assigned a single parent sampling ID to all associated data. Sampling conditions were noted, and additional weather data from the day of sampling were downloaded from the PNNL-Sequim Marine and Coastal Research Laboratory (MCRL) data site (MCRL Data, 2025).

Depending on which fractions were already processed, we then transported samples to our laboratory in Sequim WA to process the remaining filter fractions. All tubing and vessels were rinsed with 70% ethanol and/or autoclaved prior to priming with the appropriate matrix (we used each sequential filtrate to prime the subsequent collection vessels). For most samples, this involved a combination of in-line filtrations (20µm 47 mm nylon filter and a 5 µm 47 mm mixed-cellulose ester (MCE) membrane secured in either a Millipore 47 mm polypropylene nylon filter holder or a Cytiva 47mm polycarbonate in-line filter holder using a peristaltic pump, Supplementary Tables 2, 3) and single vacuum flask filtrations. In-line filtration was followed by mixing in a single-flask and filtration onto a sterile 0.22 µm polyethersulfone (PES) or MCE filter in a vacuum flask to capture the pico fraction. Flow-throughs from the 0.22 µm filter were precipitated using a method adapted from John et al. (2011); iron chloride hexahydrate was added to the filtrate at 1 mL of 10g/L FeCl · 6H₂O solution in deionized water per 20 L of seawater and incubated for 1 hour, and precipitates were captured on a sterile 0.22 µm PES filter via vacuum filtration.

Additional filtration setups we tested included a 0.05 µm filter after 0.22 µm filtration, and various alternate size combinations with the three filter types: nylon, MCE membranes, and PES. Specifically, we tested 47 mm and 90 mm 0.22 or 0.2 µm MCE or PES filters, either as individual, sterile membranes or sterile vacuum units. We tested the 5 µm and 0.22 µm MCE membranes in an in-line filtration setup using either a polycarbonate filter-holder purchased from Avantor Sciences or a polypropylene filter-holder purchased from Millipore Sigma. The long-term use of 47 and 90 mm filters for differing filter fractions was driven by limited commercial availability of each material-pore size-diameter combination and for time-efficiency. A list of the tested filters, in-line filter holders, and other product information can be found in Supplementary Tables 1-4.

In the final sampling scheme, the sampling process included separating environmental samples by serial filtration into four fractions (Figure 2), labelled by their target organism size and the nominal pore sizes of the filters which they passed through and on which they were collected: macro (> 20 µm), micro (large micro, 5-20 µm), pico (small micro, 0.22-5 µm), and viral (< 0.22 µm).​ Further details are provided in Figure 3. Note that filter subsamples for the 20 µm and 5 µm filters were obtained by changing filters as they became clogged (as estimated qualitatively by reduced flow rate), whereas most filter subsamples for the 0.22 µm filters (both the pico and the viral fractions) were obtained by physically cutting the filters into thirds or quarters after completion of filtration (Figure 3).

DNA was extracted from selected filter subsamples using Zymo’s Quick DNA Plant/Seed Miniprep Kit, and three Qiagen kits: DNeasy PowerWater Kit, PowerSoil Pro Kit, or Blood and Tissue Kit (Supplementary Table 4). Following extraction, DNA quality and quantity were assessed using a NanoDrop One spectrophotometer (ThermoFisher Scientific, Waltham, MA, USA) and Qubit4 fluorometer system (ThermoFisher Scientific, Waltham, MA, USA). Samples were then frozen, and a select few were later analyzed on an Agilent 4150 TapeStation for quality assessment (Supplementary Data Figure 1). Samples were opportunistically selected for comparative analyses and were chosen based on sample factors such as number of technical replicates, the filter size, date of sampling and extractions, and tidal stage. While we did not have sufficient sample numbers for controlled comparisons in many cases, we chose the most comparable samples from the larger exploratory set and excluded samples that were too variable.

For direct comparisons of the kit effect on DNA recovery, DNA from select sampling events were extracted using two separate kits (Blood & Tissue and PowerWater) and analyzed on an Agilent 4150 TapeStation using a D5000 ScreenTape and the relevant kit reagents (D5000 buffer and D5000 ladder) following the kit protocol with no modifications (Agilent Technologies, Palo Alto, CA, USA). To understand whether levels of variability due to kits are large compared to other sources of variation, we also make comparisons of DNA recovery without controlling for other known variables, which we refer to as ‘kitchen sink’ comparisons (a mix of conditions and replicates) to clarify the context in which they should be interpreted.

PNNL-Sequim is uniquely positioned for marine environmental research, leveraging its bench-to-bay capabilities through direct access to Sequim Bay, raw seawater access in the facilities, and a unique array of permanently deployed environmental sensors (MCRL Data, 2025). Most sampling events took place from the PNNL-Sequim floating dock at the mouth of Sequim Bay with an extensive range of environmental parameters documented from the MCRL Data, while other locations relied solely on parameters captured from a YSI ProDSS multiparameter water quality meter (YSI Xylem, Yellow Springs OH) (Supplementary Data Table 5). For PNNL-Sequim sampling events, we acquired data post-sampling using targeted scripts (MCRLdata Sandbox, 2025a; MCRLdata Sandbox, 2025b) to process environmental data from an online data repository (MCRL Data, 2025). One sampling event used for a controlled membrane test was executed in PNNL’s marine mesocosm system, which circulates unfiltered seawater through a large (ca. 62,000 L) outdoor tank (Marine and Coastal Research Laboratory, 2025). This mesocosm system was used for a single sampling event due to restricted access to the PNNL-Sequim floating dock during the scheduled sampling time. While the samples from the mesocosm and dock are not directly comparable, the mesocosm system is constantly circulated with raw seawater with a submerged intake a few meters away the floating dock.

We used data collected for tidal stage, water temperature, and salinity by the sensors described above to provide environmental context for sampling events off the PNNL-Sequim dock (Figure 4). Data collected by sensors are sent to the cloud where they undergo automated initial quality control and displayed in near real-time (MCRL Data, 2025) and are provided as part of the joint metadata summary (MCRLdata Sandbox, 2025a; MCRLdata Sandbox, 2025b). Additional quality control combining automated (e.g. de-spiking) and manual visual assessment were then applied, resulting in several data gaps originating from sensor maintenance, sensor errors, or unreliable data (identified during quality control). To fill any gaps, we used the ranger R package to construct Random Forest models (Wright and Ziegler, 2017). First, we gap-filled water levels by building a Random Forest model using water levels (NOAA Tides and Currents, 2025) pulled using the VulnToolkit R package along with hour of day, and day of year as predictors (Hill et al., 2021). The resulting model (R² = 0.96) was then used to gap-fill water levels. Those gap-filled water levels were used (along with hour of day and day of year) to predict both water temperature and salinity measured off the MCRL dock by CTD (conductivity, temperature, and depth sensor), and those models were used to gap-fill both temperature and salinity time-series. Model goodness-of-fit was high for both temperature and salinity models (R² = 0.93 and 0.95, respectively).

DNA subsamples were shipped on ice to Azenta Life Sciences (South Planfield, NJ,USA) for sequencing. Sequencing libraries were prepared using NEBNext Ultra DNA Library Prep Kit (New England Biolabs, Ipswich, MA, USA) according to manufacturer’s instructions. Briefly, genomic DNA was acoustically fragmented using a Covaris S220 instrument (Covaris, Woburn, MA, USA) and fragmented DNA was cleaned and end repaired. Sequencing adapters were ligated to 3’ ends of fragments after adenylation and enriched by limited-cycle PCR. Sequence libraries were validated using a High Sensitivity D1000 ScreenTape on an Agilent TapeStation (Agilent Technologies, Palo Alto, CA, USA) and quantified using a Qubit 3.0 fluorometer (ThermoFisher Scientific, Waltham, MA, USA). Sequence libraries were also quantified using real-time PCR (Applied Biosystems, Carlsbad, CA, USA) libraries from samples that passed all quality control steps and were multiplexed and clustered onto an Illumina NovaSeq flowcell (Illumina, San Diego, CA, USA) and sequenced using 2x150 paired-end configuration. Base calling and image analysis were conducted using the NovaSeqControl Software. Raw sequence data was converted to fastq format and de-multiplexed using Illumina’s bcl2fastq software v2.20, allowing one mismatch for index sequence identification. Sequencing datasets were processed through the sequence processing pipelines described below for evaluation of the efficacy of filter fractionation methods and to refine strategies for further sampling. Raw paired-end reads were sequentially pre-processed by first interleaving paired reads into a single sequence read file preserving read pair information using the reformat functionality of BBTools v37.62 (Bushnell et al., 2017). Interleaved reads were clustered, with optical and PCR duplicate reads removed using clumpify within BBTools. Reads derived from low-quality regions of flowcells were removed using the filterbytile functionality from BBTools. BBTools bbduk was used to sequentially remove Illumina sequencing adapters, bases from the ends of reads below a phred33 quality score of ten, any reads with an overall phred33 quality score lower than ten, and any reads mapping to PhiX spike-in. Read quality was assessed before and after preprocessing using FastQC v0.12.1 (Babraham Institute, 2023).

We generated metagenome contig assemblies using MEGAHIT v1.2.9 (Li et al., 2015) under default settings in paired-end mode. Output contigs were filtered by length as required for downstream analyses using the reformat functionality of BBTools. Assembly qualities were assessed using QUAST v5.3.0 (Gurevich et al., 2013). Putative coding regions and putative protein sequences were predicted within metagenome assemblies filtered to contigs larger than 2500 bases using Metaprokka v1.15.0 (Seemann, 2014).

Raw short read pairs, interleaved pre-processed reads, contig assemblies, and predicted genes and proteins are found in Supplementary Data 3 and were deposited to BioProject PRJNA1336290.

We assigned each sampling event a unique parent ID that was tagged to the relevant environmental data for every filter and DNA subsample, which also had its own unique subsample ID. To track sampling events, we created a standardized field form with nomenclature and categorization based on the National Microbiome Data Collaborative (NMDC) water sample template for data documentation compatibility with the other public repositories (National Microbiome Data Collaborative, 2026). The unique sample identifiers for all parent and subsamples were used to further facilitate consistent data integration. Sequences, processed sequences, analytical results and associated metadata were captured and a) stored for sharing among the team using DataFed federated data storage and Globus-enabled federated IDs to ensure that this data is searchable and well-structured (George et al., 2025; EcoKMER GitHub Software Repository, 2025) and b) harmonized and visualized using the concurrently-developed EcoKMER program (George et al., 2025) with final output stored as a CSV file (Supplementary Data 1). This integrated metadata file combined sequencing summary information with environmental and contextual data, adhering to FAIR (Findable, Accessible, Interoperable, and Reusable) principles (Wilkinson et al., 2016) and ensures that this data is searchable and well-structured.

To assess sequence diversity across environmental metagenomes, we applied Nonpareil (v3.5.5) to interleaved pre-processed reads (Rodriguez-R et al., 2018). Because the default kmer size for Nonpareil is 24, reads shorter than 24 base pairs (bp) were excluded. Nonpareil estimates average community coverage and diversity using kmer redundancy, modeling the rate at which new 24 bp sequences are encountered with increasing sequencing effort. The resulting Nonpareil diversity index (Nd) used does not cleanly reflect microbial diversity, as it reports diversity inferred from identical 24mers sampled, and does not account for variable rates along genomes, or variable genome size, but it has nevertheless been demonstrated to be useful in detecting differences in community composition (Rodriguez-R et al., 2018).

61 DNA subsamples produced in year one were sequenced, filtered, assembled into contigs, and 24mer sequence diversities were calculated. To quantify taxon-specific enrichment of organisms with varying sizes, we mapped metagenomic reads to curated reference genome sets using BBMap (v39.13) (Bushnell, 2014). Five reference genome sets were constructed from high-quality, complete NCBI genomes representing Prochlorococcus (26 reference genomes), Synechococcus (78 reference genomes), their phages (341 Synechococcus phage genomes; 37 Prochlorococcus phage genomes), and Acanthamoeba spp (5 reference genomes; accession numbers of each reference set in Supplementary Data 2). Our main reason to focus on these taxon groups was to identify the relative occurrence of the two picocyanobacterium and cyanophages that might be associated with them. In addition, we separately mapped reads to a mycetezoan amoeba to determine the effectiveness of filtering out large-celled organisms. Each set was indexed using the BBMap indexing utility. Quality-filtered reads from each sample were aligned to each reference genome set using default parameters, and coverage statistics were calculated with the BBMap pileup utility (Supplementary Data 1), which attaches each read to the reference genome in a run that matches it best, and leaves unattributed any reads that do not pass the minimum matching cutoff to any included genome(s). We repeated this analysis using the curated CyanoRak database of Prochlorococcus and Synechococcus reference samples (Supplementary Figure 3), with highly correlated results between our analysis and that using CyanoRak (Garczarek et al., 2021).

Each of the six genome groups listed above were analyzed separately. During preliminary analysis, we observed that the bulk of Acanthamoeba reads were consistently attributed to the five major chromosomes of Dictyostelium discoideum, so we opted to present the sum of plus and minus reads (relative to genome orientation) most preferentially attributed to these five chromosomes as a marker for large-celled amoeboid organisms. Dictyostelium discoideum is nominally a soil organism, but the representation of diverse amoeboid organisms in the genome references was not strong, and these results should not be construed to positively indicate an undue degree of soil contamination in estuarine waters or to identify the actual genus or species source in these aquatic samples—relative amoeba read counts were only used as a marker for removal of large-celled organisms. For other runs, we report the sum of the plus and minus reads for all genomes in the reference set. To normalize differences in sequencing depth, total read counts per sample were used to compute reads per million (RPM). For comparison, we ran Kraken/Bracken (Lu et al., 2022) using version 2.0.7-beta of the Kraken database using default parameters, outputting the total reads that were attributable to the groups Prochlorococcus, Synechococcus, and Dictyoselium discoideum. The output is reported in Supplementary Data 1. Read counts using Bracken are generally lower than for BBMap, likely due to use of a smaller reference database and differences in default parameters.

Available pumps (see Methods and Supplementary Table 1), including eDNA and peristaltic models, did not reliably generate enough power to overcome sampling challenges at height from piers or vessels, so direct surface sampling was preferred over in-situ filtration. Biomass accumulation quickly clogged the 0.22 µm inline filter, which inhibited consistent flow and led to eventual removal of 0.22 µm inline filters from our protocol in favor of flask filtration of the mixed 5 µm filtrate in the tertiary filtration step (Figures 2, 3).

Millipore’s 47 mm polypropylene nylon in-line filter holders improved flow rates over Cytiva’s 47 mm polycarbonate filter holders, reducing macro and micro filtration time by over 50% and reducing sample loss due to inhibitive filtration times and insufficient gasket seals in alternative holders. The 30 primary sampling events were ultimately segregated by filtration and technical replicates into 294 filter subsamples and 424 DNA subsamples, with filters cut and used for replicate DNA extractions. With the final protocol detailed in Figure 3, we were able to produce all four filter fractions in as little as three hours for water with relatively low biomass, though this could be greatly reduced if processing smaller volumes. Filtrations could extend up to 10 hours in some cases when the collection volume was high (up to 60 liters) and/or the sampled substrate contained higher levels of biomass, such as during an algal bloom. Compared to similar methods used for the isolation of picocyanobacteria and cyanophages, such as flow cytometry, eDNA, density-gradient centrifugations, culturing, dilution to extinction, and titer plating, this protocol can be completed the same day as sampling with relatively inexpensive equipment, moderate processing volumes, and can be easily adapted for varying target organisms based on size. Additionally, this protocol can be used with sufficient preparation in-situ as well as in a laboratory setting for filtration processing, and produced multiple filter fraction samples that reduce sequencing noise in downstream applications. The 30 primary sampling events were ultimately segregated by filtration and technical replicates into 294 filter subsamples and 424 DNA subsamples, with filters cut and used for replicate DNA extractions.

As expected, we found that water temperatures increased at the entrance of Sequim Bay through summer months, peaking from late July to early August 2024, with lows in February 2025 (Figure 4). The lowest recorded water temperature of the 26 sampling events from PNNL Sequim was 7.7 °C on 2 January, 2024, and the highest recorded temperature was 14.7 °C on 17 July, 2024. The sampling event at Cline Spit had the highest recorded temperature of 19.8 °C on 7 August, 2024. The highest practical salinity value derived from the MCRL dock data was 32.1 on 18 August, 2024 and the lowest was 30.3 on 7 April, 2025. Surface sampling events occurred across a range of tidal events, ranging from -0.2 to 2.2 meters in water depth. Of all 30 sampling events, two were sampled during a slack low tide, three during a slack high tide, 17 during a flood tide, seven during an ebb tide, and one event from the Sequim outdoor mesocosm tank did not have an assigned tidal stage.

DNA yields for various method and environmental context comparisons were statistically analyzed by using a Tukey’s post hoc test (p < 0.05) after performing an analysis of variance (ANOVA). Due to the highly limited sample sizes in some comparisons (as little as two in each group), these results can only be considered preliminary. Conclusions should be validated with further studies and greater sample sizes. The results for low sample size comparisons must be interpreted with caution, acknowledging that the impact of outliers or deviations may be amplified, and group variance and differences may not be accurately depicted.

In an exploratory controlled comparison of DNA yields, the Blood & Tissue kit outperformed the PowerWater and PowerSoil kits, yielding approximately 3900, 690, and 2100 ng DNA, respectively (Figure 5). For this controlled comparison, two pico filters (90mm PES) from a single sampling event were cut into thirds, and DNA was extracted from each third using a different kit (n = two extractions per kit). Samples labeled “ab” in Figure 5 denote an overlap between groups a and b, indicating that DNA yields in these samples are not significantly different from either group a or group b, while groups a and b are significantly different from each other, as determined by a Tukey’s post hoc test (p < 0.05) after performing an analysis of variance (ANOVA). However, a ‘kitchen sink’ comparison (a mix of conditions and replicates, controlling only for filter pore size and location) of each filter fraction across all DNA subsamples collected at PNNL-Sequim (macro, micro, pico, and viral) that also included the Zymo Plant/Seed kit (Figure 6) resulted in highly variable yields under each condition such that there were no significant differences in the micro, pico, and viral fractions, and overlapping significance groups in the macro fraction. Yields in the three direct filter fractions (macro, micro, and pico) were on average relatively high (~500–1000 ng or more), with a large proportion having acceptable yields for downstream sequencing (50 ng or greater; Figure 6). This supports our approach of running filters to saturation to obtain sufficient material for DNA extraction. Perhaps unsurprisingly, the viral fractions had generally lower DNA yields, averaging 450 ng compared to an average of 1210, 1272, and 924 ng for the pico, micro, and macro fractions, respectively. However, viral DNA yields were often still high enough for sequencing.

A similar kitchen sink comparison (controlling only for filter pore size and location) of MCE versus PES membranes also had high variation in yield. While the MCE membranes yielded slightly more DNA than the PES membranes for the pico fraction (1260 ng average for pico MCE and 1143 ng average for pico PES), the difference was insubstantial compared to the variation among samples (Figure 7a). In a more controlled, exploratory comparison of DNA yields from samples, the MCE membranes in the pico fraction were significantly greater than the PES membranes (p < 0.05, 1818 and 1673 ng DNA, respectively), but the differences were small relative to the DNA yield (Figure 7b). Differences between PES and MCE in the viral fraction were also relatively small (324 ng difference in broad comparison, and 54 ng difference in focused comparison) and either were not significantly different or had overlapping significance groups based on a Tukey’s post hoc test (p < 0.05) after performing an ANOVA. Losses of yield when extracting DNA from filters after a long freezer storage period (over 100 days) were small and not significant (Figure 8).

In-line serial filtration of the macro and micro fractions onto nylon and MCE membranes, respectively, using the Millipore Sigma polypropylene in-line filter holders (see Methods and Supplementary Data Tables) was observed to be the fastest method to remove larger-sized organisms, and facilitated large-volume processing in a laboratory while monitoring and changing out 20 µm or 5 µm filters if flow lessened as they became clogged. We note that for this reason the replicates for the 5 µm filters are taken sequentially from the same sample flow through a single 20 µm filter, and may differ based on differential coating, clogging, or disruption of the 20 µm filter over that time series (see relative read counts section below). In contrast, the single-flask sterile PES filters were observed to be more time-efficient for pico and viral fractions than in-line filtration on an MCE membrane, reducing processing time per liter from over an hour to a matter of minutes in some cases. The larger (90 mm) filter size also allowed collection of quality technical replicates of DNA extraction and sequencing by cutting PES filters into equal fractions prior to freezing. A list of all yield data for samples used in these analyses can be found in Supplementary Data 3.

Read mapping results for four sites (PNNL-Sequim, John Wayne Marina, Cline Spit, and Discovery Bay) show that amoeba read counts were dramatically reduced between 3- and 30-fold after passage through 5 µm to 0.22 µm filters (Figure 9), with the greatest fold reduction from the highest count 5 µm filter (John Wayne Marina) and the least fold reduction from the lowest count 5 µm filter (Cline Spit). In the Cline Spit sequences, there were only 287 reads from the 5 µm filter that matched amoeba sequences. Although this indicates that large-bodied cells are effectively screened out, relative Prochlorococcus reads were not substantially increased or decreased in any of the 0.22 µm filter samples, while relative Synechococcus reads were substantially reduced in some 0.22 µm filter samples by up to a third (Figure 9; Supplementary Data 1). Prochlorococcus normalized read counts were similar across sites and filters (range 718–1300 per 10 M reads), indicating that substantial numbers are trapped by the larger filters and possibly that the differential filtering of larger cells is not enough to increase the relative enrichment of small cells on the subsequent filter; other possibilities are discussed below. The coefficient of variation among Prochlorococcus mean-normalized read counts for the Figure 9 technical replicates (replicates of all steps from filtration onwards) on the 0.22 µm filters was also small at 6.8% (7.5%, 8.8%, and 5.0% for data in Figures 9a-c, respectively). The Synechococcus read counts were also similar across samples and showed slightly higher reads than those of Prochlorococcus with the exception of Discovery Bay, in which the Synechococcus read counts were approximately 10-fold higher than at Cline Spit (Figures 9c, d). The amoeba reads from the 5 µm filters at Cline Spit and Discovery Bay, were about 50%-100% greater than Synechococcus reads at PNNL-Sequim and John Wayne Marina. The number of reads mapping to amoeba and Synechococcus fluctuate much more than reads mapping to Prochlorococcus, but they are not in sync with each other across sites. It is also notable that the 24-mer nonpareil sequence diversity measures across all sites (Figure 9e) were consistently higher and more variable in the 5 µm filter fractions (range ~20.5-22.5) than the 0.22 µm filter fractions (range ~19.5-20). The 5 µm filters are thus at levels that are comparable to marine and soil samples (21-25) or freshwater and sandy soils (20-22) (Rodriguez-R et al., 2018), while the 0.22 µm diversity is comparatively reduced.

Synechococcus read counts were remarkably similar at ebb and flood tides for both the 5 µm and 0.22 µm filters (Figure 10). In contrast, amoeba reads were about 3x less frequent at flood versus ebb tide on the 5 µm filter and greatly reduced in frequency on the 0.22 µm filters (Figure 10). Prochlorococcus reads were also reduced in flood versus ebb tides and were reduced by 50%-60% on the 0.22 µm filter compared to the 5 µm filter. We note that the flood tide samples were extracted considerably later than the ebb stage samples, but time stored in freezer prior to extraction was not expected to change relative frequencies based on our exploratory analysis showing no significant impact of freeze-time on yield (Figure 8).

Relative frequencies of microbial taxa on the 0.22 µm filters over time at PNNL-Sequim, which lasted from mid-April to late October 2024, are shown in Figure 11. These comparisons unavoidably included variation in tidal stage, and also in extraction method and time between collection and extraction, and are across a single year (see Figure 11 legend), so must be interpreted in that context. As in previous figures, the amoeba reads are generally low, indicating relative exclusion by the 5 µm filter. Amoeba and Prochlorococcus reads had high relative frequencies in early spring and late fall, indicating that they were relatively more frequent than other organisms in the filter sample). In contrast, Synechococcus counts were higher from April to October, with the exception of counts in July. Counts on July 17, 2024 were clustered by processing conditions for all three taxon groups, indicating that the combination of kit and time in freezer prior to extraction can affect results at nearly the same level of variation as change over time (PSK kit and short time in freezer had much higher counts). September 9, 2024 amoeba samples were strongly clustered by condition, with high relative counts for all samples with a short freezer storage duration and extracted with the PWK kit (the Prochlorococcus samples were also clustered in the same fashion). Across all three plots, DNA samples extracted within a few days of filtration had consistently clustered higher counts than samples extracted greater than 50 days after filtration (stored frozen and thawed). This suggests that extended freezer time may contribute to decreased identifiable sequence reads, yet our focused comparisons to assess the impact of freezer storage duration on DNA yield found no significant effect (Figure 8).

The extraction of high-quality DNA from picocyanobacteria remains challenging, likely due to relatively low biomass collected in oligotrophic marine water, and potentially due to environmental contaminants, cell lysis inefficiencies due to cyanobacterial extracellular polymeric substances, and polysaccharide production (Wagner et al., 2024). Additionally, the consistently low biomass found at sampling sites—excluding Cline Spit, which had much higher biological material than the other sampling sites—necessitated high filtration volumes and extensive processing time.

The low sample size for many compared groups was a challenge to use for analyses and must be interpreted with caution. While ANOVA requires a minimum of one sample above the number of groups and can technically be performed under these conditions, it is not considered a sufficient sample size for statistical conclusions. We recognize that these results are only exploratory analyses to guide further research with larger sample sizes.

For our preliminary estimates, we used limited high-level estimates of sequence composition, and small numbers of technical replicates from identical samples (same time, location, and water sample with same filtration protocol). We also had few replicates of biological/environmental effects such as tidal state, time of year, or sampling location. Due to continuous sampling methodology adaptations throughout year one, there was high variation in sample production, which significantly decreased the number of quality technical replicates for focused comparisons. Especially for the non-controlled comparisons, results should be validated with larger sample sets and using consistent methods; however, the results presented here provide a valuable understanding of the impact that methodological variation can have.