Patterns of red tree vole distribution and habitat suitability: implications for surveys and conservation planning

Environmental regulations often require wildlife surveys prior to habitat disturbance to avoid impacts or as the basis for planning mitigation, yet project-level surveys may not provide the insights needed to guide long-term management. Management of the red tree vole (Arborimus longicaudus) has largely been based on such surveys. As an alternative approach, we evaluated distribution patterns using frequency of red tree vole occurrence and habitat suitability models to guide conservation planning. We developed a suite of models based on subsets of covariates from two previously developed models and evaluated the extent to which spatial covariates improved the models. We used presence–absence data that were collected from 364 randomly selected 1-ha Current Vegetation Survey and Forest Inventory and Analysis plots to develop models and describe occurrence patterns. The best models included a spatial covariate, maximum tree diameter, distance from suitable habitat, forest age class, and the interaction between maximum tree diameter and forest age class. We compared performance of the previously published models, our best model, and an ensemble model that used predictions from all three models. Under the ensemble model, correct classification rates were relatively high and considerably improved, suggesting that the application of all three models provided greater accuracy than any individual model. We argue that habitat models, coupled with spatial patterns of the frequency of occurrence, can provide useful tools for addressing species management and may provide more insight than project-level surveys. The use of habitat suitability models can therefore be closely tied to red tree vole management decisions and conservation strategies, as well as reducing survey costs that otherwise often make projects infeasible.


