Field-work was conducted in north-eastern Alberta and the southern Northwest Territories, Canada, predominately within the boundaries of Wood Buffalo National Park (Figure 1). The area is part of the Taiga Plains and Boreal Plains Ecoregions, featuring a distinctive landscape consisting of extensive areas of interconnected river systems, wetlands and lakes, surrounded by boreal forest (59–61); (see Figures 2A–D). The area is characterized by a sub-humid climate with low mean annual air temperatures (62), a high rate of environmental disturbances (e.g., wildfires and floods), and ecoregion specific characteristics such as (semi-) permafrost (63, 64). Field reconnaissance surveys began the first week of May in 2016 and 2017 when spring melt occurs in the area, and continued until the end of June. This timing maximized our ability to acoustically detect all amphibian species expected in the area as well as physically encounter multiple life stages of each species. To determine the species present at each site and to estimate their abundance, we conducted a total of 221 visual, auditory, and dip netting surveys (75 in 2016; 146 in 2017) of suitable amphibian habitats such as floodplains, wetlands and meadows. Opportunistic surveys were conducted along the margins of the respective wetland, and down to an approximate water depth of 1 m. Three species of amphibians were encountered: wood frogs (Rana sylvatica, WF), boreal chorus frogs (Pseudacris maculata, CF), and Canadian toads (Anaxyrus hemiophrys, CT). A species was determined to be present at a site if the species was caught, or visually/auditorily identified even if we were unable to ultimately capture and collect tissue samples from the species at a site. Tadpoles of a given species were considered to be present when either visually confirmed, or when breeding behavior in adults was observed, and eggs were found. Throughout our extensive surveys we encountered three compositions: WF only, WF/CF, or all three species. Sites with CF or CT only, WF/CT, or CF/CT were never encountered.

A total of 18 unique sites were chosen for sampling (Figure 1; see Appendix I for full names and coordinates). Tissue sampling of terrestrial individuals (fully metamorphosed young of the year, juveniles, and adults) was conducted at 13 sites between May 15 and June 23, 2016. Tissue sampling of aquatic stages [tadpoles and metamorphs; Gosner stages 26 – 42 (65)] was conducted at 16 sites (11 identified in 2016 + 5 new in 2017) between June 14 and June 29, 2017. Due to logistical constraints, two of the sites sampled in 2016 were not sampled in 2017. To allow for comparisons of amphibian abundance among sites, we sampled each wetland using the same methodology, conducting opportunistic sampling for all species potentially present in the area. We quantified relative abundance as Catch Per Unit Effort (CPUE): the number of individuals encountered per unit of search time (here Ind./hour) according to (66). To allow for sufficient sample sizes for pathogen detection and statistical analysis, we selected sites with a minimum abundance of 10 terrestrials or 30 aquatic individuals per 1 h search time. If the required sample size could not be achieved at the first sampling, sites were sampled a second time within a 2 week time-frame of the initial sampling. To characterize abiotic characteristics of our study sites, temperature, water pH, and conductivity were recorded during each sampling event using an YSI 556 MPS (YSI Inc., Yellow Springs, United States), and were averaged if a site was sampled more than once. To more completely assess temperature profiles in our study wetlands, we deployed temperature data-loggers (Onset HOBO Fresh Water Data Logger U24-001) at three wetlands that differed in terms of size and other abiotic characteristics. The loggers were deployed 10 May−1 July, 2017, and were attached to a stationary pole approximately 50 cm below the water surface. Due to the landscape in the area, with vast areas of karst topography (i.e., sinkholes), it was not feasible to determine wetland depth.

