IS IMBMG Metabarcoding & Metagenomics Research Article 8 Metabarcoding and Metagenomics 6: 305-317 DOI 10.3897/mbmg.6.85199 National eDNA-based monitoring of Batrachochytrium dendrobatidis and amphibian species in Norway Omneya Ahmed Osman’, Johan Andersson’, Pedro M. Martin-Sanchez**, Alexander Eiler’ 1 eDNA solutions AB, Bjorkasgatan 16, SE-43131 Molndal, Gothenburg, Sweden 2 Watercircle, Rantmastaregatan 23c, SE-416 58 Gothenburg, Sweden 3 Department of Biosciences, University of Oslo, P.O. Box 1066 Blindern, 0316 Oslo, Norway 4 Instituto de Recursos Naturales y Agrobiologia (IRNAS-CSIC), Avenida Reina Mercedes 10, 41012 Seville, Spain Corresponding author: Omneya Ahmed Osman (omneya@ednasolutions.se) Academic editor: Florian Leese | Received 13 April 2022 | Accepted 17 August 2022 | Published 7 October 2022 Abstract Freshwaters represent the most threatened environments with regard to biodiversity loss and, therefore, there is a need for national monitoring programs to effectively document species distribution and evaluate potential risks for vulnerable species. The monitor- ing of species for effective management practices 1s, however, challenged by insufficient data acquisition when using traditional methods. Here we present the application of environmental DNA (eDNA) metabarcoding of amphibians in combination with quan- titative PCR (qPCR) assays for an invasive pathogenic chytrid species (Batrachochytrium dendrobatidis -Bd), a potential threat to endemic and endangered amphibian species. Statistical comparison of amphibian species detection using either traditional or eDNA-based approaches showed weak correspondence. By tracking the distribution of Bd over three years, we concluded that the risk for amphibian extinction is low since Bd was only detected at five sites where multiple amphibians were present over the sampled years. Our results show that DNA-based detection can be used for simultaneous monitoring of amphibian diversity and the presence of amphibian pathogens at the national level in order to assess potential species extinction risks and establish effective management practices. As such our study represents suggestions for a national monitoring program based on eDNA. Key Words environmental monitoring, invasive species, Metabarcoding, qPCR, risk assessment Introduction Freshwater environments are threatened by biodiversity loss, and there are over 100 documented cases of extinc- tion in such environments during the 2000s (Dudgeon et al. 2006; Tickner et al. 2020). In the case of amphibi- ans, 30% of all freshwater amphibian species are threat- ened by extinction (Duefias et al. 2021). The spread of invasive species because of human activity has been proposed as one of the greatest threats to indigenous species (Beaury et al. 2020). One such invasive species is the chytrid fungus Batrachochytrium dendrobatidis (Bd) which causes the disease chytridiomycosis in am- phibians. Chytridiomycosis has led to large reductions in the populations of amphibian species, and in some cases extinction across the world (Laurance et al. 1996; Berger et al. 1998; Longcore et al. 1999; Mendelson et al. 2006). Tickner et al. (2020) recently presented an emergency recovery plan for freshwater biodiversity, which includes protecting and restoring critical habitats, and preventing the introduction and spread of non-native species. In addi- tion, international treaties such as the global biodiversity framework under the Convention on Biological Diversity and the EU’s biodiversity strategy for 2030 have been ini- tiated. The latter has been established to facilitate an EU- wide restoration plan and network of protected areas, as well as introducing measures to enable necessary transfor- mative changes and to tackle the global biodiversity chal- lenge. Here the monitoring of programs of biodiversity Copyright Omneya A. Osman et al. This is an open access article distributed under the terms of the Creative Commons Attribution > PENSUFT. License (CC BY 4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. 306 Omneya A. Osman et al.: eDNA monitoring of pathogen fungus Batrachochytrium dendrobatidis in Norway is a prerequisite to design, implement and validate these international conservation and restoration efforts. To be successful, biodiversity monitoring must be de- signed in such a way as to minimize the main sources of error (Skalski and Robson 1992; Thompson et al. 1998). A common source of error comes from the fact that all survey methods fail to detect all individual species in an ecosystem (Mathieu et al. 2020). A second source of error arises from the difficulty in efficiently investigating large areas, meaning that conclusions must be based on a low number of samples taken from a few locations. This is compounded by the fact that many collection strategies are based on subjective assessments of how represen- tative certain sampling locations are, or how easily ac- cessible sites are. Thirdly, environmental managers face major problems in identifying and counting species, since organisms can be difficult to detect or to distinguish from one another. Hence, there is a need for novel methods to survey biodiversity at the national level. Technological advances in molecular biology over the past decades now allow us to minimize some sources of error. Environmental DNA (eDNA) analysis 1s a rapidly evolving methodology used for studies of current and past biodiversity (Valentini et al. 2009; Taberlet et al. 2012). eDNA has broad applications in the analysis of biodiver- sity in microbes, plants, and animals (Eiler and Bertilsson 2004; Zinger et al. 2012; Valentini et al. 2016), analysis of diet (Deagle et al. 2005; Pompanon et al. 2012), recon- struction of past biodiversity or environmental changes (Jorgensen 2012; Parducci et al. 2013; Langenheder et al. 2016), and environmental monitoring (; Eiler et al. 2013). In principle, biodiversity across the entire tree of life present in a particular system can be assessed by DNA metabarcoding (Stat et al. 2017). Substantial advantages of eDNA-based methods are higher cost- and time- ef- fectiveness compared to many traditional survey methods (Evans et al. 2017), their noninvasive nature (Cristescu and Hebert 2018), and high specificity and sensitivity (Wilcox et al. 2013). Testing for the presence and concentration of micro- bial pathogens such as Bd using eDNA methods is ap- pealing because it relies on noninvasive sampling, and free-living stages persistent over several weeks can be detected (Brannelly et al. 2020). Similarly, amphibian diversity can be assessed without the need to find and capture animals in the environment, thus avoiding harm- ing animals and introducing sampling biases (Valentini et al. 2016). Despite the high sensitivity of eDNA methods, species detection by eDNA approaches is not free from sources of error similar to those found using conventional methods. In fact, incomplete detection is inevitable when assessing the presence/absence and number of species, regardless of the methods used (MacKenzie et al. 2006). This can be partially overcome by the high sensitivity and capacity of eDNA methods that allow more samples to be processed, leading to a higher likelihood of detection, than when conventional methods are used. In addition, https://mbmg.pensoft.net adaptive sampling can increase the likelihood of detec- tion, for example by sampling during the optimal season and adapting spatial sampling strategies to the target or- ganisms (Buxton et al. 2018). In this study, we aimed to assess the distribution of the chytrid fungus Bd from water samples by employing eDNA methods. Bd is regarded as a generalist fungal pathogen, presently known to occur in Sweden (Rosquist 2020) and Denmark (Scalera et al. 2008), but little is known about the spread and virulence of Bd in Norway (Taugbel et al. 2021). Two quantitative PCR (qPCR) assays and a metabarcoding approach were used to detect Bd and amphibian species, respectively, throughout southern Norway. The specific objectives of this study were to (1) assess the spread of Bd in the southeastern part of Norway over several years, and (11) explore potential risks for amphibian hosts using co-oc- currence analysis and by reviewing the literature. Lastly, (ii1) we also intend to provide methodological suggestions for eDNA based national monitoring programs. Materials and methods Sampling design Sampling sites located in Southern Norway (regions of Viken, Oslo, Vestfold/Telemark, Agder and Innlandet) were chosen randomly from sites with observations of am- phibian species as reported in the Norwegian species da- tabase “Artsdatabanken”. Resampling for each month was done on the complete database, meaning that sites were sampled multiple times. As the density of sites with spe- cies observations is related to human presence, our sam- pling focused on anthropogenic impacted and potential Bd positive sites. Introduction of Bd is proposed to be related to human activities (Nielsen et al. 2019). During 2019, 110 bodies of water, including forest-, agricultural- and urban ponds, were sampled. To test for the optimal sampling sea- son, we collected and analyzed samples from three differ- ent months in 2019: May or early June (37 samples), July (36 samples) and September (39 samples). These sampling times correspond to life history stages such as spawning, post-spawning and juveniles, respectively (See Suppl. ma- terial 1: Table S1). To compare the detection probabil- ities between the traditional and eDNA-based sampling, 91 bodies including forest, agricultural and urban ponds, and 10 swabs from various amphibian specimens were collected over a single season in 2020. Swabs were taken from frog and salamander skins at Hasseldalen (1 swab), Vermyr (3 swabs), Ringveien 52 (3 swabs) sites, while one swab from captured tadpoles was collected from each of Butjenna, Oyenkilen, and Sodra Skauen Gard water- bodies (See Suppl. material 2: Table S2). Swabs and tad- poles were preserved in 70% ethanol on site and stored at 4 °C until further analysis. Details about the survey design including geographical and temporal distribution of the samples are given in Suppl. material 3: Fig. S1. Metabarcoding and Metagenomics 6: e85199 Water samples were collected with a beaker attached to an expandable sampling pole (a device that can reach 4 meters). At least 5 samples (of 20-100 ml) per site were pooled together. We used cartridge filters (Sterivex filter units 0.22 um pore diameter, Millipore, Merk, Darmstadt, Germany) and sterile 50-ml syringes to filter between 75 (minimum) and 1440 ml (maximum) of pooled water samples from each site (median volume of 480 ml). Sam- pling was performed by the same two persons throughout the study. A syringe - cartridge filter setup allowed for reproducible sampling without the need for heavy equip- ment and access to electricity. Filters were immediately frozen in dry shippers cooled with liquid nitrogen and then stored at -80 °C in the laboratory until further anal- ysis. Chemical and physical parameters such as tempera- ture, conductivity, pH, and turbidity were also measured on-site (See Suppl. materials 1 and 2). DNA extraction and Bd qPCRs All lab benches and equipment were cleaned by using 5% sodium hypochlorite and 70% ethanol in the field and lab- oratory. This included the three laboratories used for DNA extraction, prePCR or postPCR. Here pipettes and con- sumables were subjected to UV exposure for 30 min prior to usage. For water samples, DNA was extracted from the Sterivex filters using a Qiagen DNeasy PowerWater Ste- rivex Kit (Qiagen, Germany) including one unused filter as a negative control in one of the DNA extraction runs. The collected swabs were removed from ethanol and left to dry in clean Eppendorf tubes for 2 hours at 50 °C to evaporate residual ethanol. The DNA extraction was sub- sequently performed using the Qiagen Blood & Tissue DNA extraction kit. The quantity and quality (280/260 ra- tio) of the extracted DNA were assessed using NanoDrop. Bd qPCR assays: We used two TaqMan assays am- plifying different regions of the ribosomal RNA op- eron (ITS1 — 5.8S rRNA gene — ITS2). First the Boyle TaqMan assay (Boyle et al. 2004), targeting the Internal Transcribed Spacer 1 (ITS1) region where qPCR mixture contained 10 wL of 2x BioRad SsoAdvanced Universal Probes Supermix, 0.8 uM of the forward primer ITS1- 3 Chytr (5’- CCTTGATATAATACAGTGTGCCATAT- GTC-3') and the reverse primer 5.8S Chytr (5'- AG- CCAAGAGATCCGTTGTCAAA-3’') and 0.24 uM of Chytr MGB2 (FAM-5S' TTCGGGACGACCC-3'-NFQ- MGB), | uL of internal positive contro (IPC, if not co-ex- tracted with the sample), 1 uwL IPC primer and probe mixture, 5 wL of environmental DNA sample and RNase/ DNase free water to reach 25 wL as final volume. The second assay was based on the commercially available product Batrachochytrium dendrobatidis 5.8S rRNA Ad- vanced Genesig Kit (Primer design Ltd, UK) which was carried out according to the manufacturer’s instructions. Both qPCR assays were conducted in a BioRad CFX96 thermocycler (USA) using the following program: 2 min at 95 °C, followed by 50 cycles of 10s at 95 °C and 1 min 307 at 60 °C. At least 8 serial standard dilutions were test- ed and validated before choosing five standard dilutions (1000, 100, 10, 1 and 0.1) to perform qPCR runs with en- vironmental samples in triplicates. Bd -positive controls were added in the post-PCR room while the qPCR master mix was prepared in pre-PCR room. For obtaining the standard curves needed to quantify Bd in environmental samples with the Boyle TaqMan assay, we used a DNA extract from Bd spores provided by Professor Anssi Laurila’s research group at Uppsala University. This extract had an expected concentration of 100 genomes (or genome equivalent) per uL. The number of rRNA operons in each spore was not known, therefore concentration was given as genome equivalent. Samples collected in 2019 were analyzed with only the commercial assay while samples collected in 2020 were analyzed using both assays in triplicates including two negative controls, in 96 well plates. Samples were de- clared Bd positive if two or more of the three qPCR runs generated a positive Bd signal within the range of serial dilutions of the standard. If only one of the three qPCR runs was positive, an additional fourth qPCR run was conducted. Sanger sequencing was performed to confirm the amplification of the target fragment for samples with Ct values > 39.5. To compare our results with known oc- currences of Bd, we included data from a previous Nor- wegian survey in 2017 (Taugbel et al. 2021). In addition, to further assess the sensitivity of this as- say and estimate the target copy numbers (ITS1 copies), we constructed an extra standard curve using a synthet- ic gBlocks double-strain ITS1 fragment (“Bd_ 26-271”; this name refers to the positions in the reference sequence NR_ 119535; 246 bp) based on eight decimal serial dilu- tions from 0.39 to 3.9 x10° ITS1 copies per reaction. The limit of detection (LOD) was estimated using a discrete approach and defined by the highest dilution where at least a single replicate of the triplicated standard dilution of 5 serial concentrations showed amplification in the respec- tive assay consistently throughout all qPCR experiments with at least 95% confidence (modified from Kylmus et al. 2019). To define the limit of quantification (LOQ), we identified the highest dilution where a R?> 0.98 and an ef- ficiency between 85 and 115% was obtained when fitting a linear regression to log-transformed Bd quantity and Cq values of the dilution series for each individual experi- mental (qPCR) run. Internal positive control: In order to determine the DNA extraction efficiency and the presence of inhibi- tors, a separate internal extraction control DNA, prim- ers and probe mixture with different target region, IPC, was provided with the Batrachochytrium dendrobatidis 5.8S rRNA Advanced Genesig Primer design kit (Unit- ed Kingdom). In 2019, the batch of samples collected in September were co-extracted with IPC by mixing 4 uL of IPC with Qiagen lysis buffer (as described in the manu- facturer’s instructions) and tested in a qPCR run. In 2020, 30 randomly selected DNA samples were co-extracted https://mbmg.pensoft.net 308 Omneya A. Osman et al.: eDNA monitoring of pathogen fungus Batrachochytrium dendrobatidis in Norway with IPC while the remaining extracted DNA samples were spiked with the IPC. To spike DNA with the internal control, the internal control was diluted 1:20 and 1 uwL was added to a subset of samples which did not under- go internal control co-extraction, prior to the qPCR. The internal control is detected through the VIC fluorophore channel. We assumed a 100% extraction efficiency if Cq value of 27 was observed with Cq values of 27+2 within the normal range. Samples with Cq values of 30 or higher were considered to have PCR inhibitors. Amphibian mock community Mock communities are now widely used in metabarcod- ing as a positive control and to validate analysis proce- dures from library preparation to bioinformatics analyses. A mock community of five amphibian species was pre- pared from samples provided by the lab of Anssi Laurila (Uppsala University) and Nordens Ark, including DNA extracts of Rana arvalis and Bufo bufo, and swabs pro- vided by Anssi Lab (preserved in 95% ethanol) of Bufotes variabilis, Epidalea calamita, and Pelophylax lessonae. DNA extraction from the swabs was conducted as de- scribed above. Quantification of the DNA concentration was performed using the PicoGreen assay (ThermoFish- er, US). Approximately 10 ng of the extracted DNA from each species were mixed in equal amounts and diluted to prepare five dilutions: 1.0, 0.5, 0.1, 0.01 and 0.001 ng/uL. Amphibian DNA metabarcoding-Library preparation Paired-end sequencing on the Illumina MiSeq plat- form was based on two steps of PCR. The first step amplified the mitochondrial 12S rRNA gene of am- phibian species using the following primers: forward batra-F 5'-ACACTCTTTCCCTACACGACGCTCTTC- CGATCTNNNNNNACACCGCCCGTCACCCT-3' and reverse batra-R 5'-AGACGTGTGCTCTTCCGATCT- NNNNNNGTAYACTTACCATGTTACGACTT-3'. Human blocking primer: 5'-TCACCCTCCTCAAG- TATACTTCAAAGGCA-SPC3I-3' was used to bind to human DNA and prevent its amplification (Valentini et al. 2016). This first PCR was performed in a final volume of 25 wL containing 5 uL of 5x Q5 reaction buffer, 0.2 uM “batra” primers, 0.2 mM dNTPs, 4 uM human blocking primer and 0.02 U/uL Q5 High-Fidelity DNA Polymerase (New England Biolabs) as well as 1 wL of positive control (mock community) or 5 uL of eDNA extract as a tem- plate. Amplification was performed in a Veritipro Ther- mal Cycler PCR Machine (Applied Biosystems, USA) under the following conditions: 30 s at 98 °C, followed by 35 cycles of 30 s at 98 °C, 30 s at 57 °C and 1 min at 72 °C and a final extension step at 72 °C for 7 min (Val- entini et al. 2016). Each sample was amplified in triplicate and then pooled before the purification step. We also used two negative PCR products to monitor contamination and the mock community as a positive control of the ampli- fication in each PCR run. 1% agarose gel electrophoresis https://mbmg.pensoft.net was run to check the successful amplification of the target band (170 bp) before purifying the PCR products using AMPure XP magnetic beads (Beckman Coulter, US). Twenty forward and reverse index primers developed based on Sinclair et al. (2015) were used for Illumina sequencing. This second PCR was performed in a final volume of 20 ul containing 4 ul of 5x Q5 reaction buf- fer, 0.2 mM dNTPs, a combination of the indexed for- ward and reverse primers (0.25 uM each), 0.02 U/pl Q5 High-Fidelity DNA Polymerase (New England Biolabs) and 2 uwL of purified first step PCR product. Amplifica- tion conditions were 30 s at 98 °C, followed by 15 cycles of 10 s at 98 °C, 30 s at 66 °C, and 30 sec at 72 °C, then a final extension step at 72 °C for 2 min. The in- dexed target band was checked using 1% agarose gel electrophoresis. The second PCR product was purified using the same reagent as above, and quantified using the PicoGreen assay. The detailed DNA metabarcoding protocol for amphibians was published on protocols.io (dx.doi.org/10.17504/protocols.10.973h9qn). The first and second PCR mixtures were prepared in pre-PCR lab- oratory while the first purified PCR product was added in the post-PCR laboratory. Equal amounts of DNA from each amplified sample were mixed to prepare the sequencing library. In addition to sequencing the five dilutions of the mock community, two or three environmental samples were spiked with one of the mock community species which was processed in the same sequence pool to confirm the presence and ab- sence of the identified species. Furthermore, one random negative PCR product (entire sample) was included in each Illumina sequence run performed in 2019 and 2020. To detect all amphibian species in positive swab samples, the three swabs collected from positive Bd sites were se- quenced while one swab was only sequenced from the other sites. The final pooled samples were visualized by gel electrophoresis to ensure that no extra unwanted frag- ments existed. The amplicon library was then sequenced on an Illumina MiSeq machine using the MiSeq v2 PE 150 Micro protocol, generating approximately 20 million paired-end reads of 150 bp length and exact matches to batra primers. Raw sequence data are available at NCBI through the SRA accession no. PRJNA821076 for 2019 and PRJNA822096 for 2020. Reference database construction Amphibian 12S rRNA gene sequences were obtained by searching NCBI and BOLD using species names of the Norwegian amphibian species and *12S’ as search criteria. Next, a homology search using Basic-Local-Alignment- Search-Tool (BLAST) was used to extract rRNA gene se- quences of high (90%) homology. This two-step process led to a 12S rRNA gene sequence database of mostly am- phibian sequences. To clean up the database, sequences were aligned and the “batra” primers were matched. This resulted in a database of 7 sequences including all species that have been reported in the wild in Norway. Metabarcoding and Metagenomics 6: e85199 Observational data mining Traditional observational data of amphibian diversity was obtained from the Norwegian public biodiversity databank “Artsdatabanken” (https://artskart.artsdatabanken.no/). Data for amphibians was downloaded on 2020-05-11 using “Am- phibia” as a query resulting in 4086 species observations in the region overlapping with our samples. Species observa- tions from the last 20 years were kept and then grouped if GPS coordinates overlapped as defined by resolution to the second decimal corresponding to a place up to 1.1 km”. Sequence analysis Raw sequences were first processed with cutadapt v1.18 to remove PCR primers, and then analyzed with the R pack- age dada2 v1.14.1 (Callahan et al. 2016) for denoising and sequence-pair assembly. After manual inspection of quality score plots, forward and reverse reads of the amplicon se- quencing run were trimmed to 50 bp length. This assured an almost complete overlap of the forward and reverse reads as the amplicons without primers were 51—54 bp long. Additional quality-filtering steps removed any sequences with unassigned base pairs and reads with a single Phred score below 20. After reads were dereplicated, forward and reverse error models were created in dada2 with a subset of the sequences (10° nucleotides). Amphibian 12S rRNA gene amplicons were assembled by merging the read pairs. Chi- meras were removed using ‘removeBimeraDenovo’ from the dada2 package, which resulted in the final Amplicon Se- quencing Variant (ASV) table for performing taxonomic as- signments. Taxonomy (to species level) was assigned using blastn and an in-house 12S rRNA gene database. Sequenc- es were assigned to specific species using cutoff criteria of identity > 96% (which corresponds to approximately two mismatches) and alignment length > 45 bp. Statistical analyses Heatmaps were constructed using the ggplot2 R pack- age, plotting the number of samples uniquely detected by either eDNA or traditional observation of species 309 deposited to the Norwegian species database “Artsdata- banken” (ABD method), as well as those samples where the species was detected using both methods. Correspon- dence of species observations between eDNA and those recorded in Artsdatabanken (using GPS data rounded to the second decimal) were evaluated by correlation analysis, fisher exact test and Pearson’s Chi-square test using R (version 3.6.2). In order to trace the outcome of the sequence analysis we used the sequenced amphibian mock community. Generalized linear (glm; base R) and additive (gam; “mgcv” package in R) models were used to explore de- tection rates of amphibian species in relation to time of sampling and pH. We also performed co-occurrence anal- ysis using the R package “cooccur’’. The latter did not re- veal any significant relationships, probably due to sparse Bd detection. Results Bd detection and comparison of two TaqMan assays We tracked the distribution of the invasive chytrid fungus Bd in the sampled water bodies. By using two TaqMan assays we were able to detect Bd in six water bodies and could confirm its presence in most cases by repeated sampling (Table 1). Still, our results indicate that spring is the season of highest detection probabil- ity for Bd, as samples taken during other seasons were all negative. We used two TaqMan assays to amplify different re- gions of the rRNA operon (ITS1 and 5.8S rRNA gene). qPCR efficiency of the Bd 5.8S rRNA Genesig Advanced Kit and the Boyle TaqMan assay was on average 98.82% and 99.34%, respectively. The limit of quantification (LOQ) of the assays depending on the run varied between 100 and 10 gene copies for the commercial Genesig kit and | and 0.1 genome equivalents for the Boyle TaqMan assay (Fig. 1; Suppl. material 3: Table S3). The limit of detection (LOD), as assessed from at least three indepen- dent runs, varied between 10, 1 and 0.1 gene copies for Table 1. List of positive Bd samples based on commercial and Boyle qPCR assays from 2019-2020. One of the three collected swabs from the site no. 38 were Bd positive using both assays. Data points from 2017 are taken from Taugbol et al. (2021), list of all samples with Bd results and amphibian species detections are represented in Suppl. materials 1 and 2. Vol: volume, Temp: temperature, Turb: Turbidity. ID Month Year Vol Temp “Turb 10 26-May 2017 600 ilo) 0 13 26-May 2017 180 15 3.61 6 26-May 2017 200 17 0.64 3 26-May 2017 190 15 31,83 1 26 May 2019 540 12.5 NA 3 26 -May 2019 190 IS 31.83 12 26-May 2019 340 14 8.63 69 01-June 2020 600 24 2.16 po 01-June 2020 360 23.4 4.47 74 01-June 2020 NA 24.6 0.81 Swab 38B 30-May 2020 540 29 6.33 10 25-May 2020 480 18 DN A7 30-May 2020 450 16 1.75 pH Type Locality Assay Be Forest Frogn Taugbol et al. 2021 6.71 Forest Nesodden Taugbel et al. 2021 7.0 Agricultural As Taugbgl et al. 2021 7.89 Agricultural Ostre Stokken Taugbel et al. 2021 NA Urban Varli Commercial 7.89 Forest QOstre Stokken Commercial FAS Urban Ringveien 52 Commercial Jeo) Agriculture QOstre Stokken Commercial & Boyle 7.09 Suburb Ringveien 52 Commercial & Boyle EF Suburb Odden Commercial & Boyle 7.8 Agriculture § Vermyr, Enningdalen Commercial & Boyle 5.96 Forest Butjenna Boyle 74 Agriculture Rokkevannet Boyle https://mbmg.pensoft.net 310 cq 29 y = - 3.6875x + 33.56 15 R®= 0.99 E=87% Omneya A. Osman et al.: eDNA monitoring of pathogen fungus Batrachochytrium dendrobatidis in Norway 45 40 X 35 - 30 25 ion kc y = - 3.4583x + 38.98 15 R? = 0,99 E=95% Log Bd Quantity (Genome Equivalents) Log Bd Quantity (ITS1 copies) Figure 1. Examples of standard curves used to determine the limits of detection (LOD) and quantification (LOQ) for the Boyle TaqMan assay using spore DNA (providing genome equivalents) (A) and the synthetic gBlocks fragment Bd_27-271 (providing ITS copy numbers) (B) as standards. Crosses represent standard dilutions above the LOD but below the LOQ while circles represent standard dilutions used for determining the corresponding standard curves and their R? and efficiency values (detailed on the panels). For the spore DNA curve (A), detailed statistics on multiple runs are provided in Suppl. material 3: Table S3. the commercial Genesig kit and 1, 0.1 and 0.01 genome equivalents for the Boyle TaqMan assay (see Suppl. mate- rial 3: Table S3). Thus it can be speculated that the Boyle TaqMan assay is slightly more sensitive than the Gene- sig kit assay. The higher sensitivity of the Boyle TaqMan assay was also confirmed using dilutions of the gBlocks fragment Bd_27-271 as standards, as the lowest concen- tration detected (in 2 of 3 triplicates) corresponds to 0.39 ITS1 copies (Fig. 1B), while the LOQ (based on a single qPCR experiment) seems to be two orders of magnitude higher (39 ITS1 copies). It is well known that rRNA oper- on copy number per genome (or spore) 1s highly variable between strains of the same fungal species, as reported by Longo et al. (2013) for Bd genomes which contain from 10 to 144 ITS1 copies. Both qPCR assays generated positive Bd signals from water samples of three locations: Ostre Stokken, Ringvei- en 52, and Odden, and from one swab from the location Vermyr (Table 1). Two additional locations, Butjenna and Rokkevannet, generated positive Bd signals when us- ing the Boyle TaqMan assay. However, as a high Cq value (as in the case of Butjenna and Rokkevannet) can repre- sent unspecific amplification and thus false positive Bd results, both PCR products were purified and Sanger-se- quenced for validation. The resulting sequences were used to query the NCBI nt-database using blastn, which resulted in 100% matches with sequences from Bd, and thus the presence of Bd was concluded. Internal positive control We used internal positive controls (IPC) consisting of a fully synthetic DNA, primers and a TaqMan probe, which are included in qPCR reactions or co-extracted with the DNA extract to detect PCR inhibition and false negative results due to inhibitory substances in the environmental samples (Hoorfar et al. 2004; Phillips 2004). In 2019 and 2020, five and six samples, respectively, showed delayed Cq values above 27+2 for the IPC. This variability in Cq results strongly suggests the presence of inhibitory com- https://mbmg.pensoft.net pounds in a minor fraction of samples, which needs to be taken into account when interpreting negative results (Suppl. material 3: Fig. S2). The IPC was detected in all co-extracted samples indicating the efficiency of the ex- traction protocol used. Three of these samples generated a normal IPC-qPCR signal after dilution, which is one way to overcome the effect of inhibitors. It should be mentioned that only half of these samples were success- fully sequenced and one amphibian species was identified in each of them, thus metabarcoding results from these samples are likely biased as well. By reporting IPC-qPCR signal, potential inhibition can be observed providing a good indicator of false-negative PCR results. A way to solve some false-negative detection might be by sample dilution or optimizing qPCR conditions (Kamal et al. 2017; Lance and Guan 2019) or the use of an inhibitor removal kit (McKee et al. 2015). Evaluation of amphibian metabarcoding Amphibian diversity was assessed using a metabarcoding approach based on amplification of the mitochondrial 12S rRNA gene of amphibian species and subsequent sequencing. Positive controls suchas the mock community including DNA from R. arvalis, B. bufo, B. variabilis, E. calamita, and P. lessonae showed amplification down to 10° ng of template DNA and resulted in sequences for most mock species (Fig. 2). The exception was B. variabilis that could not be detected when the targeted DNA was below 107 ng of template DNA. Most of the species were detected in low concentration indicating that the metabarcoding approach 1s still sensitive to very low amounts of genomic DNA. Next, we evaluated species discrimination of native Norwegian amphibian species by the DNA metabarcod- ing approach. The comparison of the available reference sequences by pairwise alignment revealed sufficient di- vergence amongst native Norwegian species for annota- tion at species level. The closest species were R. arvalis and R. temporaria with a minimum of two nucleotide Metabarcoding and Metagenomics 6: e85199 Sx = species Bufo bute Rana_arvalis © Epidalea_ecalamita © Bufotes_variabilis © Pelophylax_esculentus ix = Dilution 0% = 1000 x= or T T T — B00 25 0.50 ars 1.00 proportion of reads obl Species Bufo_bufo Rana_arvalis © Epidalea_calamita ©) Bufotes_variabilis © Pelophylax_esculentus — i= a d Dilution 1000 x = 0.00 0,25 o.50 075 1.00 proportion of reads Figure 2. Sequencing results from the dilution series of an amphibian mock community in 2019 (A) and 2020 (B). Both plots showing the five added species in the original, 5x and 10x diluted mock samples, while Bufotes variabilis disappeared in the 100x diluted sample. differences in the amplified 12S rRNA gene region while all other species differed in at least four positions. Spe- cies discrimination may, however, become an issue with this metabarcoding approach in areas of high amphibian diversity including multiple closely related species, such as in the tropics. Analyses of the 110 and 91 samples from 2019 and 2020, respectively, resulted in all of the samples being well amplified as reflected by the amplified amplicon size (approx. 50 bp, expected to the target region). No amplification was obtained from negative control filters, thus allowing us to exclude the possibility of cross-con- tamination among filters and confirming the clean DNA extraction and purification steps. Moreover, there was no amphibian species detected in any negative PCR control of the metabarcoding after the end-point (no resulting se- quences), even if a band appeared in gel electrophoresis due to the primer dimer formation, confirming a clean PCR setup. A stringent raw sequence data quality filtering resulted in a loss of approximately 60% of the reads from an av- erage of 15,736 raw paired end reads per sample (range from 830 to 87,543). From an average of 6,181 (range from 587 to 29,720) quality-filtered paired end reads per sample, denoising and merging resulted in a mean num- ber of 5,676 high quality sequences (range from 567 to 29,136). Amphibian sequences were only detected in 68 of 110 samples collected in 2019, and in 69 of 91 sam- ples from 2020. Unspecific amplification of non-targeted species could be observed in almost all amplified sam- ples. In 2019 and 2020 the mean percentages of sequenc- es per sample assigned to amphibian taxa were 9% (range 0 to 96%) and 30% (range 0 to 100%), respectively, cor- responding to 947 (range 0 to 26,327) and 1,862 (range 0 to 12,667) amphibian reads per sample. Closer inspection of the non-amphibian sequences showed various taxo- nomic assignments from bacteria to fish. This indicates unspecific amplification in the PCRs and emphasizes the need for a high sequencing depth to detect amphibian spe- cies when using the batra primers. Still, as we obtained on average more than 5,000 paired reads per sample, false negatives as a result of sequence library preparation could be minimized. Looking at the 2019 data, amphibian detection was sig- nificantly lower in September compared to May/June and July (estimate = -1.635, z-value = -3.710, p-value = 0.0021 as given by generalized linear models). This was corrob- orated by the total number of detected amphibian occur- rences among all sites (Table 2) and at six sites that were sampled during the three sampling occasions in 2019. In autumn, amphibians could not be detected at any of the six sites while in May/June and July there were multiple species detections. Thus it can be concluded that sampling during spawning and post-spawning results in higher like- lihood of detection when compared to autumn sampling. In a previous study, a comparison of the detection of the pool https://mbmg.pensoft.net 312 Omneya A. Osman et al.: eDNA monitoring of pathogen fungus Batrachochytrium dendrobatidis in Norway frog (P. lessonae) using both eDNA and traditional methods revealed that traditional methods gave a higher rate of ob- servation in June, whereas eDNA gave at least as many or more observations during other months (Eiler et al. 2018). Table 2. Comparison of species detection among the three sam- pling times (May, July and September) in 2019. Numbers rep- resent the number of sites where a species was detected, clearly showing that the number of amphibian species positive sites was lowest in September. A list of all samples, raw sequences, sampling volume and environmental parameters is presented in Suppl. materials 1 and 2. ND: not detected. Species May/June July September Rana arvalis 2 ND ND Rana temporaria 9 Z 3 Bufo bufo 10 3 1 Triturus cristatus 6 10 ND Lissotriton vulgaris 9 14 3 Amphibian species distribution in Norway Five amphibian species B. bufo, L. vulgaris, R. arvalis, R. temporaria, and T: cristatus are prevalent in Norwegian waterbodies. This limited number of species was observed by our eDNA approach with more than one amphibian spe- cies being detected in several samples. A study on the in- fluence of pH on amphibian diversity in Norway reported that R. temporaria was detected in all pH levels but there is a decrease in the frequency of B. bufo and T. vulgaris (Dolmen 1988), while R. arvalis and T. cristatus seemed to be sensitive to pH changes. Using our eDNA data and generalized additive models revealed a significant decrease in both R. arvalis (Estimate = -1.388, z-value = -2.958, p-value = 0.0031) and 7 cristatus (Estimate = -3.566, z-val- ue = -3.059, p-value = 0.0023), stronger than for L. vulgaris (Estimate = -0.979, z-value = -2.186, p-value = 0.0288) and R. temporaria (Estimate = -0.789, z-value = -2.013, p-value = 0.0441). Furthermore, our results confirm previous con- clusions that acidic areas are lower in amphibian diversity in the inland/highland region of Southern Norway as given by generalized additive models revealing a significant de- crease in overall species observations with decreasing pH (Estimate = -3.350, z-value = -3.644, p-value = 0.0003). Comparison between eDNA and traditional observa- tion monitoring methods Five amphibian species were detected by both traditional observation methods (ABD) and eDNA (Fig. 3). A much larger number of sites has been covered by ABD when compared to eDNA methods. Still, there were 133 sites where overlapping data was available for comparative analysis. An all data heatmap is shown in Fig. 3. Fisher exact and Pearson’s Chi-squared tests indicate a signif- icant association between eDNA and ABD in detecting four of the amphibian species while no significant asso- ciation could be observed for R. temporaria (Table 3). To conclude, correspondence of species observations https://mbmg.pensoft.net between data from a national species reporting ABD- based program and eDNA metabarcoding was significant but weak. This may reflect both an incomplete overlap in species detection and differences in the number of total detections of individual species, and the different time- frames (20 years vs 2 years). Lissotriton vulgaris Triturus cristatus Bufo bufo Rana temporaria Rana arvalis ND shared uniaue CDNA — uninue ADB Figure 3. Comparison of amphibian species detection between our metabarcoding results (2 years of data) and species inven- tories in the Norwegian species database “Artsdatabanken” (20 years of data). ND — number of sites where the respective species was not detected by any of either method; shared — the number of sites where the respective species was detected by both methods; unique eDNA — the number of sites where the respective species was detected only by eDNA; unique ABD - the number of sites where the respective species was found only in the database. Detailed statistics on correspondence between methods are given in Table 3. Table 3. Detailed statistics on the comparison of amphibian spe- cies detection between our metabarcoding results and species in- ventories in the Norwegian species database “Artsdatabanken”, whose results are shown in Fig. 4. Fisher exact Pearson’s test (odds ratio/ Chi-squared Correlation (R/p value) Species p-value) test (p-value) Rana arvalis 0.24/<0.001 5.37/0.017 0.016 Rana temporaria 0.13/<0.001 2.17/0.204 0.2 Bufo bufo 0.25/<0.001 3.70/0.006 0.004 Triturus cristatus 0.41/<0.001 inf/<0.001 <0.001 Lissotriton vulgaris 0.21/0.04 3.56/0.017 0.02 Discussion A major advantage of metabarcoding is that it facilitates the analysis of multiple taxonomic groups from the same sample. In principle, biodiversity across the entire tree of life present in a particular system can be assessed by DNA metabarcoding (Stat et al. 2017). In addition to amphibian Metabarcoding and Metagenomics 6: e85199 diversity, we tracked the distribution of the invasive chy- trid fungus Bd in the sampled water bodies. Similar studies have found Bd in various water bodies from other Nordic countries in Denmark and Sweden, and also in association with infected amphibians using genet- ic assays (Rosquist 2020). In this previous study, it was argued that skin swabs provided a higher likelihood of Bd detection. However, the comparison of detection proba- bility between the two methods in the Rosquist (2020) study was highly biased, as only a single water sample was taken from an individual system, while swabs were taken from multiple specimens (up to 97) during multi- ple sampling occasions spread over several years from a single system. We therefore attempted to assess detection probabilities of the two methods using comparable sam- pling efforts, i.e. using the same amount of time to take water and swab samples. Since it takes a much longer time to obtain swab samples than taking single water sam- ple, we only obtained very few swab samples. In terms of both time-efficiency and non-destructive handling, water sampling should be the preferred method. However, the low number of Bd positive waterbodies in Norway did not allow us to perform such a comparison. Still, con- firming results with diverse methods such as swab- and water-based detection, as well as multiple molecular as- says, 1S recommended since this can increase detection probability and confidence in positive detection. Lessons learned for a national monitoring program based on eDNA Screening amphibian occurrence by field techniques usually requires a suite of tools, which introduce several sources of errors. One of the errors is the high variability of amphibian detection in different seasons throughout the year (Heyer et al. 1993) and another is error in the sampling technique. For example, different types of am- phibian traps affect detection probabilities (Dodd 2010; Petitot et al. 2014). This has been proposed by previ- ous observations such as in neotropical streams in Bra- zil where both ABD and eDNA metabarcoding methods were used (Sasso et al. 2017). Our results confirm pre- vious observation that amphibians should be monitored in late spring/early summer in Scandinavia when using eDNA methods (Eiler et al. 2018). Thus, national moni- toring programs for amphibians based on eDNA are best performed from late May to early July. The high dispersion and low detection probability represents a challenge for monitoring programs in risk assessment, as large areas including high numbers of sys- tems need to be monitored. There is also a challenge to integrate the monitoring results in management strategies to protect amphibian diversity. To map large areas, such as across Norway, a random sampling strategy should be set up. In our case for the risk assessment of Bd, we focused on systems impacted by humans since the main pathway of introduction of Bd is related to human activ- ities (VKM report). Since sampling of a vast amount of SiS systems is necessary, efficient sampling techniques with low likelihood of contamination need to be applied. We chose a low-tech approach using a sterile beaker attached to an expandable (easy to clean) sampling pole, sterile syringes and cartridge filters. This kept the time spent per site to a minimum (20-45 minutes) including the assess- ment of water parameters with online sondes. For sample preservation, we recommend flash freezing of filters in liquid nitrogen on-site using a dry shipper and transfer to -80 °C freezers in the lab as this allows the long- term preservation of all kinds of environmental and tissue samples. This has been the “gold standard” for down- stream DNA and RNA analysis for at least three decades (Seutin et al. 1991; Reiss et al1995; Kilpatrick 2002; An- chordoquy and Molina 2007; Frampton and Sam 2008; Wong et al. 2012). For DNA extraction, we recommend commercial kits that are optimized for filter cartridges and remove potential PCR inhibitors while still keeping DNA yield high, i.e. the DNeasy PowerWater Sterivex kit. This is an important initial step of the analysis work- flow that needs to be kept constant as extraction protocols determine to a large extent the outcome of metabarcoding studies (Majaneva et al. 2018). For library preparation in the case of metabarcoding we strongly argue for a two- step PCR. There are several reasons: 1) Short adapters with undefined nucleotides between target primer and sequencing primer minimize biased PCR amplification when compared to a 1-step (fusion) PCR adding adapters and index, 2) a two-step PCR minimizes PCR cycles in the first step limiting PCR primer bias as cycle number can be kept low (also in the case of using a ligation based library preparation) (Thompson et al. 2002). 3) Indexing by PCR, which adds adapter sequence, minimizes index jumping when compared to ligation based library prepa- ration. 4) Compared to a ligation-based protocol the two- step PCR protocol is cheaper. Bioinformatic workflows for national monitoring pro- grams should denoise reads based on quality scores and optimize for each sequencing run, for example as imple- mented in the dada2 algorithm (Callahan et al. 2016) in- stead of clustering to improve species resolution. Using such methods allows discrimination of species with two nucleotide mismatches. Reference databases should be developed specifically for the national monitoring proj- ect taken from national species inventories including potential invasive species. This was rather simple in our case, with only seven amphibian species being present in Norway, but might be more challenging in species-rich taxonomic groups or locations. Species resolutions of the metabarcoding assay can be tested in silico using for ex- ample database-against-database homology searches and lowest common ancestor analysis (Huson et al. 2007). In the case of qPCR, we recommend the use of IPC to identify samples that exhibit PCR inhibition. These sam- ples can be flagged as problematic and run again after dilution which can lead to positive detection of targets. There is also a need to decrease LOD and LOQ, keeping amplification efficiency between 90-110%. The lowest https://mbmg.pensoft.net 314 number of replicate PCR reactions should be three, but this can be increase if allowed by budget constraints (Pig- gott et al. 2016). In addition, high Ct samples should al- ways be evaluated by sequencing to detect false positives. Besides these recommendations there are still many obstacles standing in the way of the implementation of eDNA-based biodiversity monitoring. International coor- dination is paramount to ensure comparability of data in time and space and to avoid unnecessary duplication of efforts and resulting delays in wide application. Aware- ness and knowledge of the new methods among practi- tioners 1s essential for implementation and to provide platforms for critical discussion on when, where and how the eDNA-based methods should be applied. Standardiza- tion efforts need to advance and cover all steps from sam- pling design to sequence data interpretation. Additionally, in order to allow for comprehensive mapping of genet- ic diversity, the number of published genomes (at least mitochondrial genomes) should cover all known species. Modelling and analysis tools should be applied and de- veloped in tandem to facilitate efficient decision making. Putting our results in the context of previous risk as- sessment for amphibians due to Bd Our study shows that eDNA surveys provide detailed in- sights into the distribution of amphibians and an associated pathogen, thus allowing the assessment of risks associated with the invasive chytrid fungus Bd. Reproducible detec- tion of Bd over multiple years strongly suggests that this pathogen is established at multiple sites in Southern Nor- way. Bd positive samples were obtained throughout an area of at least 1000 km? south of Oslo with Bd detection highly dispersed within the area, which is comparable to similar climatic regions in Sweden and the UK (Rosquist 2020). Multiple methods have been proposed to assess risks for biodiversity including Species Distribution Models (SDMs) combined with species specific biotic indices (Rodder et al. 2009), expert opinion (Non-native Risk Assessment scheme; Roy et al. 2014) and an integrative analytical approach combining ecological traits and an- thropogenic risk factors (Munstermann et al. 2021). Bd is known to be a particularly temperature and moisture dependent species (Lips et al. 2006, Woodhams et al. 2008). The in vitro growth optimum range of Bd is at 17— 25 °C, whereas temperatures higher than 29 °C, freezing and desiccation are lethal (Piotrowski et al. 2004). These findings are supported by observations in the field (Kriger et al. 2007). A previous study on the geographic extent of Bd’s climatic niche based on SDM suggested the pres- ence but low risk of Bd for anuran amphibian species in Norway (Rodder et al. 2009). In another recent study it was argued that transmission from an environmental Bd reservoir and climate change could increase the ability of Bd to invade new amphibian populations and increase their extinction risk. Still, a recent study concludes that Bd-induced extinction dynamics were far more sensitive to host resistance and tolerance than to Bd transmission https://mbmg.pensoft.net Omneya A. Osman et al.: eDNA monitoring of pathogen fungus Batrachochytrium dendrobatidis in Norway (Wilber et al. 2017). Infection-tolerant host species seem to dominate in Norway, as suggested by previous Europe- an studies (Smith 2014). Hence, the overall likelihood of a negative impact of Bd on indigenous amphibian species has been suggested (with moderate confidence) to be min- imal, while for an introduced frog species (P. /essonae) there was a moderate risk identified (Nielsen et al. 2019. Considering previous outbreaks of chytridiomycosis in Europe (Bosch and Martinez-Solano 2006) and extensive studies on Bd in similar climate regions, as well as its low prevalence in our study coupled with co-occurrences of various indigenous species over multiple years, it can be assumed that chytridomycosis outbreaks in Norway are unlikely. Still, an annual follow-up of Bd presence at least for the positive sites and nearby areas may be advisable, to keep an eye on the spread of this fungal pathogen. Acknowledgements We would like to thank Mats Topel and Thomas Larsson for their help in planning this work. Funding was provid- ed by the Norwegian Environment Agency. References Anchordoquy TJ, Molina MC (2007) Preservation of DNA. Cell Preservation Technology 4(4): 180-188. https://doi.org/10.1089/ cpt.2007.0511 Beaury EM, Finn JT, Corbin JD, Barr V, Bradley BA (2020) Biotic re- sistance to invasion is ubiquitous across ecosystems of the United States. Ecology Letters 23(3): 476-482. https://doi.org/10.1111/ ele.13446 Berger L, Speare R, Daszak P, Green DE, Cunningham A, Goggin CL, Slocombe R, Ragan MA, Hyatt AD, McDonald KR, Hines HB, Lips KR, Marantelli G, Parkes H (1998) Chytridiomycosis causes amphibian mortality associated with population declines in the rain forests of Australia and Central America. Proceedings of the Na- tional Academy of Sciences of the United States of America 95(15): 9031-9036. https://doi.org/10.1073/pnas.95.15.9031 Bosch J, Martinez-Solano I (2006) Chytrid fungus infection related to unusual mortalities of Salamandra salamandra and Bufo bufo in the Pefialara Natural Park, Spain. Oryx 40(1): 84—86. https://doi. org/10.1017/S0030605306000093 Boyle DB, Olsen V, Morgan JA, Hyatt AD (2004) Rapid quantitative detection of chytridiomycosis (Batrachochytrium dendrobatidis) in amphibian samples using real-time Taqman PCR assay. Diseases of Aquatic Organisms 60: 141-148. https://do1.org/10.3354/da0060141 Brannelly LA, Wetzel DP, Ohmer MEB, Zimmerman L, Saenz V, Rich- ards-Zawacki CL (2020) Evaluating environmental DNA as a tool for detecting an amphibian pathogen using an optimized extraction method. Oecologia 194(1-2): 267-281. https://doi.org/10.1007/ s00442-020-04743-4 Buxton AS, Groombridge JJ, Griffiths RA (2018) Seasonal variation in environmental DNA detection in sediment and water samples. PLoS ONE 13(1): e0191737. https://do1.org/10.1371/journal.pone.0191737 Metabarcoding and Metagenomics 6: e85199 Callahan BJ, McMurdie PJ, Rosen MJ, Han AW, Johnson AJ, Holmes SP (2016) DADA2: High-resolution sample inference from Illu- mina amplicon data. Nature Methods 13(7): 581-583. https://doi. org/10.1038/nmeth.3869 Cristescu ME, Hebert PDN (2018) Uses and misuses of environmental DNA in biodiversity science and conservation. Annual Review of Ecology, Evolution, and Systematics 49(1): 209-230. https://doi. org/10.1146/annurev-ecolsys-110617-062306 Deagle BE, Tollit DJ, Jarman SN, Hindell MA, Trites AW, Gales NJ (2005) Molecular scatology as a tool to study diet: Analysis of prey DNA in scats from captive Steller sea lions. Molecular Ecology 14: 1831-1842. https://doi.org/10.1111/).1365-294X.2005.0253 1.x Dodd CKJ (2010) Amphibian Ecology and Conservation: A handbook of techniques. Oxford University Press, Oxford. Dolmen D (1988) Coexistence and niche segregation in the newts 7rit- urus vulgaris (L.) and Triturus. cristatus (Laurenti). Amphibia-Rep- tilia 9(4): 365-374. https://doi.org/10.1163/156853888X00044 Dudgeon D, Arthington AH, Gessner MO, Kawabata Z, Knowler DJ, Lévéque C, Naiman RJ, Prieur-Richard AH, Soto D, Stiassny ML, Sullivan CA (2006) Freshwater biodiversity: Importance, threats, status and conservation challenges. Biological Reviews of the Cambridge Philosophical Society 81(02): 163-182. https://doi. org/10.1017/S 1464793 105006950 Duefias MA, Hemming DJ, Roberts A, Diaz-Soltero H (2021) The threat of invasive species to IUCN-listed critically endangered species: A systematic review, Global Ecology and Conservation 26: e01476. https://doi.org/10.1016/j.gecco.2021.e01476 Eiler A, Bertilsson S (2004) Composition of freshwater bacterial com- munities associated with cyanobacterial blooms in four Swedish lakes. Environmental Microbiology 6(12): 1228-1243. https://doi. org/10.1111/).1462-2920.2004.00657.x Eiler A, Drakare S, Bertilsson S, Pernthaler J, Peura S, Rofner C, Simek K, Yang Y, Znachor P, Lindstrém ES (2013) Unveiling dis- tribution patterns of freshwater phytoplankton by a next generation sequencing based approach. PLoS ONE 8(1): e53516. https://doi. org/10.1371/journal.pone.0053516 Eiler A, Lofgren A, Hjerne O, Norden S, Saetre P (2018) Environmental DNA (eDNA) detects the pool frog (Pelophylax lessonae) at times when traditional monitoring methods are insensitive. Scientific Re- ports 8(1): e5452. https://doi.org/10.1038/s41598-018-23740-5 Evans NT, Shirey PD, Wieringa JG, Mahon AR, Lamberti GA (2017) Comparative cost and effort of fish distribution detection via envi- ronmental DNA analysis and electrofishing. Fisheries 42(2): 90-99. https://doi.org/10.1080/03632415.2017.1276329 Frampton M, Sam D (2008) Evaluation of Specimen Preservatives for DNA Analyses of Bees. Journal of Hymenoptera Research 17: 195— 200. https://doi.org/10.1016/S0262-4079(08)62477-X Heyer R, Maureen AD, Mercedes F, Roy M (1993) Measuring and monitoring biological biversity: Standard Methods for Amphibians. Smithsonian Institution press. Hoorfar J, Malorny B, Abdulmawjood A, Cook N, Wagner M, Fach P (2004) Practical considerations in design of internal amplification controls for diagnostic PCR assays. Journal of Clinical Microbi- ology 42(5): 1863-1868. https://doi.org/10.1128/JCM.42.5.1863- 1868.2004 Huson DH, Auch AF, Qi J, Schuster SC (2007) MEGAN analysis of metagenomic data. Genome Research 7(3): 377-386. https://doi. org/10.1101/gr.5969107 eke Jorgensen D (2014) Not by Human Hands: Five Technological Te- nets for Environmental History in the Anthropocene. Environment and History 20(4): 479-489. https://doi.org/10.3197/09673401 4X14091313617163 Kamal R, Dhand Navneet K, Whittington Richard J, Plain Karren M (2017) PCR Inhibition of a Quantitative PCR for Detection of Mycobacterium avium Subspecies Paratuberculosis DNA in Feces: Diagnostic Implications and Potential Solutions. Frontiers in Micro- biology 8: e115. https://doi.org/10.3389/fmicb.2017.00115 Kilpatrick CW (2002) Noncryogenic preservation of mammalian tissues for DNA extraction: An assessment of storage methods. Biochemical Genetics 40(1/2): 53-62. https://doi.org/10.1023/A:1014541222816 Kylmus KE, Merkes CM, Allison MJ, Goldberg CS, Helbing CC, Hunt- er ME, Jackson C, Lance RF, Mangan AM, Monroe EM, Piaggio AJ, Stokdyk JP, Wilson CC, Richter CA(2019) Reporting the limits of detection and quantification for environmental DNA assays. En- vironmental DNA 2(3): 271-282. https://doi.org/10.1002/edn3.29 Kriger KM, Pereoglou F, Hero JM (2007) Latitudinal variation in the prevalence and intensity of chytrid (Batrachochytrium dendro- batidis) infection in Eastern Australia. Conservation Biology 21(5): 1280-1290. https://doi.org/10.1111/j.1523-1739.2007.00777.x Lance RF, Guan X (2019) Variation in inhibitor effects on qPCR as- says and implications for eDNA surveys. Canadian Journal of Fish- eries and Aquatic Sciences 77(1): 23-33. https://doi.org/10.1139/ cjfas-20 18-0263 Langenheder S, Comte J, Zha Y, Samad MS, Sinclair L, Eiler A, Lindstrom ES (2016) Remnants of marine bacterial communities can be retrieved from deep sediments in lakes of marine origin. Environmental Microbiology Reports 8(4): 479-485. https://doi. org/10.1111/1758-2229.12392 Laurance WF, McDonald KR, Speare R (1996) Epidemic disease and the catastrophic decline of Australian rainforest frogs. Conser- vation Biology 10(2): 406-413. https://doi.org/10.1046/j.1523- 1739.1996.10020406.x Lips KR, Brem F, Brenes R, Reeve JD, Alford RA, Voyles J, Carey C, Livo L, Pessier AP, Collins JP (2006) Emerging infectious disease and the loss of biodiversity in a neotropical amphibian communi- ty. Proceedings of the National Academy of Sciences of the Unit- ed States of America 103(9): 3165-3170. https://doi.org/10.1073/ pnas.0506889103 Longcore J, Pessier A, Nichols D (1999) Batrachochytrium dend- robatidis gen. et sp. nov., a Chytrid Pathogenic to Amphibians. Mycologia 91(2): 219-227. https://doi.org/10.1080/00275514.199 9.12061011 Longo AV, Rodriguez D, da Silva Leite D, Toledo LF, Mendoza Almer- alla C, Burrowes PA, Zamudio KR (2013) Correction: ITS1 Copy Number Varies among Batrachochytrium dendrobatidis Strains: Im- plications for qPCR Estimates of Infection Intensity from Field-Col- lected Amphibian Skin Swabs. PLoS ONE 8(10). https://doi. org/10.1371/annotation/e4c527e6-4b6f-4c 12-970c-a9e2e6ed8f79 MacKenzie DI, Nichols JD, Royle JA, Pollock KH, Hines JE, Bailey LL (2006) Occupancy Estimation and Modeling: Inferring Patterns and Dynamics of Species Occurrence. Elsevier, San Diego, California, USA. Mathieu C, Hermans SM, Lear G, Buckley TR, Lee KC, Buckley HL (2020) A Systematic review of sources of variability and uncertainty in eDNA data for environmental monitoring. Frontiers in Ecology and Evolution 8: 1-14. https://doi.org/10.3389/fevo.2020.00135 https://mbmg.pensoft.net 316 Omneya A. Osman et al.: eDNA monitoring of pathogen fungus Batrachochytrium dendrobatidis in Norway Majaneva M, Diserud OH, Eagle SHC, Hayjibabaei M, Ekrem T (2018) Choice of DNA extraction method affects DNA metabarcoding of unsorted invertebrate bulk samples. Metabarcoding and Metage- nomics 2: e26664. https://doi.org/10.3897/mbmg.2.26664 McKee AM, Spear SF, Pierson TW (2015) The effect of dilution and the use of a post-extraction nucleic acid purification column on the accuracy, precision, and inhibition of environmental DNA samples. Biological Conservation 183: 70—76. https://doi.org/10.1016/j.bio0- con.2014.11.031 Mendelson JR, Lips KR, Gagliardo RW, Rabb GB, Collins JP, Diffen- dorfer JE, Daszak P, Ibafiez DR, Zippel KC, Lawson DP, et al. (2006) Biodiversity. Confronting amphibian declines and extinctions. Sci- ence 313(5783): 48. https://do1.org/10.1126/science. 1128396 Munstermann MJ, Heim NA, McCauley DJ, Payne JL, Upham NS, Wang SC, Knope ML (2021) A global ecological signal of extinction risk in terrestrial vertebrates. Conservation Biology 36(3): e13852. https://doi.org/10.1111/cobi. 13852 Nielsen A, Dolmen D, Héglund J, Kausrud K, Malmstrom M, Taugbol A, Vralstad T, Ytrehus B, de Boer H, Hindar K, Kirkendall L, Nilsen E B, Rueness EK, Velle G (2019) Assessment of the risk to Nor- wegian biodiversity from the pathogenic fungi Batrachochytrium dendrobatidis (Bd) and Batrachochytrium salamandrivorans (Bsa). Opinion of the Panel on Alien Organisms and Trade in Endangered Species (CITES) of the Norwegian Scientific Committee for Food and Environment. VKM report 2019: 4. Parducci L, Matetovici I, Fontana SL, Bennett KD, Suyama Y, Haile J, Kjaer KH, Larsen NK, Drouzas AD, Willerslev E (2013) Molec- ular- and pollen-based vegetation analysis in lake sediments from central Scandinavia. Molecular Ecology 22(13): 3511-3524. https:// doi.org/10.1111/mec.12298 Petitot M, Manceau N, Geniez P, Besnard A (2014) Optimizing occu- pancy surveys by maximizing detection probability: Application to amphibian monitoring in the Mediterranean region. Trends in Ecology & Evolution 4(18): 3538-3549. https://doi.org/10.1002/ ece3.1207 Phillips JM (2004) Chapter 2 Real-time RT-PCR: What lies beneath the surface. In: Bustin SA (Ed.) A-Z of Quantitative PCR. International University Line, La Jolla, CA, USA, 47-85. Piggott MP, Wilson R, Banks SC, Marks CA, Gigliotti F, Taylor AC (2016) Evaluating the effects of laboratory protocols on eDNA de- tection probability for an endangered freshwater fish. Wildlife Re- search 35: 617-624. https://doi.org/10.1071/WR08040 Piotrowski JS, Annis SL, Longcore JE (2004) Physiology of Batrachoch- ytrium dendrobatidis, a chytrid pathogen of amphibians. Mycologia 96(1): 9-15. https://doi.org/10.1080/15572536.2005. 11832990 Pompanon F, Deagle BE, Symondson WO, Brown DS, Jarman SN, Taberlet P (2012) Who is eating what: Diet assessment using next generation sequencing. Molecular Ecology 21: 1931-1950. https:// doi.org/10.1111/).1365-294X.2011.05403.x Reiss RA, Schwert DP, Ashworth AC (1995) Field Preservation of Cole- optera for Molecular Genetic Analyses. Environmental Entomology 24(3): 716-719. https://doi.org/10.1093/ee/24.3.716 Rodder D, Kielgast J, Bielby J, Schmidtlein S, Bosch J, Garner TWJ, Veith M, Walker SF, Fisher MC, Lotters S (2009) Global Amphib- ian Extinction Risk Assessment for the Panzootic Chytrid Fungus. Diversity 2009(1): 52-66. https://doi.org/10.3390/d1010052 Rosquist G (2020) Coordinated actions against chytridiomycosis in the Nordic countries. The County Administrative Board Skane. https://mbmg.pensoft.net Roy HE, Preston CD, Harrower CA, Rorke SL, Noble D, Sewell J, Walk- er K, Marchant J, Seeley B, Bishop J, Jukes A, Musgrove A, Pearman D, Booy O (2013) GB Non-native Species Information Portal: Docu- menting the arrival of non-native species in Britain. Biological Inva- sions 16(12): 2495-2505. https://doi.org/10.1007/s10530-014-0687-0 Sasso T, Lopes CM, Valentini A, Dejean T, Zamudio KR, Haddad CFB, Martins M (2017) Environmental DNA characterization of amphibi- an communities in the Brazilian Atlantic forest: Potential application for conservation of a rich and threatened fauna. Biological Conser- vation 215: 225—232. https://doi.org/10.1016/).biocon.2017.09.015 Scalera R, Adams M, Galvan SK (2008) Occurence of Batrachochytri- um dendrobatidis in amphibian population in Denmark. Herpetolog- ical Review 2: 199-200. Seutin G, White BN, Boag PT (1991) Preservation of avian blood and tissue samples for DNA analyses. Canadian Journal of Zoology 69(1): 82-90. https://do1.org/10.1139/z91-013 Sinclair L, Osman OA, Bertilsson S, Elier A (2015) Microbial communi- ty composition and diversity: Evaluating the Illumina platform. PLoS ONE 10(2): e0116955. https://doi.org/10.1371/journal.pone.0116955 Skalski JR, Robson DS (1992) APPENDIX 1 - General Variance Com- ponent Formula. In: Skalski JR, Robson DS (Eds) Techniques in Wildlife Investigations. Academic Press, 215-216. https://doi. org/10.1016/B978-0-12-647675-0.50012-9 Smith F (2014) The epidemiology of the amphibian pathogen. Biology, 1-189. Stat M, Huggett MJ, Bernasconi R, DiBattista JD, Berry TE, Newman SJ, Harvey ES, Bunce M (2017) Ecosystem biomonitoring with eDNA: Metabarcoding across the tree of life in a tropical marine en- vironment. Scientific Reports 7(1): e12240. https://doi.org/10.1038/ $41598-017-12501-5 Taberlet P, Coissac E, Hajibabaei M, Rieseberg LH (2012) Environmen- tal DNA in ecology. Molecular Ecology 21(8): 1789-1793. https:// doi.org/10.1111/).1365-294X.2012.05542.x Taugbol A, Berum KM, Dervo BK, Fossey F (2021) The first detection of the fungal pathogen Batrachochytrium dendrobatidis in Norway with no evidence of population declines for great crested and smooth newts based on modelling on traditional trapping data. Environmen- tal DNA 3(4): 760-768. https://doi.org/10.1002/edn3.180 Thompson WL, White G, Gowan C (1998) Monitoring vertebrate pop- ulations. Academic Press, San Diego, CA. Thompson JR, Marcelino LA, Polz MF (2002) Heteroduplexes in mixed-template amplifications: Formation, consequence and elim- ination by ‘reconditioning PCR’. Nucleic Acids Research 30(9): 2083-2088. https://doi.org/10.1093/nar/30.9.2083 Tickner D, Opperman JJ, Abell R, Acreman M, Arthington AH, Bunn SE, Cooke SJ, Dalton J, Darwall JW, Edwards G, Harrison I, Hughes K, Jones T, Leclere D, Lynch AJ, Leonard P, McClain ME, Muru- ven D, Olden JD, Ormerod SJ, Robinson J, Tharme RE, Thieme M, Tockner K, Wright M, Young L (2020) Bending the Curve of Glob- al Freshwater Biodiversity Loss: An Emergency Recovery Plan. Bioscience 70(4): 330-342. https://doi.org/10.1093/biosci/biaa002 Valentini A, Pompanon F, Taberlet P (2009) DNA barcoding for ecol- ogists. Trends in Ecology & Evolution 24(2): 110-117. https://doi. org/10.1016/j.tree.2008.09.011 Valentini A, Taberlet P, Miaud C, Civade R, Herder J, Thomsen PF, Bellemain E, Besnard A, Coissac E, Boyer F, Gaboriaud C, Jean P, Poulet N, Roset N, Copp GH, Geniez P, Pont D, Argillier C, Baudoin JM, Peroux T, Crivelli AJ, Olivier A, Acqueberge M, Le Brun M, Metabarcoding and Metagenomics 6: e85199 Moller PR, Willerslev E, Dejean T (2016) Next-generation moni- toring of aquatic biodiversity using environmental DNA metabar- coding. Molecular Ecology 25(4): 929-942. https://doi.org/10.1111/ mec. 13428 Wilber MQ, Knapp RA, Toothman MH, Briggs CJ (2017) Resistance, tolerance and environmental transmission dynamics determine host extinction risk in a load-dependent amphibian disease. Ecology Let- ters 9(9): 1169-1181. https://doi.org/10.1111/ele.12814 Wilcox TM, McKelvey KS, Young MK, Jane SF, Lowe WH, Whiteley AR, Schwartz MK (2013) Robust detection of rare species using en- vironmental DNA: The importance of primer specificity. PLoS ONE 8(3): e59520. https://doi.org/10.1371/journal.pone.0059520 Wong PB, Wiley EO, Johnson WE, Ryder OA, O’Brien SJ, Haussler D, Koepfli K-P, Houck ML, Perelman P, Mastromonaco G, Bentley AC, Venkatesh B, Zhang Y, Murphy RW (2012) Tissue sampling methods and standards for vertebrate genomics. GigaScience 1(1): 8-20. https://doi.org/10.1186/2047-217X-1-8 Woodhams DC, Alford RA, Briggs CJ, Johnson M, Rollins-Smith LA (2008) Life-history trade-offs influence disease in changing cli- mates: Strategies of an amphibian pathogen. Ecology 89(6): 1627— 1639. https://doi.org/10.1890/06-1842.1 Woodhams DC, Kenyon N, Bell SC, Alford RA, Chen S, Billheimer D, Shyr Y, Rollins-Smith LA (2010) Adaptations of skin peptide de- fenses and possible response to the amphibian chytrid fungus tn pop- ulations of Australian green-eyed treefrogs, Litoria genimaculata. Diversity & Distributions 16(4): 703-712. https://doi.org/10.1111/ j.1472-4642.2010.00666.x Zinger L, Gobet A, Pommier T (2012) Two decades of describ- ing the unseen majority of aquatic microbial diversity. Molec- ular Ecology 21(8): 1878-1896. https://doi.org/10.1111/j.1365- 294X.2011.05362.x Supplementary material 1 Table S1 Authors: Omneya A. Osman, Johan Andersson, Pedro M. Mar- tin-Sanchez, Alexander Eiler Data type: Excel file. Explanation note: List of samples from 2019 with raw se- quences, sampling volume and environmental parameters. Samples with positive (+) and negative (-) Bd results and amphibian species detected are represented. Copyright notice: This dataset is made available under the Open Database License (http://opendatacommons.org/li- censes/odbl/1.0/). The Open Database License (ODbL) is a license agreement intended to allow users to freely share, modify, and use this Dataset while maintaining this same freedom for others, provided that the original source and author(s) are credited. Link: https://do1.org/10.3897/mbmg.6.85199.suppl1 317 Supplementary material 2 Table S2 Authors: Omneya A. Osman, Johan Andersson, Pedro M. Mar- tin-Sanchez, Alexander Eiler Data type: Excel file. Explanation note: List of samples from 2020 with raw sequenc- es, sampling volume and environmental parameters. Samples with positive (+) and negative (-) Bd results and amphibian species detected are represented. All swabs collected from positive Bd sites were sequenced while one swab was only sequenced from the negative &d sites. Copyright notice: This dataset is made available under the Open Database License (http://opendatacommons.org/li- censes/odbl/1.0/). The Open Database License (ODbL) is a license agreement intended to allow users to freely share, modify, and use this Dataset while maintaining this same freedom for others, provided that the original source and author(s) are credited. Link: https://do1.org/10.3897/mbmg.6.85199.suppl2 Supplementary material 3 Table S3, Figures S1, 82 Authors: Omneya A. Osman, Johan Andersson, Pedro M. Mar- tin-Sanchez, Alexander Eiler Data type: MS Word file. Explanation note: Table S3. Summary statistics on the limit of detection (LOD) and the limit of quantification (LOQ) of the Bd assays. LOD and LOQ: 10, and 1 gene copies for commercial kit and 1 and 0.1 genome equivalent for Boyle TaqMan assay. Figure S1. Sampling scheme of the sam- pling in 2019 showing the distribution of sampling locations during the three sampling occasions (A). Distribution of Bd as assessed by the analysis of water samples, swabs, and tadpole DNA around Oslo (B) and Southern Norway and Central Sweden (C). In the latter map the distribution of all samples obtained in this study (data from 2019 and 2020) are shown. Additional data was obtained from Rosquist (2020) and Taugbol et al. (2021). Full symbols indicate detection of Bd by the TaqMan assays while open symbols indicate nega- tive results. Figure S2. Scatter plot of Cq values of standard and samples multiplexed with the internal extraction control. The figure only shows three samples with Cq values more than 27 and less than 50, while the other three samples had Cq value more than 50, not detected by qPCR. Copyright notice: This dataset is made available under the Open Database License (http://opendatacommons.org/li- censes/odbl/1.0/). The Open Database License (ODbL) is a license agreement intended to allow users to freely share, modify, and use this Dataset while maintaining this same freedom for others, provided that the original source and author(s) are credited. Link: https://do1.org/10.3897/mbmg.6.85199.suppl3 https://mbmg.pensoft.net