Evaluating Amphibian Declines with Site Revisits and Occupancy Models: Status of Montane Anurans in the Pacific Northwest Usa

Amphibian declines have been reported in mountainous areas around the western USA. Few data quantify the extent of population losses in the Pacific Northwest, a region in which amphibian declines have received much attention. From 2001–2004, we resurveyed historical breeding sites of two species of conservation concern, the Western Toad (Bufo [=Anaxyrus] boreas) and Cascades Frog (Rana cascadae). We detected B. boreas breeding at 75.9% and R. cascadae breeding at 66.6% of historical sites. When we analyzed the data using occupancy models that accounted for detection probability, we estimated the current use of historically occupied sites in our study area was 84.9% (SE = 4.9) for B. boreas and 72.4% (SE = 6.6) for R. cascadae. Our ability to detect B. boreas at sites where they were present was lower in the first year of surveys (a low snowpack year) and higher at sites with introduced fish. Our ability to detect R. cascadae was lower at sites with fish. The probability that B. boreas still uses a historical site for breeding was related to the easting of the site (+) and the age of record (None one of the variables we analyzed was strongly related to R. cascadae occupancy. Both species had increased odds of occupancy with higher latitude, but model support for this variable was modest. Our analysis suggests that while local losses are possible, these two amphibians have not experienced recent, broad population losses in the Oregon Cascades. Historical site revisitation studies such as ours cannot distinguish between population losses and site switching, and do not account for 167 colonization of new habitats, so our analysis may overestimate declines in occupancy within our study area.


Introduction
Amphibian declines have been reported from many parts of the world, with concentrations in Australia, Central America, and western North America [1].Declines in species from mountainous areas such as the western USA have been difficult to quantify and their causes can be complex [2][3][4].The Cascade Range spans the latitudinal length of the northwestern USA, and is an example of a montane region that has not experienced the direct anthropogenic disturbance present in most western lowlands, but where amphibian declines have been reported [2,3,5,6].Most of the research on putative amphibian declines in the Cascade Range has focused on experimental investigations of hypothesized stressors such as UV radiation (e.g., [3]).There has been relatively little work analyzing amphibian status or trends over the region.
One approach to quantifying amphibian declines is to resurvey sites of known historical occurrence to assess whether the species remains present [7,8].Comparisons of extant sites with those where the species is lost allows development of hypotheses about the potential causes of any declines.Studies of historical sites have helped clarify the scope of declines in some other western amphibians [2,4,[9][10][11].However, many of these studies use single surveys at each site and have not attempted to account for imperfect detection of target species (but see [11]).Occupancy modeling is a relatively new technique that accounts for imperfect and variable detection probabilities [12].The procedure can relate occupancy and detectability to site-and observation-level covariates and provides a means to produce unbiased estimates of the proportion of sites occupied.
Western Toad (Bufo [=Anaxyrus] boreas) and Cascades Frog (Rana cascadae) are pond-breeding anurans that were historically widespread in the Cascade Range in Oregon and the Mt.Lassen and Trinity Alps areas of northern California [13,14].Declines of both species have been reported and both species have figured prominently in recent literature on amphibian declines [2,3,5,15,16].Population losses have been documented in portions of both species' ranges in California [2,4,10,17] and, for toads, in Colorado [18] and eastern Oregon [11].However, quantitative field data addressing potential declines are sparse in the Oregon Cascades.
Our objective was to use occupancy modeling and resurveys of historical breeding sites to assess the status of B. boreas and R. cascadae in the Oregon Cascade Range.We also examined potential links between occupancy and site characteristics that are often associated with population persistence.

Study Region and Amphibians
The Cascade Range extends over 1,100 km from its intersection with the Rocky Mountain chain near the Canadian border to its southern limit near Mt.Lassen in northern California.The main axis of the Cascades runs north-south and has several volcanic peaks >3,000 m elevation.Our study sites are within the Oregon section (ca.350 km north to south) of the Cascade Range.This region often gets heavy snow between November and April, and has short dry summers (June to September).Most of our survey sites are within relatively closed forest of Douglas fir (Pseudotsuga menzeisii), western hemlock (Tsuga heterophylla), mountain hemlock (Tsuga mertensiana), and lodgepole pine (Pinus contorta).Some sites are in open subalpine forest of the latter two tree species.Many of the larger lakes and streams in our study area have introduced brook (Salvelinus fontinalis) and rainbow trout (Oncorhynchus mykiss) [19].Upper elevations of the Oregon Cascades were glaciated in the Pleistocene and are populated with lakes, ponds and wetlands [20].Today much of the crest of the Cascades is in federal Wilderness areas, and the lower and middle elevations are used mainly for timber production.The Willamette Valley and the Umpqua and Rogue valleys occur along the western foot of the Cascades and support intensive agriculture and the largest population centers in Oregon.
Bufo boreas and Rana cascadae generally breed in ponds and lakes at mid-and upper elevations in the Oregon Cascades (ca.600 to >2,100 m) [14].Breeding by both species takes place soon after snow melt, and larvae transform in two to four months in our study area [14].Adults of both species can be long lived (>7 yr) [21] and are capable of extensive movements away from breeding ponds [14,22,23].Rana cascadae ranges from just south of the Canadian border to Mt. Lassen in northern California.Disjunct populations occur west of the Cascade Range in the Olympic Peninsula (Washington) and the Trinity Alps (California).Bufo boreas (which recent genetic data suggests may be a complex of taxa rather than single species [24]) is found from Alaska south along the Coast Ranges, Rocky Mountains, and Cascade Range, with the southernmost populations in New Mexico and California [22].

Data Collection
Our objective was to estimate the proportion of historical breeding sites that were currently used for breeding by B. boreas and R. cascadae.We also examined the influence of selected habitat attributes and stressors that are hypothesized to contribute to population status: vegetation cover [25], non-native fish [26], hydrologic permanence [27], pond origin [11], roads [28], and proximal forest cover [29].We included UTM location data (EASTING, NORTHING in meters) as occupancy covariates since declines appear to be more severe in southern or windward portions of ranges for other western amphibians [2,4,9,30].
We compiled historical site records from museum collections, graduate theses, publications, our own studies in prior years, and personal communications with herpetologists and agency biologists with direct experience with the target species.We limited the dataset to records of breeding (eggs, larvae, metamorphic individuals).Field surveys were conducted by 2-4 person crews in summers of 2001-2004.Sites were surveyed 1-4 times but were never surveyed >1 time in a single year.Field variables were assessed during the first survey at each site (Table 1) and are indicated in the text by capital letters.Crews used visual searches and dip-netting along shore and throughout wadeable habitat to detect amphibians and FISH.We consulted US Forest Service and Oregon Department of Fish and Wildlife fishery biologists with knowledge of study sites to complement our field assessments of FISH.We visually estimated the percentage of site perimeter that had FOREST within 50 m of shoreline.We used field visits, aerial photography, and USGS topographic maps to quantify whether sites had PAVED or ANY ROAD (paved or unpaved) within 500 m.We visually assessed whether sites were MANMADE (e.g., impoundments or excavations), and estimated the percentage of each site that had aquatic VEGETATION.We called sites PERMANENT if they had aquatic vegetation and depth (>1 m) suggesting they would not dry during most low snow years.We considered AGE of each record to be the number of years between the end date of the historical record and 2001, which was the first year of our study.Sites were scored '1' for each survey on which eggs, larvae, or metamorphs of the target amphibian were found, and '0' when none of those life stages were encountered.

Analysis
Besides reporting the naï ve occupancy rate (the proportion of sites where we detected the target species at least once), we used the single-season model in program PRESENCE [12] to estimate the proportion of historical sites that were used for breeding at least one year over the course of our study.A non-detection for a target species at a surveyed site can mean the species was not present at the site or we failed to detect it when it was present.PRESENCE uses a joint likelihood model to estimate the probability of missing a species when it is present at the site (p, hereafter 'detectability') and the probability that a site is occupied (Ψ).This is typically done using data from repeat observations made over a period of time during which site occupancy is assumed to be constant (e.g., summer months between breeding and metamorphosis).In this way, a non-detection from a site with at least one detection can be treated as a false negative and the probability of failing to detect a species, conditional on presence, over multiple observations can be estimated.In our analysis, repeat observations were made over multiple seasons and we cannot assume that site occupancy was constant during that period.Because of this, the definitions of the parameters we estimate differ somewhat from the norm.Our estimate of Ψ represents the probability that a species used a site for breeding in any year of the study.Our estimate of p it has two multiplicative components and represents the probability we detected breeding at a site, given that breeding occurred at that site in that year.The first is the probability that a species that bred at site i during any year of our study, bred at a site i in year t.The second is the probability that the species was detected at site i in year t conditional on presence at site i in year t.The practical implication of the multiplicative components of p is that our estimate of the proportion of sites used by a species takes into account irregular use of sites for breeding.To the extent that this sporadic breeding use occurs, our estimate of occupancy will be higher than the proportion of sites where breeding actually occurred in any one year.
Prior to analysis, we screened predictor variables for Pearson correlation >0.70 [31].We normalized all continuous variables to have a mean of zero and a standard deviation of 1 for analysis in PRESENCE.Occupancy models assume that detectability is either constant or that any variation is modeled using covariates.We suspected that detectability could be affected by 2 site (FISH, VEGETATION) and 1 survey (YEAR1) characteristics.We used the term YEAR1 to test for effects of the 2001 drought on detectability.Precipitation during the 2001 water year (October 2000-June 2001) at 5 SNOTEL stations in our study area was 60-68% of respective 25-36 yr station averages (http://www.wcc.nrcs.usda.gov/cgibin/state-site.pl?state=OR&report=precsnotelmon).Drought is thought to affect breeding in B. boreas and several other pond-breeding amphibians [8,10].YEAR1 was coded 1 for surveys in 2001 and 0 for surveys in all other years.Because of this hypothesized variation in detectability, we began our analysis by using AIC (Akaike's Information Criterion) to select the most parsimonious combination of covariates to account for variation in detectability among sites and surveys.To assess detectability covariates for both our amphibians, we compared the null model with models for direct effects of habitat (effects of FISH, VEGETATION) and drought (YEAR1).For R. cascadae, we added FISH ×VEGETATION because, in contrast to B. boreas, R. cascadae is palatable to fish and can modify risky behaviors such as cover use in presence of fish [26,32,33].For each species, the model with the lowest AIC was considered the best detectability model and the covariates it contained were included in all subsequent analyses to account for variation in detectability.
Next, we used AIC to rank models comprised of all combinations of factors that could affect occupancy at 3 spatial scales: pond, local landscape, and regional.We analyzed each scale separately for both species: pond scale with covariates FISH, VEGETATION, PERMANENT and MANMADE; landscape scale with FOREST, PAVEDROAD, and ANYROAD for both species; regional scale with ELEVATION, NORTHING, EASTING, and ELEVATION × NORTHING.We included the ELEVATION × NORTHING model as a proxy for relative present exposure to ultraviolet-b (UVb) radiation.Other things being equal (e.g., dissolved organic carbon, horizon shading), sites at higher elevation and more southern latitude receive greater doses of UVb [4].Due to minimal variation in FOREST among R. cascadae sites, we converted it to a 1 if the sites' perimeter was >80% forested and a 0 for all other values.
In the final step of our analysis, we compared all models from the 3 spatial scales that had ΔAIC > 2 with each other and with a model that only contained AGE as a site covariate.We included AGE because there is evidence that habitat changes associated with succession can influence pond breeding amphibians (e.g., [8,34]).We used the best ranked model from this step to estimate occupancy.Because the models are logit-linear, we report effect sizes as odds ratios.These were calculated as exp(β) for binomial covariates and as exp(β/SD) for normalized covariates where β is the coefficient for a given covariate and SD is the standard deviation of the covariate prior to normalization.We assessed model fit as indicated by c-hat for the best supported model at each scale and the final models [12].

Bufo Boreas
We surveyed 54 historical B. boreas breeding sites in 1 year (n = 21 sites), 2 years (n = 20), 3 years (n = 12), or 4 years (n = 1).Eleven of the 33 sites (33.0%) that we surveyed in > 2 years had at least 1 year where we detected B. boreas and 1 year where we did not.We detected evidence of B. boreas breeding during at least 1 survey at 41 sites (Figure 1), yielding a naï ve estimate of occupancy of 75.9%.Our estimate of occupancy was 84.9 % (SE = 4.9).The average detection probability given breeding that year was 0.84 (SE = 0.05).We found evidence for effects of YEAR1 and FISH on detectability of B. boreas (Table 2).The odds of detecting B. boreas given breeding that year were lower by a factor of 0.21 (SE = 0.13) in 2001 (the driest year) than in the other years (2002)(2003)(2004).The odds of detecting B. boreas in sites where we detected FISH were 6.8 (SE = 4.7) times the odds of detecting B. boreas in sites where we did not find FISH.Among all predictor variables, AGE and EASTING had the greatest effects on the probability that a site was used for breeding during the study period.For each 10-year increase in the age of the record, holding other factors constant, the odds that B. boreas still used a site for breeding decreased by a factor of 0.61.For every km east within our study area, the odds that B. boreas used a site for breeding increased by a factor of 1.04.

Rana Cascadae
We surveyed 66 historical R. cascadae breeding sites in 1 year (n = 32), 2 years (n = 19), or 3 years (n = 15).Seven of the 34 sites (20.6 %) that we surveyed in >2 years had at least 1 year where we detected R. cascadae breeding and 1 year where we did not.We detected evidence of R. cascadae breeding in at least 1 survey at 44 of the sites (Figure 2), yielding a naï ve estimate of occupancy of 66.6%.The average detection probability given breeding in that year was 0.85 (SE = 0.05).The odds of detecting R. cascadae breeding at sites where they bred and we detected FISH were 0.25 times the odds of detecting R. cascadae breeding at sites where they bred and we did not detect fish.Our final estimate of occupancy adjusted for detectability was 72.4% (SE = 6.6) of historical sites.None of the predictor variables we analyzed were strongly related to occupancy for R. cascadae (Table 3).

Discussion
The western USA has been an area of high concern regarding amphibian declines.Much experimental work has been directed at potential causes.Quantifying the scope and severity of declines can be difficult due to many factors, such as lack of historical information on amphibian occurrence, difficulty accessing sites, and challenges accounting for imperfect detection of target species.Resurveys of historical breeding sites along with use of occupancy models can help address some of these limitations.Our estimates of occupancy among historically used sites (B.boreas 84.9 ± 4.9%; R. cascadae 72.4 ± 6.6%) were higher than in the only study that used comparable methods (B.boreas 49.5 ± 15.4%, R. luteiventris 52.9 ± 9.5%; [11]).Our naï ve estimates of occupancy (75.9% for B. boreas, 66.6% for R. cascadae) exceed rates unadjusted for detectability from historical revisit studies elsewhere in the western USA: B. boreas 17%, R. pipiens 12%, R. sylvatica 69% in Colorado [7]; B. boreas 13% [9], and B. canorus 48%, R. cascadae 34%, R. draytonii 32%, R. muscosa 17% in California [4].These results from around the western USA suggest substantial geographic variation in status among closely related western anurans.Such variation underscores the need to include local and regional factors along with global ones when investigating potential causes of declines.We conducted surveys over multiple years and limited our historical records to breeding sites rather than including records of more mobile post-metamorphic stages.These practices should help provide an accurate and conservative estimate of declines, and may in part explain why our estimates of decline were lower than shorter field studies that included adult records and did not account for imperfect detection.
Our analysis indicates that the probability of B. boreas using a site for breeding decreases with the AGE of the record and increases with EASTING.A variety of evidence indicates B. boreas is well-suited for breeding in disturbed or early successional habitats [23,35,36], so it is possible that this association relates to a decline in breeding site suitability over time.However, the lack of association with proximal FOREST suggests no strong effect of encroaching trees, which we only noted at a minority of the sites.The effect of EASTING on B. boreas occupancy was largely due to our finding no toads in 5 of 6 historical sites in the Klamath Mountains which run from the Cascade Range to the Pacific Coast.The Klamath Mountains have a more complex geology, contain fewer lentic wetlands, and have hotter, drier summers than most of the Cascade Range in Oregon [37].These were also the southernmost records in our study area.
After accounting for detectability, the odds of both B. boreas and R. cascadae breeding at historically used sites were 1.01 times greater for each km north in our study area (Table 4).While evidence for a latitudinal trend in occupancy was modest within our study area, other evidence suggests declines of both species are more pronounced in southern portions of their overall range.For example, R. cascadae appears to be near extirpation at its southern limit around Mt. Lassen [2,17].Some studies refer to R. cascadae declines in Oregon, but field data are not presented [5,14,16,38].The best documented and potentially largest declines of B. boreas are reported from southern portions of its range in California [9,10], New Mexico [39], and Colorado [7,18].Populations farther north in the Rocky Mountains appear to be faring better [30].Both species are common in parts of the Cascade Range in Washington [35,40], and we are not aware of reports of declines in that region.Rana cascadae appears to be relatively common in the two mountain areas that are separated from the main axis of the Cascade Range outside of Oregon: Olympic Mountains in northwest Washington [41] and Trinity Alps in northwestern The presence of trout affected the probability of detection for both amphibians.Our detection of R. cascadae given breeding was negatively affected by fish presence.Larval stages of R. cascadae and most western ranids are palatable to fish and can reduce risky behaviors when fish are present [32,33].Increased use of refuges by R. cascadae larvae with fish could reduce the likelihood of detection.California [26].Alternatively, R. cascadae could breed less regularly at sites with fish, resulting in our detection of the species one year and not another.In contrast to R. cascadae, we found a positive relationship between B. boreas detection and fish presence.A similar positive relationship was reported between toads and trout presence in northern California [26].Most fish find B. boreas unpalatable, but toad larvae are consumed by aquatic insect predators [32,42].Fish may thus benefit B. boreas by reducing their palatable amphibian competitors or predators such as predaceous macroinvertebrates (e.g., [43]).If that results in greater abundance of B. boreas or more conspicuous behavior, toad larvae might be more easily detected at sites that have fish.
We also found evidence of lower detectability of B. boreas breeding in 2001 (a low snow year) than in the subsequent years, which may be due to a lower probability of detection given breeding or a lower probability of breeding in that year compared to other years.The importance of this factor for B. boreas and several natural history observations imply that B. boreas populations have less regular breeding than R. cascadae.First, B. boreas larvae are more conspicuous (e.g., they are black and often aggregate close to shore) than R. cascadae larvae, and thus are less likely to escape visual detection during surveys.Second, B. boreas breeding is likely to be influenced by drought conditions.Adult B. boreas of both genders at intensively monitored sites in Oregon [44] and Colorado [45] are known to skip breeding years and this may be related to snow conditions.In California's Central Valley, B. boreas breeding was more variable than three syntopic amphibians and related to winter precipitation [10].Shallowly sloped littoral habitats preferred by B. boreas for oviposition sometimes do not flood or can strand eggs after winters of low snowmelt ( [46]; CAP, pers.obs.).Some of these observations relate to individuals skipping breeding but our findings suggest that entire populations either did not breed, bred at other sites (e.g., Bufo americanus [8]), or bred but were more difficult to detect in 2001 than in the other years of our study.An alternative but not exclusive explanation is that some populations go extinct and are recolonized within several years [47].Regardless of the mechanism, surveys over multiple years seem advisable when trying to assess site occupancy or population declines in variable breeding species such as B. boreas.
Several aspects of historical site studies temper our conclusions: (1) this type of study is biased toward detecting declines because they do not account for site switching or expansions to nearby ponds [8,11]; (2) we only estimated occupancy, which may not directly address trends in abundance; (3) our historical records were mainly within the last several decades, so we had a limited ability to document older declines; and (4) historical sites do not constitute a probabilistic sample, so our results may not reflect patterns at sites outside our sample.Despite these limitations, our relatively high estimates of occupancy are not consistent with the hypothesis that B. boreas and R. cascadae have experienced broad declines in the Oregon Cascade Range.Further investigation of B. boreas status in the Klamath Mountains of southwestern Oregon is warranted.Other studies with our target species, along with our own finding of a modest latitudinal trend in occupancy, suggest that further study of north-south gradients in declines may help understand broad patterns in amphibian status around the western USA.Imperfect detection for both amphibians in this study argues for the use of occupancy models to estimate occupancy rates and evaluate the effects of covariates.

Table 1 .
Site variables assessed as covariates of detectability or occupancy for B. boreas (N = 54) and R. cascadae (N = 66) historical breeding sites in the Cascade Range, Oregon, USA.Continuous variables are expressed as x (SD) and Range.Categorical variables are expressed as percentages of sites where respective variable was Found or Not Found (binomials) or in each 10-yr class relative to the first year of surveys, 2001 (AGE).
1For B. boreas, Forest is % site shoreline Forested; for R. cascadae, Found is >80% site shoreline Forested and Not Found is < 80% site shoreline Forested.

Table 2 .
Ranking for models of detectability (p) and occupancy (ψ) for historical B. boreas breeding sites (N = 54) in Cascade Range, Oregon, USA.First models are detectability alone; subsequent groups are models of occupancy that incorporate top-ranked detectability covariate (see Methods).

Table 3 .
Ranking for models of detectability (p) and occupancy (ψ) for historical R. cascadae breeding sites (N = 66) in Cascade Range, Oregon, USA.First models are detectability alone; subsequent groups are models of occupancy that incorporate top-ranked detectability covariate (see Methods).

Table 4 .
Percentages of historical breeding sites where B. boreas or R. cascadae were found breeding > 1 year, from 2001 to 2004 (Oregon, USA).Latitudinal range and main river basins are also shown for sample sites in each region (per Figures1 and 2).