Individual amphibians were caught by hand or by dip-netting and then held individually in new, single use plastic bags containing air and pond water (for aquatic stages), or wet vegetation (for terrestrial stages), until processed (<30 min). Care was taken to replenish air in the bags frequently and to protect the animals from direct sun and overheating. Aquatic individuals were measured (snout-vent-length; SVL) and staged (65) in their respective bag. Terrestrial individuals were weighed (g) and their SVL was measured. To test for RV infections, we non-lethally collected a ~1 mm toe-clip (terrestrials) or tail-clip (aquatic stages) from each individual. To prevent cross-contamination and spread of pathogens among individuals during handling, a new pair of nitrile gloves and a new scalpel blade was used for processing each animal, and all surfaces and tools were regularly disinfected with 99% isopropyl. Tissue samples were stored individually in 1.5 ml microcentrifuge tubes (filled with 95% ethanol, Eppendorf AG, Hamburg, Germany). To prevent accidental transfer of pathogens among sites, all gear and equipment that came into physical contact with the processed individuals, including nets, boots and waders were disinfected for 5 min with a 10% bleach solution, rinsed with clean tap water, and allowed to dry before being used at the next site [see (67)]. We solely used non-lethal tissue samples for ethical reasons and because they are easy to collect and store under remote field conditions. Tail and toe clips have been shown to be a viable alternative to lethal sampling and are in a generally high agreement with results from corresponding organ samples, but often slightly underestimate infection prevalence (44, 68–71). There are no known long-lasting negative effects on the amphibians when a single phalanx or tail clip is collected aseptically (68–70), and non-lethal sampling significantly decreases the impact on biodiversity, particularly when rare or threatened species are studied (44, 69, 72).

Ranavirus detection and quantification was performed using quantitative PCR (qPCR) following the cycling conditions described by Leung et al. (73), and quantification described by Hoverman et al. (11). DNA was extracted from toe and tail clips using Qiagen DNEasy® Blood and Tissue kits according to manufacturer specifications (QIAGEN Inc., Valencia, CA, USA). A Synergy H1 Hybrid Multi-Mode Reader was used to quantify concentration of genomic DNA in each sample (BioTek, Winooski, VT, USA). Briefly, the qPCR mixture contained 10 μl TaqMan Universal PCR Master Mix 2X (Thermofisher Scientific, Waltham, MA, USA), 1 μl forward primer MCPRV_F-5GTCCTTTAACACGGCATACCT3 (10 μM), 1 μl reverse primer MCPRV_R-5ATCGCTGGTGTTGCCTATC3 (10 μM), and 0.05 μl TaqMan probe MCP_NFQ-5TTATAGTAGCCTRTGCGCTTGGCC3 (100 μM). Subsequently, we added 250 ng of template DNA and DNA grade water to a final reaction volume of 20 μl. The qPCR was performed using an Mx3005P QPCR System (Agilent Technologies, Santa Clara, CA, USA). Cycling conditions were 50°C for 2 min, 95°C for 10 min, and 50 cycles of 95°C for 15 s, and 60 C for 30 s. Samples were run in duplicate in 96-well-plates. Each 96-well-plate included a no-template control (molecular grade water), and a serial dilution of a known quantity of cultured FV3 (1 x 10⁶-10¹ copies/μl) to create a standard curve with precise fit (R² >0.95). Individuals were considered positive if the target DNA was clearly amplified in both of the duplicates (i.e., surpassed the respective cycle threshold). If only one of the two wells amplified, a third reaction was conducted to either confirm or dismiss the previous positive result. Standard curves from each plate were used to calculate viral load following Yuan et al. (74) and reported as viral copies/250 ng of gDNA, as recommended by Gray et al. (69). Due to large standard deviations, viral loads are reported as log10 copies/250 ng gDNA (hereafter “logVL”). Primers used in our study target a consensus Major Capsid Protein (MCP) sequence of the major amphibian associated RV isolates: FV3: GenBank No AY548484, (75); TFV: AF389451, (76); CMTV: JQ231222, (77); EHNV, FJ433873, (78), and allow for a high analytical sensitivity when used in combination with a TaqMan probe (73). Pathogen identification and infection was confirmed through immunohistochemistry and in-situ hybridization on a subset of collected tissues (79).