INTRODUCTION
Plant and animal surveys are often conducted to estimate trends in abundance (McComb et al. 2010), to evaluate effects of management (Nichols and Williams 2006), to develop habitat models (Scott et al. 2002), and to document species presence to facilitate mitigation (USDA and USDI 2001) and conservation design (Groves 2003). Many federal and state agencies require field surveys prior to habitat disturbance to identify locations of plants or animals of conservation concern to meet policy requirements. However, in many of these cases, it is not clear how projectlevel surveys can guide management decisions because the data do not answer the questions they are intended to address. In addition, even if surveys can address the management issues, the spatial scale of sampling is often inadequate because of costs, especially for species that are difficult to detect.
Habitat suitability models have the potential to help guide the need for surveys that are intended to document species presence to avoid harm, aid in developing conservation plans, and address mitigation by identifying potential habitat that may be disturbed, protected, or restored (Dunk et al. 2004, USFWS 2012. Habitat suitability models should also be able to help managers decide if and where surveys are most needed to address management questions and provide a broader perspective for species management than surveys alone. The use of habitat suitability models can therefore be closely tied to survey requirements that are part of many mitigation strategies. Because suitable habitat is not always occupied (Capen et al. 1986, Hanski 1999, Fielding 2002, predictions of habitat suitability from models may be more useful than the recognition of whether particular locations are occupied or not at a single point in time. Identifying suitable but unoccupied habitat is critical for conservation (Camaclang et al. 2015). If surveys fail to detect the target species, the area may no longer be considered of conservation value and habitat modification may be allowed. This is especially important for wide-ranging species that are unlikely to be present in a small portion (e.g., a small sampling unit) of their home range, even in highquality habitat. Removing or modifying suitable but unoccupied habitat eventually results in population declines. Such understanding is the cornerstone of metapopulation dynamics but can be operative even if populations are not structured as such (Hanski 1999:158). More effective conservation strategies may be best served when long-term occupancy or species persistence is addressed rather than the short-term view obtained by surveys conducted at a single point in time (Camaclang et al. 2015).
Wildlife surveys are a key requirement of mitigation approaches in the Northwest Forest Plan (NWFP). The NWFP provides one of the most comprehensive strategies to protect plant and animal species in the world (USDA and USDI 1994a, DellaSala and Williams 2006, Carroll et al. 2010a, and was conceived as an ecosystem management plan to protect species associated with old forests on 9.8 million ha of lands managed by the Bureau of Land Management and U.S. Forest Service in northwestern California, western Oregon, and western Washington (USDA and USDI 1994a).
The NWFP was a compromise between protecting habitat for species associated with old forests and timber production. New reserves created by the NWFP brought the protected area of federal land to approximately 80% of the total federal forest land within the planning area and included approximately 86% of the remaining old forest (USDA and USDI 2000a, Thomas et al. 2006). The reserve strategy was thought to ensure the viability of most old-forest-associated species, but for some species, there was concern for their continued persistence (USDA andUSDI 1994b, Molina et al. 2006). The Forest Service and Bureau of Land Management implemented an extensive program of surveys in which they attempted to better document the distribution and habitat associations of over 400 species that were thought to be vulnerable (Molina et al. 2006, Olson et al. 2007). These surveys, collectively referred to as the Survey and Manage Standards and Guidelines (USDA and USDI 2001;hereafter SM Program), were intended as a "fine filter" approach for managing specific species to augment the "coarse filter" strategy of reserve designation (USDA and USDI 1994b, Edwards et al. 2004, Molina et al. 2006. For some of the rarer species, the SM Program may have achieved its mitigation goals: Known sites were protected and conservation planning for these species did not greatly affect timber harvest (Molina et al. 2006). The greatest controversy arose over species that were located frequently during surveys in areas proposed for timber harvest (Molina et al. 2006, Thomas et al. 2006. The red tree vole (Arborimus longicaudus), a small arboreal mammal endemic to coniferous forests in western Oregon and northwest California, has been the most controversial species included in the SM Program and provides an excellent case study in the use of habitat suitability models in conservation planning. Federal agencies are required to avoid harvest within 4 ha of occupied or recently occupied tree vole nests in old forest (USDA and USDI 2001). In lieu of this approach, federal agencies can conduct site assessments to determine whether a site is non-high priority to the species' persistence, or develop conservation strategies at a watershed scale (Huff 2016). Surveys to find nests and document occupancy by climbing trees are timeconsuming and expensive, detection rates are extremely low (Marks-Fife 2016), and where the species is frequently found, the no-harvest buffers may conflict with timber harvest or forest restoration goals. Because of high cost, surveys conducted prior to management activities (hereafter, "pre-disturbance surveys") occur only within project areas (USDA and USDI 2000b) and thus often cannot provide the broader landscape perspective that is desired for addressing species persistence mandates. These characteristics have made tree voles difficult to manage and have led to lawsuits regarding tree vole management in proposed federal timber sales.
Although surveys can address the management goals of providing disturbance buffers surrounding occupied sites, and in part, preparation of assessments as non-high-priority or high-priority sites (Huff 2016), they are often neither effective nor efficient. Because tree vole nests are extremely difficult to detect, identification of areas in which to place buffers is incomplete. Similarly, because surveys underestimate the number of voles, surveys fail to fully document areas where the species is not in need of management focus. More importantly, project-scale surveys may not address patterns of long-term occupancy. Finally, surveys alone may not provide sufficient information because conservation planning for tree voles is typically conducted at the fifth-field watershed scale, a scale too cost prohibitive for surveys. Rather than relying entirely on surveys, incorporating data on geographic patterns of occurrence along with habitat suitability models provides a scientifically credible approach for addressing tree vole management and conservation of species on managed landscapes more generally.
We developed and evaluated habitat suitability models and identified patterns of geographic occurrence to inform management of tree voles. We discuss the efficacy of using these tools to inform management to an extent that projectlevel surveys by themselves do not. If models perform sufficiently well, their use could provide a cost-effective method of predicting potential site occupancy by tree voles and informing conservation strategies. Similarly, patterns of occurrence and habitat models could inform managers of the likelihood of tree voles being present at areas proposed for disturbance. Nonetheless, such ideas are not currently institutionalized and implemented by land management agencies.
There have been three previous large-scale models developed to estimate habitat suitability or habitat associations for tree voles. One of these used presence-only data (Forsman et al. 2016), the second used logistic regression under a generalized additive modeling approach (Dunk and Hawley 2009), and the third used a Bayesian approach to add spatial patterning as a random effect to the Dunk and Hawley model (Carroll et al. 2010b). We did not explore Carroll et al.'s (2010b) model because we modeled spatial pattern as a fixed effect. We developed a suite of models based on covariates from the first two models and then evaluated the extent to which additional spatial covariates improved models. Because of the high spatial variation of red tree vole abundance (Forsman et al. 2004(Forsman et al. , 2016, we hypothesized that models with spatial covariates would be better predictors of occurrence and that these covariates would allow greater insight for management. Our goals were to (1) evaluate the frequency of tree vole occurrence among different spatial scales and geographic locations to identify geographic patterns of occurrence that can be used in conjunction with habitat models to inform management and (2) develop models of tree vole occurrence based on vegetative and spatial variables to reliably predict the presence or absence of tree voles and potentially serve as a surrogate to field surveys in developing conservation strategies for tree voles.
Columbia River to the Klamath River in northern California, with the exclusion of the Willamette Valley (Fig. 1). The study area was dominated by mountainous terrain and included a diversity of forest types. In the northern and central portion, coniferous forests of Douglas-fir (Pseudotsuga menziesii) and western hemlock (Tsuga heterophylla) dominated, whereas in the southern portion of the study area, drier conditions favored mixed conifer-hardwood forests (Franklin andDryness 1973, Forsman et al. 2016

Red tree vole surveys
Our models, those from Hawley (2009) and, in part, Forsman et al. (2016), were fit with data from a stratified random survey conducted on federal lands in 2001(Rittenhouse et al. 2002. We used the same plots as Dunk and Hawley (2009), except we omitted one plot because geographic information systems (GIS) data were unavailable for the appropriate time frame. This resulted in 364 1-ha plots subsampled from the grid of 1003 Current Vegetation Survey and Forest Inventory and Analysis plots located in the study area on lands managed by the Forest Service and Bureau of Land Management ( Fig. 1; Rittenhouse et al. 2002). Plot selection was stratified by forest age class and land-use allocation. Age stratification was based on two age classes: (1) old forest (>80 years) and (2) young forest (≤80 years), the transition age when old forest characteristics generally begin to develop (Old-Growth Definition Task Group 1986, Thomas et al. 2006). Land-use allocation strata were based on two groupings: (1) reserved lands where management is focused on maintenance and restoration of old forests (Reserve) and (2) non-reserved lands where timber harvest is emphasized (Matrix; USDA and USDI 1994a). The sample was stratified as 60% reserve/old forest; 20% reserve/young forest; 10% matrix/old forest; and 10% matrix/young forests, but was slightly modified because some plots were removed due to logistical and safety issues (Dunk and Hawley 2009). Mean age of trees in young-forest plots ranged from 0 (recent harvest) to 80 years old, whereas trees in old-forest plots ranged from 81 to 384 years old.
Tree vole surveys were conducted by at least two observers who searched for potential nests along four 100-m transects spaced 25 m apart within a 1-ha square plot (Rittenhouse et al. 2002). All trees observed with potential nests were climbed, and each nest was searched for evidence of tree voles, consisting of conifer cuttings, resin ducts, fecal pellets, and debarked twigs. If nests were not found during ground surveys or in plots with complex canopies, tree climbers searched for nests in five randomly chosen Douglas-fir trees that were >61 cm diameter at breast height (dbh), or the largest trees if none met this size criterion (Rittenhouse et al. 2002). Tree vole nests were classified as "active" if they contained evidence of recent occupancy or "inactive" if nest material was old (Rittenhouse et al. 2002).
We defined tree vole presence as a plot with ≥1 active or inactive nest observed, consistent with Hawley (2009) andForsman et al. (2016) because tree vole nests only persist for a couple of years (Thompson and Diller 2002). Using nests as a surrogate for tree voles probably caused us to overestimate the proportion of plots occupied by voles, whereas imperfect survey detection resulted in an underestimate. Dunk and Hawley (2009) estimated that 6% of the plots had false negatives based on finding nests from the random tree searches in plots in which nests were not otherwise detected. Our analyses assume absence of tree voles from plots where they were not detected.

Model development
We extended models by Forsman et al. (2016) and Dunk and Hawley (2009) by adding a suite of spatially explicit predictor variables (Fig. 2). We included spatial covariates related to geographic regions and to the distribution of suitable habitat. Suitable habitat was estimated from the previously published models and this required extrapolation of the Dunk and Hawley (2009) model to the entire study area, whereas the Forsman et al.'s (2016) developed estimated suitability throughout the study area (Fig. 2).
Previously published models.- Dunk and Hawley (2009) fit generalized additive models with the random plot data. We included only their best model, which we refer to as DH-GAM. Explanatory variables in DH-GAM included basal area of trees 45-90 cm dbh, maximum dbh, standard deviation of dbh for conifer trees >2.5 cm dbh, ❖ www.esajournals.org slope, and the Universal Transect Mercator (UTM) coordinates of the center of the plot. Forsman et al. (2016) created their models with program MaxEnt (Phillips et al. 2006) using presence-only data from random plots, predisturbance surveys, and other location data. Their model included forest structure and tree species composition covariates from the Gradient Nearest Neighbor (GNN) vegetation database and models Gregory 2002, Moeur et al. 2011). GNN data linked ground-based data with remote sensing data (Moeur et al. 2011) and were based on the average 30 9 30 m pixel GNN predictor variables resampled to 1 ha. Forest structure and species composition covariates that were most influential in their models included (1) large trees per ha (≥75 cm dbh); (2) percentage of hardwood cover; (3) diameter diversity index, a measure of canopy structural diversity (Ohmann et al. 2012); (4) percentage of conifer cover; (5) percentage of total basal area of tree species groups based on region; and (6) a tree vole food source variable estimated as the percentage of total basal area of Douglas-fir, western hemlock, Sitka spruce, and grand fir (Abies grandis). The most influential abiotic covariate was the frequency of summer fog. Forsman et al. (2016) constructed separate models for each of the four  (Forsman et al. 2016). We incorporated environmental and habitat suitability-related covariates from these models into our generalized linear models (GLMs). After evaluation of the models with red tree vole survey data from the random plot and the pre-disturbance survey areas, we developed the ensemble model that consisted of only consistent predictions from the three models as suitable or non-suitable habitat.
GNN modeling regions, which were delineated based on physiographic province and ecoregions. They spliced together their regional models to provide a single model to estimate habitat suitability for the entire study area. It is this model that we refer to as DA-MAX and from which comparisons were made to the other models.
GLMs.-We estimated the probability that tree voles occurred within a random plot using logistic regression with generalized linear models (GLMs) in PROC GENMOD (SAS 2000) and used Akaike's information criterion (AIC, Akaike 1973) to rank models (Burnham and Anderson 2002). We developed a set of models within 10 unique themes and a set of models incorporating multiple themes ("hybrid models"; Table 1; Appendices S1-S11). We included linear, logistic, and quadratic relationships of the probability of tree vole occurrence based upon hypothesized relationships with the covariates. We then developed hybrid models that incorporated covariates from the best models of one or more themes. To the best hybrid models, we added forest age class and land-use allocation to account for the stratified sampling design and included the interaction term between maximum tree diameter and forest age class. The interaction term evaluated the hypothesis that in young forest large trees have a greater effect on tree vole occurrence than in old forests (Thompson and Diller 2002). We evaluated a total of 118 unique models with AIC. To illustrate the estimated effects from our best model, we created histograms showing the number of plots with and without tree voles in relation to the parameter value (Smart et al. 2004).
1. Covariates from previous models-We used environmental covariates from DH-GAM and DA-MAX in constructing our GLMs. We did not include UTM plot locations because we could not replicate the DH-GAM approach of modeling plot locations into our GLM (Yee and Mitchell 1991, Table 1. Highest ranking logistic regression models of the probability of red tree vole occurrence in 1-ha random plots within each set of related models (theme) and among all 118 models. Notes: K is the number of parameters in the model, ÀL is thelog likelihood, AIC is Akaike's Information Criterion, W t is the AIC weight for the model within its own theme, ΔAIC all is the difference of the model's AIC from the model with the lowest AIC from the comparison using all models, and W all is the associated AIC weight. Hybrid models were developed by combining models from more than one theme. Models are ranked from "best" to "worst" based on AIC.
† Sets of models were constructed for covariates related to suitability of areas surrounding 1-ha plots (Surrounding Hab), distance to suitable habitat (Dist to Hab), environmental covariates included in DH-GAM (Dunk and Hawley, Dunk and Hawley 2009), area of suitable habitat in which a plot was embedded (Patch Size), covariates included in DA-MAX (Davis and Andrews, Forsman et al. 2016), distance of plots to the center or periphery of the study area, or high-density areas (Distance Metrics), geographic boundaries (Boundaries), mean habitat suitability within fifth-field watersheds or sub-basins (Watersheds Hab), forest age (AGE), covariates related to selection of plots (Stratification), and more than one theme (Hybrid).
‡ Covariates in the models include age class of forest (AGECLASS), basal area of trees with a diameter between 45 and 90 cm (BA4590), density contours of tree vole occurrence (CON), canopy cover of dominant conifer trees (CONCOV), diameter diversity index (DDI), frequency of summer fog (FOG), percentage of total basal area of conifer species fed upon by red tree voles (FOOD), maximum tree diameter (MAXDBH), standard deviation of the diameter of conifer trees (SDCONdbh), slope (SLP), and log-transformed values of mean habitat suitability estimated from DH-GAM at the scale of the sub-basin (L_basin-hab_GAM), mean habitat suitability estimated from DH-GAM of 4 ha surrounding each plot (L_buff_4 ha_GAM), distance from the center of the study area (L_meancenter), distance from suitable habitat estimated from DH-GAM with a threshold of 0.50 (L_dist_GAM_050), the number of trees >75 ch dbh per ha (L_LTPH), and the patch size of surrounding suitable habitat estimated from DH-GAM (L_size_GAM_025). Guisan et al. 2002). We included forest age as a continuous covariate and categorized as young or old. Age was included in some of the Dunk and Hawley (2009) models, but not their best model. From DA-MAX, we considered all of the covariates that were most influential in their model (Appendix S2).
2. Spatial variation-We identified areas where the density of plots with tree vole nests varied the most using a fixed kernel density estimator to delineate five contours (10, 20, 50, 80, and ≥95%) of density (Worton 1989, Beyer 2012). These contours were the smallest areas containing the stated level of the probability distribution. Each inner contour (e.g., 10%) contained the highest density of plots for the given area. We specified a bandwidth smoothing parameter using smoothed cross-validation and a cell size of 1 km for estimating the kernel density function, which we estimated in Geospatial Modeling Environment (Beyer 2012). The contours were included as a categorical covariate in the models.
We also developed covariates using a variety of geographic zones: (1) physiographic provinces  Forsman et al. (2004Forsman et al. ( , 2016, but the few plots east of the Cascade Crest were included in the North Cascades subregion (Fig. 3c), (4) survey zones which Huff et al. (2012) delineated as priority zones for conducting surveys to which we added a zone outside of where surveys would be required (Fig. 3d), (5) federal forest administrative units ( Fig. 4), (6) reserved or non-reserved land-use allocations (USDA and USDI 1994a), (7) designated critical and non-critical habitat for northern spotted owls (Strix occidentalis caurina; USFWS 2012) under the hypothesis that spotted owl critical habitat provides "umbrella" protection for tree voles, (8) delineation of the Coast Ranges north of the Siuslaw River where a distinct population segment for the red tree vole was identified (USFWS 2011), (9) sub-basins, and (10) fifth-field watersheds.
We also evaluated hypotheses that habitat suitability of a species decreases (1) toward the edge of the range, (2) away from areas of high densities, (3) with increasing distance from the center of the range (Brown 1995, Osborne andSu arez-Seoane 2002), and (4) with distance from suitable habitat. For distance from the perimeter of the range, we computed the distance from the edge of each plot to the nearest non-coastal perimeter of the study area. To evaluate proximity to an area of high density, we computed the distance to the center of the 10% contour defined by the kernel density estimator. To evaluate the association of tree vole occurrence with distance to the center of the species' range, we computed the distance from the spatial center of the study area (ArcGIS 10.0 [Environmental Systems Research Institute, Redlands, California, USA]: Mean Center tool).
For distance to suitable habitat, we developed several covariates based on the suitability of habitat as predicted from DH-GAM and DA-MAX. For DH-GAM, we first had to extrapolate the model algorithm onto GNN data for coverage of the entire study area (Appendix S12). The covariates we explored were as follows: (1) DH-GAM with suitability defined as predicted probability of occurrence ≥0.25; (2) DH-GAM with a predicted probability ≥0.50; (3) DA-MAX with a relative habitat suitability ≥0.34; and (4) the same covariates as above but with 4 ha of suitable habitat surrounding the plot, equivalent in size to mitigation buffers (USDA and USDI 2001). We used these threshold values that denoted habitat as suitable based on their previous use in Dunk and Hawley 2009), the often used midpoint, 0.50 (Hosmer and Lemeshow 1989), and the average threshold among modeling regions for DA-MAX (Forsman et al. 2016).
We examined a suite of covariates related to the amount of continuous suitable habitat in which the plot was embedded (patch size) and the average suitability surrounding plots at multiple spatial scales (neighborhood effects). For patch size, we defined three different criteria of what constituted suitable habitat based on the habitat model (DH-GAM or DA-MAX) and suitability threshold (DH-GAM: 0.25, 0.50; DA-MAX: 0.34). For neighborhood effects, we used DH-GAM and DA-MAX to estimate mean habitat suitability within circles of 4 ha and 40 ha from the center of the plot, and within fifth-field watersheds and sub-basins.

Extrapolation of model to the study area
We were interested in (1) estimating habitat suitability covariates based on DH-GAM and ❖ www.esajournals.org DA-MAX in our GLMs, (2) evaluating models with pre-disturbance survey data, and (3) providing a model to estimate habitat suitability throughout the study area. Estimates of habitat suitability from our GLMs and those from DH-GAM were limited to the random plots because one or more of the covariates were derived from plot-based measurements. Therefore, we extrapolated our best model and DH-GAM onto GNN data (http://lemma.forestry.oregonstate. edu/data) for the entire study area. We used Marine Geospatial Ecology Tools version 0.8a50 (MGET; http://mgel.env.duke.edu/mget; Roberts et al. 2010) and ArcGIS 10.0 to estimate probabilities of tree vole occurrence for each 30 9 30 m pixel within the study area (Appendix S12). DA-MAX modeled habitat suitability across the entire range, so extrapolation was unnecessary.

Performance of models
We explored the performance of DH-GAM, DA-MAX, and our best GLMs with correct Fig. 3. Percentage of plots with red tree vole nests detected among (a) physiographic provinces, (b) modeling regions described by Forsman et al. (2016), (c) subregions as depicted by Forsman et al. (2004Forsman et al. ( , 2016, and (d) survey zones described by Huff et al. (2012). We allocated the southwestern-most plots in California to Klamath Mountains Province (a) and the California plots to the South Coast or Interior Southwest subregions (c). We pooled the five survey plots that were within the northwest corner of the East Cascades subregion (circled) with the North Cascades subregion because of the similarity of vegetation at these plots (a and c). Plots with voles detected (triangle) and not detected (circle) are shown. classification rates of both absence and presence. For our GLMs, we used a threshold value of 0.25 because this threshold resulted in >70% correct classification for plots with and without tree voles. We used the same value for DH-GAM, consistent with Dunk and Hawley (2009), and 0.34 for DA-MAX, the mean of the thresholds from their four modeling regions. We interpreted estimates of probability of occurrence (our GLMs and DH-GAM) or relative habitat suitability (DA-MAX) that were greater than or equal to these thresholds as a prediction of "presence," whereas we interpreted estimates less than the thresholds as "absence." We then evaluated models using Cohen's kappa, a measure of the proportional improvement over chance for correct classification. Cohen's kappa ranges from 0 (agreement no better than random assignment) to 1 (perfect agreement ;Agresti 1990:366). In addition, we evaluated the relationship, via linear regression, between the proportion of sites with tree voles detected within each 0.05 bin of predicted probabilities (20 bins) to estimates of habitat suitability from DH-GAM, DA-MAX, and our best GLM, with the expectation that well-performing habitat models should have a positive correlation.
We evaluated model performance using two data sets. We used the random plot data that we used to develop our GLMs, DH-GAM, and in part, DA-MAX and an independent data set (not used in developing the models) from predisturbance surveys. Although using the same data that were used to fit models can lead to positively biased estimates of correct classification rates, the independent data from pre-disturbance surveys were collected with less rigorous methods. Typically, these surveys covered larger areas and trees were often not climbed (Biswell et al. 2000, 2002, Huff et al. 2012), leading to a larger percentage of tree vole nests that are not detected than in the random plots. Thus, we expected correct classification rates from the independent data to be negatively biased. Using both data sets provides what we view as reasonable bounds on estimates of correct classification rates.
We used pre-disturbance surveys conducted in 1999-2013. Observers visually searched for signs of tree voles along transects (Biswell et al. 2000, 2002, Huff et al. 2012. For survey areas with ≥ 1 tree vole detected (hereafter, "presence polygon"), we used only the location of each tree vole nest. Survey areas without tree vole detections (hereafter, "absence polygon") were searched throughout the polygon, and thus, we assumed tree voles were absent throughout, and therefore, we used randomly selected 1-ha cells, as described below.
To compare performance among models from the pre-disturbance surveys, we sampled from the presence and absence polygons and estimated habitat suitability among selected 1-ha cells. We used ArcGIS to create a grid of 1-ha square cells across the study area. We omitted cells that (1) were modified by timber harvest or fire between the time that surveys were conducted and the date of the environmental data we used in the models (LandTrendr, Kennedy et al. 2010), (2) were previously used to develop models for DA-MAX, (3) did not have their centroid or majority area within survey polygons, and (4) where GNN data were not available. For analysis of presence, we used all cells with ≥ 1 tree vole nest detected and that met the above criteria. For analyses of absence, we used a random set of cells equal in number to the sample of presence cells. For each cell, we used the maximum estimate of the probability of occurrence from DH-GAM and our best GLM for the 30 9 30 m map pixels within the 1ha cell. DA-MAX estimates were based on 1-ha resolution, and therefore, we used the relative habitat suitability of that cell. There were a total of 1014 survey polygons with tree voles (presence polygons) and 2172 polygons where tree voles were not detected (absence polygons) that met our criteria for inclusion in the analyses. There were 3667 and 18,887 cells within the presence and absence polygons, respectively. We included all 3667 cells with detections of tree voles from the presence polygons, and a random selection of 3667 cells from absence polygons.
We also evaluated classification rates for predictions that were consistent among all three models (hereafter, ensemble model) using both random plots and pre-disturbance survey data. This model can be considered a special case of a consensus ensemble model where more than a single model is used for prediction (Ara ujo and New 2007). To compute classification rates, we used only random plots or cells from survey polygons with consistent predictions from the three models as either present or absent. For pre-disturbance surveys, we used all of the cells with consistent predictions from presence polygons, but selected a new set of random cells from absence polygons equal to the number of cells from presence polygons.

Geographic boundaries
The 364 random plots were distributed across most of the study area with the exception that only a few plots were sampled in the northern Coast Ranges (Fig. 1), where federal lands were uncommon. Tree vole nests were located in 95 (26.1%) of the plots, but their frequency of occurrence varied tremendously across the study area.
The location of a plot within the study area was one of the most predictive elements of the models we evaluated. Of those based on physiographic provinces, the Forsman et al. (2004) geographic subregions had the greatest number of different strata and the lowest AIC ( Fig. 3c; Appendix S8). The South Coast and Interior Southwest subregions had the highest (44%) and lowest (7%) occurrence rate, respectively (Fig. 3c). The North Cascades and North Coast subregions had considerably lower rates of occurrence relative to the central portions of the range (Fig. 3c). These general patterns were evident with the other delineations of physiographic provinces (Fig. 3). Survey zones, which were delineated based on broad physiographic provinces and tree vole distribution patterns (Huff et al. 2012), had the most distinct differences of occurrence rates (Fig. 3d). Tree voles were detected in <5% of plots outside of all survey zones, 37% of plots within the mesic zone, 28% within the north mesic zone, and 13% in the xeric zone (Fig. 3d).
Density contours was the best-performing covariate of those in the geographic theme, and carried through to hybrid models (Table 1). The 20% contour had the highest occurrence rate (60%), followed by the 10% contour (47%; Fig. 5). There was strong support for differences between each contour compared to the ≥95% contour, and there was weaker support for differences among the other contours ( Table 2). The predictive ability of contours was greatest when combined with other variables and was an important covariate in models in which it was included (Table 1; Appendix S11). Fig. 5. Percentage of plots with red tree vole nests detected among five contours that best described patterns of density of red tree vole occurrence based on kernel density estimators. Because of the lack of plots with tree vole nests detected in the outermost contour, we pooled data within the 95% and 100% contours in the habitat suitability models. Notes: Density contours (CONTOURS) were estimated from a fixed kernel density estimator; estimates for each probability distribution are relative to the ≥0.95 contour. MAXDBH is the maximum tree size (cm), DDI is the diameter diversity index, L_dist_dunk050 is the log of the distance (m) to suitable habitat as defined by model DH-GAM with a 0.50 threshold, AGE-CLASS is the forest age class comparing young (≤80 years) relative to old (>80 years) forests, and MAXDBH 9 AGECLASS is the interaction.
Results of other geographic covariates further described the high spatial variation of tree vole occurrence. Variation among administrative units was high, ranked second in geographic models, and had frequencies of occurrence ranging from 3% to 88% of plots (Fig. 4). Although watersheds described spatial variation in tree vole occurrence well, the large numbers of watersheds led to low precision and were poor predictors of tree vole occurrence (Appendix S8). Distinct population segment, spotted owl critical habitat, and land-use allocations were all very poor predictors (Appendix S8).

Environmental covariates from previous models
Of the DH-GAM covariates, maximum tree diameter was the single best predictor of tree vole occurrence (Appendix S1) and was on average 41 cm greater in plots with tree voles (Fig. 6), with an even greater difference in young forest ( Table 3). The ground-based measurement of maximum tree diameter preformed much better than that derived from GNN (ΔAIC = 14.7). The standard deviation of the diameter of conifer trees was on average 48% higher in plots with tree voles than those without (Fig. 6) and was the second best predictor of the DH-GAM covariates (Appendix S1). None of the other vegetation covariates in DH-GAM explained much of the variation between plots with and without tree voles (Appendix S1). We found that linear forms of vegetation covariates ranked higher than nonlinear relationships in our GLMs (Appendix S1). Slope had low predictive ability (Appendix S1) and little difference between plots with and without tree voles (Fig. 6).
In our GLMs, the environmental variables in DA-MAX were less predictive than those in DH-GAM (Table 1; Appendix S2). Diameter diversity index and the number of large trees per ha had the greatest difference in mean values between plots with and without tree voles, but there was high overlap in their distributions for these and the other DA-MAX covariates (Fig. 7). Diameter diversity index and the logarithmic relationship of large trees per ha were essentially equal as best DA-MAX covariates in our GLMs, followed by fog index, food source as a logarithmic response, and as the least predictive covariate, conifer cover (Appendix S2).
Tree vole occurrence was strongly associated with forest age. Occurrence rates were over three times higher in plots in old forest (30.9%) than in young forest (8.9%). Although there was a positive relationship with the probability of occurrence of tree voles and forest age for old-forest plots (b = 0.0045 AE 0.0020), age class was a better predictor of tree vole occurrence than forest age as a continuous variable (Appendix S9).

Distance-based covariates
Of all of the distance factors we modeled, only distance to suitable habitat performed reasonably well. The logarithmic relationship of distance to suitable habitat, derived from DH-GAM with a 0.50 threshold, had the greatest predictive ability (Table 1; Appendix S6). There were relatively large differences in the distribution of distances from suitable habitat between plots with and without tree voles (Fig. 8). However, the large differences existed in part because plots that were within suitable habitat (i.e., distance = 0) were more likely occupied (Fig. 8). Furthermore, areas considered suitable habitat were sufficiently common that a majority (55.2%) of plots were within 100 m and almost all (92.3%) were within 1 km of suitable habitat. Distance from areas defined as suitable habitat based on DH-GAM with a threshold value of 0.25 had reasonable predictive ability. The other models, which included those derived from DA-MAX and those that included distance to larger patches of suitable habitat, ranked much lower (Appendix S6).
Of covariates related to plot location relative to the study area boundaries, the logarithmic relationship of distance from the center of the study area performed best (Appendix S7). Plots with tree voles tended to be closer to the center of the study area (96.6 AE 6.2 km) than those without tree voles detected (140.7 AE 3.7 km). Distance to the perimeter of the study area was least predictive of all the distance-based covariates (Appendix S7). However, none of these models were competitive with the other best-of-theme models (Table 1).

Patch size and neighborhood effects
The best patch size model was the logarithmic relationship of patch size as estimated with DH-GAM, with a threshold ≥0.25 (Table 1; Appendix S5). The effect of patch size was largely due to the estimate of suitability at the scale of the 1-ha plot. Patch size was similar for plots with (1985 AE 633 ha) and without (1574 AE 722 ha) tree voles when only plots defined as suitable habitat were compared. Thus, the effect of patch size was largely a function of the existence of suitable habitat and not the extent of such habitat.
However, we detected an association of the suitability of areas surrounding plots with tree vole occurrence. Mean suitability of habitat within 4-and 40-ha areas surrounding plots differed in plots with and without tree voles, whereas there was little difference at the scale of the fifth-field watershed and sub-basin. The logarithmic relationship of suitable habitat within 4 ha estimated with DH-GAM had the lowest AIC (Table 1; Appendix S3). The linear form of this model was competitive, whereas models estimated with DA-MAX had much lower predictive ability for both spatial scales (Appendix S3). Plots with tree voles

Hybrid models
Our best GLMs resulted from combining topranking models from various themes. The bestperforming model was BEST9, which included covariates from DH-GAM (maximum tree diameter) and DA-MAX (diameter diversity index), in addition to density contours, distance from suitable habitat, forest age class, and the interaction between maximum tree diameter and forest age class (Table 1, Fig. 9). The second best GLM (BEST10) included the same covariates as the best-performing model, including the interaction term, but included all of the non-geographic covariates from DH-GAM and DA-MAX. The third ranking GLM had the same covariates as BEST9, but included administrative units rather than contours (Appendix S11). These three GLMs accounted for 90% of the relative model weights. All of these models included a spatial covariate (contours or administrative units), tree diameter, and the interaction between maximum tree diameter and forest age class (Appendix S11). Identical models without spatial covariates were not nearly Fig. 9. Effects of covariates from the logistic regression model BEST9. For each covariate examined, we held the other continuous covariates at their mean values and explored changes in the probability of red tree vole occurrence as the covariate increased in value and as different levels of the categorical variables were examined. The histograms represent the number of plots for the given continuous covariate for plots where tree vole nests were detected (top) or not detected (bottom histogram). Fig. 8. Comparison of distance (m, log-transformed) from suitable habitat between plots with and without detections of red tree vole nests for all plots and for only those plots that did not contain suitable habitat (i.e., >0 distance). Suitable habitat for this covariate was estimated by model DH-GAM with a probability threshold of 0.50 that delineates suitable from unsuitable habitat. Values (ln distance) in the box represent the middle 50% of the measurements, the uppermost and lowermost boundaries of the box represent the upper and lower quartiles of the measurements, and the tails represent the 90th and 10th percentiles. Filled circles outside the tails represent extreme values. The solid line through the box represents the median, and the dotted line the mean.
as competitive, and those without the interaction between age and diameter had much lower AIC weights, but were among the top six models (Appendix S11). Mean maximum tree diameter was similar between young-and old-forest plots with tree voles, but was almost 50% lower in young-forest plots without tree voles (Table 3). This strong association explains the inclusion of the interaction term in the top three models.

Performance of models
Although the distribution of predicted probabilities of tree vole occurrence from the top GLMs, DA-MAX, and DH-GAM overlapped between plots with and without tree voles, mean probabilities were greater in plots with tree voles than in those without voles. DA-MAX had the least difference in mean predicted habitat suitability values between plots with (0.43) and without (0.26) tree voles, whereas the top three GLMs and DH-GAM had an approximately threefold difference in mean predicted probabilities. DA-MAX resulted in predicted probabilities that were most dissimilar to the results from the other models we evaluated (r ≤ 0.56), whereas estimates for DH-GAM were similar to BEST9 (r = 0.86) and were similar to all of the other models that included maximum tree diameter either as a GNN variable (r = 0.75) or as a groundbased measurement (r ≥ 0.85). Overall correct classification rates (~80%) and Cohen's kappa (0.52) were highest for DH-GAM and BEST9 (Table 4).
When models were applied to pre-disturbance survey data, correct classification rates and Cohen's kappa were generally high for cells where tree voles were detected (~75%) for DH-GAM and BEST9, but low for DA-MAX (Table 4). Correct classification rates from absence polygons were low and almost identical for the three models (Table 4). The frequency of tree vole occurrence was positively related to the predicted probability (Fig. 10). Although this relationship was similar among the three models, DA-MAX had much greater precision than DH-GAM and BEST9 (Fig. 10). At the threshold value where predictions were classified as suitable habitat, 46% of the 1-ha cells were predicted to have tree voles, and this was consistent among the three models (Fig. 10).
The ensemble model improved classification rates considerably (Table 4). Of the 364 random plots, 228 (62.6%) of the plots had consistent outcomes (predicted as tree voles present or absent) among the three models. Of these, 94.6% and 80.2% were classified correctly for plots with and without tree voles detected, respectively, resulting in a 63% improvement over chance (Table 4). Of the 7334 cells in the pre-disturbance survey data that we evaluated, 4078 (55.6%) had consistent outcomes, with 87.2% and 67.2% correctly classified for cells with and without tree voles detected, respectively, a 50% improvement over chance (Table 4). Notes: Tree vole presence was predicted where the probability of tree vole occurrence was ≥0.25 (models DH-GAM and BEST9) or estimates of relative habitat suitability were ≥0.34 (DA-MAX); absence was predicted when values were less than these. The ensemble model includes only plots and cells that had consistent predictions of presence and absence among the three models. We defined absence (specificity) as plots or cells where tree voles were neither detected nor predicted, and presence (sensitivity) as areas where they were both detected and predicted to occur. Cohen's kappa is a measure of the proportional improvement over chance, ranging from 0 to 1.

DISCUSSION
Public agencies use wildlife surveys as a primary tool to address regulatory requirements for species of conservation concern. For most species, especially those that either live in difficult terrain to survey or are hard to detect, surveys are a major expense and may fail to address management questions. Management of the red tree vole provides an excellent case study of how project-level surveys have been required as a way to inform management, while other tools such as habitat models and patterns of occurrence have rarely been implemented. Our development and evaluation of these tools provides compelling support for their use in meeting mitigation requirements and reducing the need for project-level surveys of tree voles. By using predictions from three models developed with different statistical approaches and different environmental and spatial covariates and assumptions, the ensemble model provides a tool to predict tree vole occurrence with relatively high accuracy, and may have greater application than surveys conducted at project-level scales.
Our analysis of frequency of occurrence and our top habitat suitability models point to only a few important covariates associated with predicting tree vole presence. Geographic covariates were in all of the top models reflecting the high spatial variation in the occurrence of tree voles, consistent with earlier studies. Forsman et al. (2016) concluded that tree voles were widely distributed in western Oregon, but that abundance of the species was highly variable, ranging from rare or uncommon in the northern Coast Ranges and northern Cascades to fairly common in old forests in the southern Coast Ranges and central Cascades. Our findings, using randomly collected samples, were consistent with their results, but the level of spatial variation that we estimated was much greater, occurring at finer spatial scales than the subregions used by Forsman et al. (2016). Although we were able to capture spatial variation by using geographic covariates, we were not able to elucidate the mechanisms, be they environmental, biogeographic history (Miller et al. 2006), or patterns of timber harvest and wildfire (Price et al. 2015) or other human uses. The spatial variation of tree vole occurrence throughout their range likely represents many of these factors. Fig. 10. Relationship between the frequency of occurrence of red tree vole nests and the predicted probabilities of relative habitat suitability of (A) DA-MAX, and occurrence from (B) DH-GAM and (C) BEST9 within 0.05 suitability bins. Occurrence rates are the proportion of 1-ha cells from the pre-disturbance surveys (n = 7552 cells) in which red tree vole nests were detected. Linear regression (solid line) was based on the midpoints of each probability/relative suitability bin. Dashed line represents the threshold probability we used to delineate predicted absence or presence of tree vole nests for each model.
In addition to spatial covariates, the only environmental covariates included in the top models were based on forest age, as reflected by both frequency of occurrence and that tree size and forest age class were important covariates in all of the top models. Red tree voles occur in Douglas-fir forests as young as 22 years (Swingle and Forsman 2009), but have a higher frequency of occurrence and greater abundance in old, but not necessarily old-growth forests (Corn and Bury 1986, Gillesberg and Carey 1991, Dunk and Hawley 2009, Forsman et al. 2016. Furthermore, relative abundance of tree voles was likely biased low in old stands because of the difficulty of detecting nests in large trees Forsman 2009, Marks-Fife 2016). Our finding that the occurrence of tree voles was associated with an interaction between maximum tree diameter and forest age class supports the consistent findings of an association of tree voles with the largest trees (Gillesberg and Carey 1991), particularly within young forests (Thompson andDiller 2002, Forsman et al. 2016). The arboreal nature of tree voles also adds a third dimension in understanding abundance patterns: We would expect greater numbers of voles in older forests simply because there is a larger volume of habitat to occupy than in stands with younger trees, but this is insufficient to explain the very large differences in occurrence and abundance. Tree voles may have higher densities in larger trees in part because there may be more suitable foundations upon which they can build their nests (Gillesberg and Carey 1991, Meiselman and Doyle 1996, Thompson and Diller 2002. Our finding of a negative relationship of the probability of tree vole occurrence with distance to suitable habitat is consistent with the hypothesis that for species with limited dispersal, proximity to occupied habitat affects occupancy. Although there are only two studies that reported on distance to neighboring nest trees as a factor affecting the distribution of tree vole nests (Vrieze 1980, Meiselman andDoyle 1996, but see Thompson and Diller 2002), the challenge in obtaining reliable data to test this makes such evaluations difficult. The random plot data do not allow a direct test because of the limited 1-ha extent of each plot and the large distances between nearest plots. Although we were unable to evaluate distance to habitat that was known to be occupied, our top models all included distance to areas estimated to be suitable habitat. A more direct approach would be to use distance to large trees rather than model-based estimates of suitable habitat. Light detecting and ranging (LiDAR) data may allow reasonably accurate estimates of distances to large trees and may improve the performance of models (Johnston and Moskal, in press) and provide greater biological insight by using covariates directly related to the factors believed responsible for tree vole occurrence.

Performance of models
Overall correct classification rates of our best GLMs and DH-GAM were high with the random plot data and performed reasonably well with the independent data. We predicted that accuracy of the habitat models would be lower using the independent data set, comprising pre-disturbance surveys, than data from the random plots. Our prediction was based on four points: (1) Test performance is almost always better when using data that were used to fit the model (Johnson 2001, Boone and Krohn 2002, Boyce et al. 2002, Fielding 2002; (2) pre-disturbance survey data were not collected with the same rigor as the random plot data, thus increasing the likelihood of false absences; (3) for DH-GAM and BEST9, remotely sensed GNN data were used rather than actual field measurements; and (4) pre-disturbance surveys were conducted only in areas considered potential habitat for tree voles, thus reducing the breadth of environmental conditions that facilitates distinguishing between presence and absence sites (Osborne andSu arez-Seoane 2002, Dunk et al. 2004). As expected, performance was lower with the pre-disturbance survey data than with the random plot data. However, all three models (DH-GAM, DA-MAX, and BEST9) indicated that the proportion of sites with voles detected increased as estimates of habitat suitability increased. This is a measure of model calibration and is expected with habitat suitability models that perform well (Boyce et al. 2002).
Using all three models in an ensemble model improved classification rates substantially. Predictions from our ensemble model are more likely to result in reliable predictions when the model is applied to novel areas than any of the singular models, consistent with the findings from Latif et al. (2013). The ensemble model could only be applied to the portion of the study area (~63%) that had consistent predictions of tree vole presence and absence among the three models. Predictions from the ensemble model may be sufficiently accurate to guide management, and areas where predictions from the three models were not consistent may warrant further evaluation or data collection (Latif et al. 2013). In these cases, it may be best to use one of the singular models, such as BEST9, or to rely on data on the primary drivers of all models and field studies: older forest, large trees, and geographic location. Furthermore, when feasible, using ground-based vegetation data rather than remotely sensed data should improve estimates, as our comparison of models incorporating GNN and ground-based measurements of maximum dbh demonstrated. Collecting ground-based forest structure data would be much more practical in terms of cost and time than surveying for red tree voles with sufficient effort to obtain a high detection probablity. LiDAR also offers potential for greater accuracy in assessing tree vole habitat (Johnston and Moskal, in press).
What is "good enough?" Although there are many approaches for evaluating predictive models (Fielding and Bell 1997), there has been little guidance in assessing how good is good enough for use as a management tool. The nature of natural resource models requires the manager to identify how well a model should perform for its intended use (Dunk et al. 2004). One cannot "validate" such a model because there is no gold standard upon which to make a definitive conclusion (Johnson 2001, McDonald andManly 2001). Some models have been so regularly applied that their use is almost unquestioned (e.g., timber growth and yield models). We suggest that predictive models of a species' presence could be used in the same way: to make informed decisions about potential actions, or to inform the choice among a series of alternatives. Hurley (1986) warned that models that may be valuable management tools are not used because the criteria for acceptance are too restrictive. These are the challenging questions that managers must decide. Our perspective is that the results of our evaluation, using a suite of management-oriented models and several extensive data sources, resulted in habitat models and an understanding of the patterns of occurrence that can be applied as tools for management of tree voles in a manner that surveys cannot. Predictions of habitat suitability from these models may be more useful than the simple recognition of whether a particular location might be occupied or not at the time of a single survey.

Management implications
The SM Program uses two approaches to mitigate against potential harm to tree voles. The first is to protect occupied areas by placing 4-ha nodisturbance buffers centered on nest trees (USDA andUSDI 2000b, Huff et al. 2012). The second is to conduct site assessments to determine whether an area is a "non-high-priority site" or to develop a high-priority site conservation strategy at the scale of a fifth-field watershed to avoid habitat disturbance in those areas (Huff 2016). The criteria for site assessments consist of the quantity and distribution of occupied sites or suitable habitat, and the measures that are in place to ensure persistence of the species (USDA and USDI 2012, Huff 2016). Surveys have been the primary tool to address these requirements. Our best-performing habitat models, particularly the ensemble model, should be able to help managers decide if and where surveys are most needed to address species management questions. Dunk et al. (2004) proposed using habitat suitability models to help guide the need for predisturbance surveys of terrestrial mollusks. They suggested that surveys could be avoided in areas where well-performing habitat models predicted very low or very high suitability by assuming the area was either not occupied or occupied, respectively. While well-performing habitat models can potentially be used to fulfill the objectives of field surveys, some initial survey, such as the red tree vole random plot surveys, is needed to develop and test the validity of those models. Developing broader conservation plans for a species, assessing the need for a species to remain on the Survey and Manage list, and developing site assessments all require modeling approaches that can predict species occurrence patterns accurately enough to satisfy policy requirements, guide management to provide for population persistence, and avoid litigation. As with any model, conditions beyond those initially sampled need to be considered prior to application of the model (Johnson 2001).
Habitat suitability modeling has grown into a relatively mature and robust ecological discipline (Elith andLeathwick 2009, Franklin 2009). The advent of better (e.g., high spatial accuracy) field data, remotely sensed data, and increased computing power has all added to this. However, researchers' development and use of these modeling approaches has not translated to their use by land managers. Responsibility lies on both sides. Researchers need to ensure that their applied research projects are really that, and able to be applied. Managers need to be open to the use of new approaches, and to ask critical questions about their utility. Collaboration at the beginning of projects between researchers and managers is most likely to lead to appropriate development and subsequent use of models to address management-related questions.