Validating metabarcoding-based biodiversity assessments with multi-species occupancy models: A case study using coastal marine eDNA

k is present in the metacommunity. Here, the random variable zik describes the occurrence (zik = 1) or non-occurrence (zik = 0) of species k at site i having parameter ψik, which describes the occu- pancy probability for species present in the metacommunity. zik � BernoulliðcikakÞ Species may be present in the metacommunity but be unobserved during surveys due to imperfect capture and detection. Summing ak provides an estimate of species richness for the metacommu- nity, including observed and unobserved species. Case study: Conception Bay, Newfoundland Sample collection, processing and sequencing. Triplicate 250 mL water samples were collected from coastal surface water at eight sites along two transects in Conception Bay, New- foundland and Labrador, Canada, on October 13–14, 2017 (see [5] for sampling details). No permits were required to collect samples since there are no regulations on collecting seawater. Water samples were filtered using 0.22 μm PVDF Sterivex filters (MilliporeSigma) and DNA was extracted from filter membranes using the DNeasy PowerWater Kit (Qiagen). Five target amplicons in the cytochrome c oxidase I (COI) region were amplified by PCR from each sam- ple. Table 1 details the primer sets used to target these amplicons. Three PCR replicates were performed for each amplicon from each sample and then PCR replicates were pooled into a single PCR cleanup for each of the five amplicons with the QIAquick 96 PCR purification kit (Qiagen). Amplicons were then indexed using unique dual Nextera indexes (IDT). All ampli-

cons were pooled into one library to normalize DNA concentration and the library was sequenced with a 300-cycle S4 kit on the NovaSeq 6000 following the NovaSeq XP workflow. Raw sequence reads are available in NCBI’s sequence read archive under accession number PRJNA574050. Primers were trimmed from sequences and then DADA2 v1.8.015 [32] was used for quality filtering, joining paired end reads and denoising to produce exact sequence variants (ESVs). Taxonomy was assigned using NCBI’s blastn tool v2.6.026 [33] to compare ESV sequences against the nt database. See [5] for detailed sampling, sequencing, and bioinfor- matic methodology. Occupancy model implementation. Under the occupancy modelling framework described above, each collection site along each transect in Conception Bay was considered a different site in the occupancy model. Replicate bottles collected at a site were considered sam- ples. Each of the five amplicons sequenced from each bottle was considered a technical repli- cate. While we conducted three replicate PCRs for each amplicon, PCR products were pooled for each amplicon prior to sequencing so we did not include PCR replicates for each amplicon separately in our models. However, PCR replicates can easily be accommodated in a multi- scale, multi-species occupancy modelling framework, such as the model described here. Table 1. Primer pairs used to amplify five target amplicons in the COI region of the mitochondrial genome from water samples collected in Conception Bay, New- foundland, Canada. Marker Target Length (bp)

Forward Primer Reverse Primer Reference Fishe (Mini_SH-E) 226 5'-CACGACGTTGTAAAACGACACYAAICAYAAAGAYATIGGCAC-3' 5'-GGATAACAATTTCACACAGGCTTATRTTRTTTATICGIGGRAAIGC-3' [34] Fishc (Mini_SH-C) 127 5'-CACGACGTTGTAAAACGACACYAAICAYAAAGAYATIGGCAC-3' 5'-GGATAACAATTTCACACAGGGAARATCATAATGAAGGCATGIGC-3' [34] F230 235 5'-GGTCAACAAATCATAAAGATATTGG-3' 5'-CTTATRTTRTTTATNCGNGGRAANGC-3' [35] Leray 330 5'-GGWACWGGWTGAACWGTWTAYCCYCC-3' 5'-TAAACTTCAGGGTGACCAAAAAATCA-3' [36] BR5 310 5'-CCIGAYATRGCITTYCCICG-3' 5'-GTRATIGCICCIGCIARIACIGG-3' [37] https://doi.org/10.1371/journal.pone.0224119.t001 PLOS ONE Multi-species occupancy models for metabarcoding PLOS ONE | https://doi.org/10.1371/journal.pone.0224119 March 19, 2020 5 / 17