Statistical analyses were conducted using R 3.3.3 (80) in RStudio Version 1.1.383 (81). Sample sizes allowed for reasonable prevalence estimates of RV at each site [95% confident that prevalence was < ~5% if no animals tested positive at a given site; (82, 83)]. We conducted analyses on the community level including all three species. Analyses at the individual and population level were conducted for WF only because single species sites were found only for this species. Aquatic and terrestrial life stages were analyzed separately. We used Pearson's Chi-squared tests to compare prevalence (percentage of positive individuals in a community) among the three different community compositions and used post-hoc pairwise comparisons to determine significance. We conducted one-way analyses-of-variance (ANOVAs) to test for differences in mean viral load among communities. We used Tukey's HSD test to detect differences among means for ANOVAs where the overall null hypothesis was rejected. To analyze RV prevalence among populations and communities, we performed beta regressions using the package “betareg” (84) for beta-distributed, dependent percentage data. To assess the influence of predictor variables on viral loads at the individual, population and community level, we performed generalized linear mixed model regressions fitted with a Poisson family and log link, using the package “lme4” (85). Each level of analysis included site-specific conductivity, pH, and species richness (WF, WF/CF, WF/CF/CT), as well as total relative community abundance (number of individuals caught per person-hour) as fixed effects. Due to the wide range of measurements for conductivity, it was scaled into seven categories (low to high conductivity). Site was treated as a random effect in all mixed model analyses. We used the dredge function in the R package “MuMIn” to create a sub-set of best-supported models of all possible sub-models from ranavirus presence and infection global models (86). The best supported models were compared separately using Akaike's Information Criterion [AIC; (87, 88)]. Due to the overall small sample sizes, we used AIC corrected for small sample sizes (AICC) in all our analyses (87, 89). Model comparison and averaging were performed using R packages ‘MuMIn’ (86) and ‘AICcmodavg’ (88), and can be found in the Supplementary Material (see Appendices II, III).

Three amphibian species were encountered in three types of community composition (Figure 1, Tables 1, 2), and we sampled a total of 1,244 individuals across all species and sites (Table 3). Wood frogs occurred at all sites where amphibians were found and were the only species we observed by itself (n = 2 sites; Table 4). The most commonly encountered community consisted of WF and CF (n = 10 sites; Table 4). If the sites featured large, shallow areas and pools, and provided nearby suitable hibernacula sites (e.g., sand dunes and pits), we often found CT as well (n = 6 sites). Abundance varied greatly among the three amphibian species: in terrestrial phases (2016), an average of 12 WF were caught per hour (range 3–37 individuals/hour) across 13 sites, whereas an average of 4 CF were caught per hour (range 1–9 individuals/hour) across nine sites (Table 1). For terrestrial stages, CT were caught at three sites at an extremely low average rate (1 individual/hour) and were observed at only four sites overall (Tables 1, 4). For aquatic stages (2017), an average of 40 WF were caught per hour across 16 sites (range 4–60 individuals/hour) whereas an average of 14 CF were caught per hour across 11 sites (range 2–41 individuals/hour; Tables 2, 4). Aquatic stages of Canadian toads were caught at two sites at an average rate of 9 individuals per hour (range 8–10 individuals/hour). We observed CT breeding activities at four additional sites, but could not find any aquatic stages during sampling. Although abundance varied among species, CF abundance did not vary significantly across the communities in which they were found (WF/CF and WF/CF/CT) in either terrestrial or aquatic life stages. Wood frog abundance however was lower in WF/CF/CT communities as compared to the two less complex communities in both terrestrial and aquatic stages. Across all our study sites, the pH ranged from 6.18 to 7.98 in 2016 and 6.90 to 8.71 in 2017 (Tables 1, 2). Conductivity ranged from 73 μS/cm to 8,290 μS/cm in 2016 and 70 μS/cm to 7,132 μS/cm across all sites (Tables 1, 2). Although we collected water temperature at the time of tissue sampling, it was not included in our analyses. Single-point measurements of temperature under field conditions are highly biased by sampling time and environmental conditions. This observation was further supported by the marked daily variability in our temperature data collected at three selected wetlands, with the means not differing from one another (see Appendix IV for data-logger recordings). Individual measures such as snout-vent length, weight, and Gosner stage for aquatic individuals, did not have any statistical significance in preliminary model comparisons and were excluded from the final analyses.

