When protected areas are not enough: low-traffic roads projected to cause a decline in a northern viper population

Animal mortality resulting from collisions with vehicles has emerged as a major human-caused threat to wildlife. While direct mortality of wildlife from vehicles has been well documented, fewer studies have focussed on the population-level effects of road mortality, particularly due to low-traffic volume roads. We conducted a population viability analysis (PVA) on western rattlesnakes Crotalus oreganus occupying a protected area with low road density and an average traffic volume of ≈350 vehicles d−1, near the northern periphery of the species’ range. We used the program Vortex with a field-derived database on road mortality, population demography, and extent of occurrence. The model showed that although the population had a high likelihood of persistence over the next 100 yr (extinction probability <0.01), a substantial decline was projected (stochastic growth rate −0.035, 97% decrease in mean population size, from 2131 to 72) under the current road mortality rate (6.6% of population yr−1); any increases in road mortality rates were projected to cause extirpation in under 100 yr. Our study provides strong evidence that road mortality is and will continue to be a significant contributor to the decline of this threatened species, even without higher traffic volumes and other significant anthropogenic impacts.


INTRODUCTION
Roads and their associated effects, from mortality to habitat loss, have emerged as a major threat to biological diversity (Forman & Alexander 1998, Forman et al. 2003, Fahrig & Rytwinski 2009, van der Ree et al. 2015. The fragmentation of ecosystems by roads inevitably results in higher-risk movements of animals to reach habitats or other resources. As a result, populations living in proximity to roads often experience higher mortality rates than those in remote areas (Mumme et al. 2000), and direct mortality due to wildlife−vehicle collisions may exceed that experienced by natural causes (Fahrig & Rytwinski 2009), thereby limiting population size and growth (Kushlan 1988, Coffin 2007, Chambers & Bencini 2010. Rapid and widespread declines of various animal populations caused by road mortality (Fahrig et al. 1995, Gibbs & Shriver 2002 may have additional consequences such as genetic isolation (Epps et al. 2005), and in extreme cases, extirpation (Jones 2000). Roads are of concern in protected areas, such as parks, and although they may constitute the only significant source of anthropogenic mortality in these areas (Andrews et al. 2008), the continual and additive mortality from traffic can have devastating consequences on local populations. Populations of atrisk species likely already suffer from depressed numbers, reduced genetic diversity, and multiple threats (Mace et al. 2008), making them even more susceptible to road mortality effects.
Over 10 yr ago, a consortium of road ecologists outlined the need for studies on population-level effects of road mortality (Roedenbeck et al. 2007), yet studies of roadkill impacts on populations remain rare in comparison to more direct assessments of road mortality rates and patterns (Forman et al. 2003, Andrews et al. 2008, van der Ree et al. 2011. Those studies that do exist are often fraught with issues (see examples of problematic studies in Fahrig & Rytwinski 2009). The scarcity of population studies reflects the difficulty in quantifying road effects: to understand fully the impact of road mortality on populations, studies must be intensive, long-term, and accompany assessments of the surrounding wildlife populations rather than just collect evidence from the road itself. Roedenbeck et al. (2007) further exhorted study designs that increase inferential strength; in road ecology, however, manipulative and beforeafter-control-impact studies often are not feasible, especially when assessing impacts of pre-existing roads. Additionally, the urgent need to make management decisions without further compromising the health of populations prevents scientists from acquiring evidence of long-term effects.
In the absence of experimental or long-term data, population viability analysis (PVA) provides a powerful tool for estimating the persistence of populations experiencing road mortality. However, like all modelling exercises, the quality of data being used to establish parameters dictates the usefulness of the outputs (Boyce 1992). If used cautiously (see Beis singer & Westphal 1998 for recommendations), PVA projections of population responses to road mortality and/or mitigation efforts can be highly accurate (Brook et al. 2000) and valuable for conservation practitioners, especially when limited empirical information is available. The flexibility and breadth of PVA allows for assessment of direct and indirect ef fects; for example, incorporating genetic analysis into a model can strengthen the assessment and further illuminate the degree of isolation of populations (Frankham et al. 2014). However, to serve any of these purposes, particularly the assessment of threats from pervasive disturbances like roads, the construction of PVA requires baseline research that generates defendable and accurate parameters to inform the models.
Here, we present a rare case study combining field-collected measurements and PVA to assess the long-term viability of a 'protected' snake population experiencing mortality along a low-traffic road. We coupled our assessments of roadkill rates (that incorporate observer and scavenger corrections; Winton et al. 2018) with intense study of population demography and extent of occurrence on the landscape. A PVA was then constructed using mark−recapture and radio-telemetry data to provide estimates of population size, movements, home range sizes, and annual percent mortality rate due to roads. Increased or decreased rates of mortality, stemming from possible traffic increases or successful mitigation measures, were simulated. Models were also manipulated to assess the sensitivity of the population to variable adult female mortality, decreased life expectancy, and variable initial population sizes. Our subject, the western rattlesnake Crotalus oreganus, is federally listed in Canada as Threatened (COSEWIC 2015); it exhibits traits characteristic of many late-maturing temperate snakes (including the families Viperidae and Colubridae) such as high adult survivorship, prolonged longevity and generation time, and low fecundity and infrequent reproduction (Parker & Plummer 1987, Macartney & Gregory 1988, Macartney et al. 1990). We posit that these characteristics combined with anthropogenic disturbances, such as roads, will limit population growth and make it difficult for northern rattlesnake populations to compensate for continual mortality, thus resulting in population decline.

Study site
Field research was conducted within the White Lake basin (49°N, 119°W) in the South Okanagan region of British Columbia, Canada. This arid ecosystem consists of open shrub-steppe and grassland habitat, predominantly characterized by bluebunch wheatgrass Agropyron spicatum and big sagebrush Artemisia tridentata surrounded by ponderosa pine Pinus ponderosa forest, rolling hills, and steep bluffs (Meidinger & Pojar 1991). Situated near the northern limit of rattlesnakes ( Fig. 1), the White Lake basin contains all critical habitat elements for Crotalus oreganus (COSEWIC 2015), including hibernacula, seasonal migratory corridors, and summer foraging habitat. Twenty-six identified rattlesnake hiber nacula were located within 1.4 km of 2 roads that bisect the basin bottom. Traffic counters showed these roads supported a mean maximum of 350 vehicles d −1 from April to October during the study years − 2017Winton 2018). At the time of our study, there were no road mortality mitigation measures in place at the study site. The rattlesnakes and their habitat within the study site are protected by both federal and provincial designation as part of an 80 km 2 biodiversity preserve. While low-intensity ranching and recreational activities (e.g. hiking) are permitted in some portions of the protected area, land development and natural resource use are restricted. See Winton (2018) for additional description and photos of the site.

Population parameters
We calculated the total area and population of rattlesnakes affected by road mortality in the study site in 2 steps. First, we estimated population size and density at 6 roadside hibernacula using mark−recapture and telemetry data. We then extrapolated the calculated density for the 6 roadside hibernacula to include the extent of occurrence for all 26 known hibernacula within the study site (i.e. the population range size) for an estimate of total population size. Six focal hibernacula within 400 m of the road were sampled (following methodology in Lomas et al. 2015) multiple times (10−30 visits) throughout spring egress and fall ingress of 2015, 2016, and 2017. All captured rattlesnakes were given a unique identifier (passive integrated transponder tag or scale clip), and total captures ranged from 172 to 257 ind. yr −1 (new captures and recaptures). Mark−recapture data collected at the 6 roadside hibernacula were used to estimate population size using the Jolly-Seber model in the 'Rcapture' package (Baillargeon & Rivest 2007) for R (version 3.4.2; R Core Team 2017). Each of the 3 years of sampling was considered a capture period, and the capture histories of each individual snake were used across all 3 periods. The fit of the model was determined by examining the Pearson residuals, all of which were <1 in our model, thus indicating high quality of fit.
To determine the portion of the landscape that supported this population (and thus permit density calculations), we measured the extent of summer movements away from the hibernacula by conducting telemetry on 27 rattlesnakes (Table S1 in the Supplement at www. int-res. com/ articles/ suppl/ n041 p131_ supp. pdf). In addition, we incorporated all locations of rattlesnakes (n = 25) marked at any of the 6 focal hibernacula and subsequently en countered over the summer months. We calculated a minimumconvex polygon (MCP) area around each hibernaculum using the aggregated snake movements and locations in Garmin BaseCamp (version 4.6.2, 2016), based on the recommendation by Row & Blouin-Demers (2006). We then calculated rattlesnake density by dividing population size by MCP area used.
To estimate the population range, we determined the maximum distance travelled from the respective hibernacula by all telemetered rattlesnakes (including those from non-focal hibernacula, n = 32 snakes; Table S1). We then followed a process similar to Weir et al. (2011) to determine the extent of occurrence, by using the mean movement distance from the hibernaculum (1.3 km) and applying buffers of this radius centred around all known hibernacula in the study site (n = 26) using ArcGIS (version 10.2.2, 2014). If an MCP area for a hibernaculum with tracked rattlesnakes exceeded the buffer area, the MCP was used. The resulting boundary of the hibernaculum buffers and MCPs established the population range. With 1 exception, the population range encompassed all locations of rattlesnake encounters (n = 176) over the summer months, including road mortalities, and contained 13.3 km of road. We used the mean road mortality rate for rattlesnakes calculated for our population (0.058 km −1 d −1 ; Winton et al. 2018) to determine the percentage of the population lost annually. The rattlesnake density previously calculated for the 6 roadside hibernacula was then applied to the entire population range for an estimate of the total rattlesnake population size. This method of estimating population size may inflate the population range area and therefore population size, thus resulting in a relatively conservative assessment of the probability of extinction of the population.

Population viability analysis
We modelled the population and assessed persistence using VORTEX (version 10.2.7.0), a program suitable for long-lived animals with low reproductive rates, and for the application of road mortality data (Lacy 1993). Population models were simulated 500 times and did not include dispersal or inbreeding effects. A quasi-extinction threshold of 10 individuals (irrespective of sex) was set, beyond which other rattlesnake populations have been considered no longer viable (COSEWIC 2012); hereafter 'extinction probability' denotes the probability of reaching this quasi-extinction value. We classified input parameters for the population model into 2 categories: population-specific parameters (as detailed in Section 2.2) and species-specific parameters. Species parameters were derived from published research pertaining to life history and demographic traits of C. oreganus in the northern portion of their range, with emphasis on previous studies involving populations in British Columbia (Table 1).
Road mortality was modelled within the PVA as the annual removal of a constant percent of individuals from the population, after the adults had bred but before they had aged; thus, road mortality constituted an additive source of mortality to natural mortality. In addition to PVA models examining the effect of variable road mortality rates on population persistence, we also simulated the impacts of varying longevity and maximum and minimum adult female road mortality rates. Proportions of the total rattlesnake roadkill for sex and age classes were determined from specimens collected during road surveys ) and incorporated into the models (Fig. S1). Our PVA models were used to determine the probability of extinction, stochastic population growth rate, and mean size of extant populations.

Population size, density, range, and road mortality rate
The estimated population size of the 6 focal hibernacula (combined) from 2015−2017 was 452 ± 59 (mean ± SE) rattlesnakes. These rattlesnakes occupied an estimated area of 8.35 km 2 around the hiber-  Table 1. Species-specific parameters and sources used for population viability analysis of western rattlesnakes Crotalus oreganus in the White Lake basin in British Columbia (BC), Canada. All estimates were taken from other studies conducted on the species within the same region (Okanagan Valley, BC) nacula ( Fig. 2), yielding a density of 54.1 rattlesnakes km −2 . The larger population range area occupied by rattle snakes from all 26 known hibernacula was 39.4 km 2 (Fig. 2), and after extrapolation from the density estimate, this area was considered to support a population of 2131 ± 279 (SE) rattlesnakes. Based on the mean annual road mortality rate of 0.058 km −1 d −1 (measured value; Winton et al. 2018), 141 rattlesnake deaths occurred each year along the portion of road encompassed by the population range (13.3 km), equivalent to 6.6% of the population. Adult female road mortality ranged from 8 to 17% of total roadkilled rattlesnakes.

Population viability analysis
Using the measured road mortality rate of 6.6%, our estimated probability of extinction for this population was < 0.01 in 100 yr and 0.0 in 50 yr. We then modelled extinction probabilities for annual road mortality rates of 0 to 20% (2% increments) of the population (Fig. 3). Above a simulated road mortality rate of 6%, the probability of extinction in 100 yr increased directly with road mortality up to a rate of 14%, at which point the probability of extinction reached 1.0. For a forecast period of 50 yr, simulated road mortality rates of ≤14% had zero probability of extinction with a similar increasing trend in extinction probability observed for higher road mortality rates (Fig. 3).
The stochastic population growth rate (r) was consistently negative through the range of simulated road mortality rates, indicating that the population would decline over time ( Table 2). Even though the population was projected to persist for the next 100 yr (0.4% probability of extinction), the PVA projected a constant decrease in population size. For example, at the current road mortality rate of 6.6%, the stochastic growth rate was −0.035, and the mean population size was estimated to decrease to 72 individuals (97% de crease in 100 yr). Even at a low simulated road mortality rate of 2%, the population would experience a drop of 54% in 100 yr. In contrast, the growth rate in the absence of road mortality was 0.005 and the population was projected to increase by 60% over 100 yr (Table 2). 135 Fig. 2. Population range of the western rattlesnake Crotalus oreganus in the White Lake basin, British Columbia, Canada. Minimum convex polygons (MCPs) of the focal hibernacula (grey) were used to determine the population density, while 1.3 km buffers centred over hibernacula in combination with MCPs were used to establish the population range area. Field-collected data used to demarcate these areas included rattlesnake hibernacula locations, encounter locations, and telemetry movements. Within the area shown on the map, 26 rattlesnake hibernacula occur within 0.0−1.4 km from the roads; the precise locations of these sites are omitted due to the sensitive nature of the data

Adult female mortality
A higher proportion of adult females dying on the road (with total adult mortality held constant) resulted in a higher probability of extinction for simulated road mortality rates between 6 and 14% (Fig. 4) as well as consistently lower population growth rates for all simulated road mortality rates (with greater differences between growth rates at higher road mortality rates). When we used the measured rate of road mortality (6.6%), the modelled probability of extinction over 100 yr shifted from 0.0 to 0.1, depending on the proportion of adult females dying (minimum 8% of roadkill to maximum 17%).

Longevity
The probability of population persistence increased when we adjusted rattlesnake longevity upwards (Fig. 5). At a simulated level of low road mortality (4%), and with snakes reaching an age of 14 yr, the population had a high probability of persistence (1.0) over 100 yr; however, this probability dropped to 0.71 when the road mortality rate was moved upwards to 6% and dropped even further to 0.23 at a road mortality rate of 8%. In fact, at a road mortality rate of 8%, even if snakes were simulated at a maximum age of 21 yr (extremely high, and improbable for a population in proximity to roads; J. Petersen et al. unpubl. data), the probability of persistence was less than 1.0. If the snakes were simulated to live to only 10 yr of age, the probability of persistence for the population approached zero (< 0.02), even under a simulated rate of low road mortality (4%), well below the measured rate of 6.6%.

DISCUSSION
Our PVA strongly suggests that the study population of snakes, residing in a wildlife reserve with minimal development or disturbance, and a traffic volume of only 350 vehicles d −1 (Winton 2018) Table 2. Stochastic growth rate, mean population size, and mean time to extinction for the population of western rattlesnakes Crotalus oreganus in the White Lake basin, British Columbia, Canada, over 100 yr at different rates of road mortality. Extinction was considered when the population reached a quasi-extinction threshold of 10 individuals. clining due to road mortality. Consequently, the impact of snake mortality from vehicle collisions is considerable at the population level, and particularly disconcerting if the population faces an increase in vehicular traffic. Our analysis is built upon reasonably extensive empirical data, along with adjusted road mortality rates calculated from intensive surveying and corrections for bias ). Thus, the foundation for the PVA is relatively robust, allowing for meaningful projections, although given that our assessment of the current population size is likely an overestimate, the situation is concerning.
In most of our scenarios, a simulated road mortality rate of 6% appeared to be the maximum that the population could tolerate without an appreciable risk of extinction in 100 yr. Although this extinction probability was not high for simulated road mortality rates ≤8%, the rate of population decline was severe. In fact, growth rates were negative for all simulated road mortality rates, including those well below the measured value (6.6%) during our study. Therefore, road mortality will reduce the population to a size that, in addition to loss of genetic diversity, could be sensitive to other anthropogenic stressors, disease outbreaks, or catastrophes that cause further population declines. Our results highlight the difficulty of maintaining a viable population of this species over time with any degree of road mortality and identify a critical need to reduce this anthropogenic source of mortality.
Projections of population decline and extirpation in our PVA models showed high sensitivity to the proportion of adult females killed on the road, as well as the longevity of individuals in the population; this is not surprising given the importance of breeding females to the persistence of wildlife populations (Hebblewhite et al. 2003, Ramp & Ben-Ami 2006, especially long-lived reptile species (e.g. Congdon et al. 2001, Row et al. 2007, Bulté et al. 2010). In our simulations, the probability of population persistence to 100 yr was always <1.0 with road mortality rates set at ≥8%, even if these simulations allowed snakes to live to extreme ages. Our results strengthen the argument that road mortality rates > 6% pose a serious risk to the viability of the population by making it less likely that individuals will survive to enter older age classes. Taken in tandem, the observations that populations are sensitive to the loss of reproductive females and that roads likely reduce life expectancy are concerning, and exemplify the compounding effects that road mortality can have on a population through the continuous removal of reproductively mature individuals, especially when natural limiting factors also are at play. Many reptiles exhibit reproductive strategies similar to Crotalus oreganus where the loss of adults is a significant threat to populations (Brooks et al. 1991, Enneson & Litzgus 2008, Hyslop et al. 2012, Brehme et al. 2018. Generally, populations of long-lived animals with low reproductive rates that experience roadkill or any other source of additive mortality tend to lack compensatory mechanisms (Congdon et al. 1993, Webb et al. 2002, Andrews et al. 2015. While road mortality often has detrimental effects on the abundance of wide-ranging mammals (e.g. Hoodicoff et al. 2009), our study and several others (Gibbs & Shriver 2002, Row et al. 2007 show that at a finer scale this negative impact is also true for reptiles with smaller home ranges. As intact habitat patches become smaller (Ibisch et al. 2016) and roads separate habitats and/or resources, more local populations will be threatened by roadkill, and conservation actions will be required to prevent extirpation and maintain contiguous species ranges.
Our study demonstrates that habitat protection (i.e. wildlife reserves) alone may be insufficient to prevent a continued decline in populations experiencing road mortality, even under relatively low traffic volumes. In turn this suggests that a reduction in road mortality rates is paramount to sustaining the population. Although a complete elimination of road mortality would go the furthest in preventing population decline and enabling recovery, such a target is not feasible without a major decrease in vehicle traffic and/ or highly effective mitigation measures. Numerous mitigation measures including fences paired with underpasses are being studied and used to successfully reduce mortality for various wildlife species (Rytwinski et al. 2016, Colley et al. 2017, but see Baxter-Gilbert et al. 2015. Other suggested methods (less rigorously studied) for reducing reptile road mortality include reduced speed limits (Farmer & Brooks 2012), dynamic signs (Hardy et al. 2006), population supplementation or translocation (Madsen et al. 1999), and road closures (Kline & Swann 1998, Palis 2016.
The results of our PVA combined with movement data from radio-telemetry (S. A. Winton unpubl. data) further suggest that protection/enhancement of connective and other critical habitat would be beneficial, particularly for reproductively mature females (e.g. breeding or gestation sites; D. Eye unpubl. data) living away from roads. The protection of migration corridors across roads paired with the establishment of breeding or refugia sites in proximity to core habitat could improve rattlesnake survival while enhancing reproductive success. Intuitively, the expansion of existing roads or other activities that would en-courage increased traffic in areas where populations are already impacted by road mortality should be avoided, as wider roads and increased traffic volumes are associated with higher road mortality rates and increased barrier effects (Forman & Alexander 1998, Valero et al. 2015. Our analysis provides a precise and detailed understanding of population-level consequences of road mortality on a northern snake population; this can be used to assess, through continued monitoring of populations and roadkill rates, the efficacy of conservation efforts.