Ranavirus was detected in individuals at all sites sampled (2016, n = 13 and 2017, n = 16; Tables 1, 2) and in all species (Table 3). The overall prevalence of RV across all species and sites was 16% in terrestrial individuals (2016) and 32% in aquatic stages (2017; Table 4). The mean infection prevalence in WF was 19% in terrestrial stages, and 37% in aquatic stages (Table 3). In CF, the mean infection prevalence was 9% in terrestrial stages, and 49% in aquatic stages (Table 3). We did not detect any infected terrestrial CT and the mean infection prevalence in CT aquatic stages was 16% (Table 3). It is notable that all RV-positive CT aquatic stages were detected at a single site during a RV die-off involving WF and CF [Tables 2–4; see Forzán et al. (79) for histological data on a subset of collected individuals]. Among all terrestrial individuals sampled, only 1% exhibited typical signs of ranavirosis such as hemorrhages at the toe-tips and ventral surfaces of the abdomen (44); see Table 1 for population level data). By contrast, 29% of all aquatic individuals sampled exhibited typical signs of RV infection such as hemorrhages in developing hind legs, at the base of the tail, and on the abdomen (44); (98 WF, 63 CF, 1 CT; see Table 2 for population level data), In addition to the aforementioned signs, aquatic individuals (n = 72) sampled during two die-off events in 2017 commonly exhibited edema. We found no effects of specific sampling dates on prevalence data for aquatic individuals (Z = 1.61, p = 0.11) or terrestrial individuals (Z = 0.03, p = 0.98), and we did not detect any temporal trends (Appendices V–VII).

The mean viral loads among positive terrestrial WF and CF individuals were low (logVL < 2.83), and RV was not detected in terrestrial CT (Tables 1, 3). Viral loads in aquatic stages varied greatly, with more than half of the individuals (61%) exhibiting higher mean viral loads in comparison to terrestrials (logVL < 3.65, Tables 2, 3). Wood frog and CF aquatic stages sampled during mass die-off events (n = 72) exceeded the commonly observed viral loads (38 individuals with logVL > 3.65 and 34 individuals with logVL > 4.65).

The percentage of RV positive individuals in a community increased with the number of species, but not significantly (χ² = 2.19, df = 2, p = 0.335; Table 4, Figure 3A). At sites where WF occurred alone, we observed a 9% mean infection prevalence, increasing to 15% when CF were also present in the community. When three species were present, the mean community prevalence increased to 23% per site on average with regression analyses showing a significant effect of species richness in explaining the observed pattern, but no effects of other variables (Table 5, Appendix II). We also observed a significant increase in mean community viral loads linked to species richness (F_{2, 7} = 6.82, p < 0.05), while other variables had no explanatory power (Tables 4, 5, Figure 3B).

At the population level (WF only), we observed a significant increase in infection prevalence from 9 (WF) to 16 (WF/CF), to 28% (WF/CF/CT; χ² = 8.5, df = 2, p < 0.05; Table 4), but did not find any significant relationships of population prevalence and abiotic variables (Table 5). We also found a significant positive relationship of mean viral loads in terrestrial WF stages with the number of syntopic species in the community, at the individual and population levels (Table 5).

In the aquatic amphibian life stages we observed that at the community level, ranavirus prevalence nearly doubled with each additional species (Table 4, Figure 3C). Wood frog aquatic stages in sites with no other amphibian species exhibited a low RV infection prevalence of 13% while in communities consisting of WF and CF we found significantly higher infection prevalence across the two species (28%; Table 4). When CT were also present in the communities, the mean community infection prevalence increased significantly to 55% across the three species (χ2 = 45.85, df = 2, p < 0.001; Table 4). Regression analyses showed significant effect of species richness on RV prevalence at the community level, and no effects of other variables (Table 6, Appendix III). We also observed a significant increase in mean community viral loads (F_{2, 13} = 1.963, p < 0.05) linked to species richness, while other variables had no explanatory power (Tables 4, 6, Figure 3D).

At the population level (WF only), we observed an increase in prevalence from 13% (WF) to 29% (WF/CF), and to 53% (WF/CF/CT; χ² = 22.437, df = 2, p < 0.001). Regression analyses showed a significant positive correlation of population prevalence in WF aquatic stages with the number of syntopic species in the community, with no significant effect of other variables (Table 6). The same patterns were observed for viral loads in WF aquatic stages at individual and population levels (Table 6). When data collected during mass die-off events were excluded from the statistical analysis, species richness remained the main predictor for prevalence and viral loads on community level and WF population level (see Appendix VIII).