Linking surface and subterranean climate: implications for the study of hibernating bats and other cave dwellers

Cavesandothersubterraneanfeaturesprovideuniqueenvironmentsformanyspecies.Theimportanceof cavemicroclimate isparticularly relevantat temperate latitudeswherebatsmakeseasonaluseof caves forhibernation.White-nose syndrome (WNS), a fungaldisease that hasdevastatedpopulationsof hibernating bats across easternandcentralNorthAmerica, hasbrought renewed interest inbathibernationandhibernaculumconditions.A recent reviewsynthesized currentunderstandingof caveclimatology, exploring thequalitative relationship between cave and surface climatewith implications for hibernaculum suitability.However, a morequantitativeunderstandingoftheconditions inwhichbatshibernateandhowtheymaypromoteormediateWNSimpacts isrequired.Wecompiledsubterraneantemperaturesfromcavesandminesacrossthewestern UnitedStatesandCanadato (1)quantify thehypothesizedrelationshipbetweenmeanannual surface temperature (MAST) and subterranean temperature and how it is influenced bymeasurable site attributes, and (2) use readilyavailablegriddeddata topredictandcontinuouslymaptherangeof temperatures thatmaybeavailable incavesandmines.OuranalysissupportsqualitativepredictionsthatsubterraneanwintertemperaturesarecorrelatedwithMAST, that temperaturesarewarmerandlessvariable farther fromthesurface, andthatevendeep within-cave temperatures tend tobe lower thanMAST.Effects of other site attributes (e.g., topography,vegetation, and precipitation) on subterranean temperatures were not detected.We then assessed the plausibility of model-predicted temperatures using knowledge ofwinter bat distributions andpreferred hibernaculum temperatures.Ourmodelunavoidablysimplifiescomplexsubterraneanenvironmentsandisnotintendedtoexplain all variability in subterranean temperatures. Rather, our results offer researchers and managers improved broad-scaleestimatesofthegeographicdistributionofpotentialhibernaculumconditionscomparedtoreliance onMASTalone.Weexpect this information tobetter support range-scale estimationofwinter batdistributions and projection of likelyWNS impacts across the west.We suggest that our model predictions should serve as hypothesestobefurthertestedandrefinedasadditionaldatabecomeavailable.


INTRODUCTION
Subterranean environments provide unique and important habitat for many organisms (Thurgate et al. 2001, Clements et al. 2006, Silva and Ferreira 2016, yet are poorly studied, in part because accessing and working in these environments can be difficult. Cave-dwelling fauna often displays unique physiological adaptations to both take advantage of the benefits of underground environments and cope with their challenges, particularly when it comes to maintaining a balanced energy budget. Cave environments generally have limited influx of energy and nutrients compared to surface environments due to the absence of sunlight and photosynthetic organisms (Engel 2007). Many cave obligates are able to survive in these low-energy environments by maintaining low metabolic rates (Hȕppop 1985(Hȕppop , 2005. This strategy may be aided by the highly stable climate conditions offered by many deep cave environments, which are likely to help overcome some of the energetic challenges associated with fluctuations in environmental conditions (Tattersall et al. 2012, Stawski et al. 2014. Furthermore, this stability protects freeze-intolerant species against seasonal freezing temperatures experienced at the surface (Novak et al. 2014) and supports specialists adapted to narrow ranges of climate conditions (Peck 1976, Howarth 1980. Some species, however, bridge the subterranean and aboveground environments. In temperate regions, the cool and stable abovefreezing temperatures of caves and mines provide seasonal refuge for hibernating bats. Like many permanent residents of cave systems, hibernating bats are adapted to a relatively narrow range of microclimate conditions necessary for successful hibernation, when they may go months without food. As heterothermic endotherms, hibernating bats incur high thermoregulatory costs while defending normal body temperature at colder temperatures, but depend on the availability of cool, stable microclimates to promote energy-saving torpor (reduction in body temperature and metabolic rate; Speakman and Thomas 2003). Caves, mines, and even deep crevices and talus slopes offer safety (e.g., from predators) as well as stable, abovefreezing temperatures that enable hibernating bats to survive cold winters (Speakman and Thomas 2003). However, these environments can also provide ideal growing conditions for the cold-loving fungus Pseudogymnoascus destructans (P. destructans) that causes white-nose syndrome (WNS; Blehert et al. 2009, Verant et al. 2012), a disease that has decimated eastern and central North American bat populations since its introduction in 2006 (Frick et al. 2015, Leopardi et al. 2015. Many studies have demonstrated that hibernaculum temperature and humidity are critical determinants of bats' ability to survive winter while also maintaining above-freezing body temperature (Thomas et al. 1990). More recently, others have suggested that hibernaculum microclimate is also a key determinant of mortality due to WNS (Langwig et al. 2012, Johnson et al. 2014. As P. destructans spreads across North America (USFWS 2019), it continues to encounter and subsequently infect new populations and species. This expanding distribution of P. destructans is particularly concerning as bat biodiversity is considerably higher in western North America than in the east (Harvey et al. 2013; though not all species hibernate). Furthermore, western bats are understudied due to the rugged terrain and remoteness characterizing much of their habitat in the west. Climate change is likely to further complicate the outlook for bats exposed to P. destructans as it is expected to influence multiple aspects of the pathogen-host relationship in ways that may compound or counter one another (Frick et al. 2010, Algeo et al. 2014. Warming temperatures may increase or decrease bats' energy use during torpor as well as the growth rate of P. destructans, and these effects are likely to vary spatially and among species. Furthermore, changes in the duration of winter will alter the cumulative impacts of these shifts on winter survival. These complex interactions make WNS survivorship in a changing climate difficult to predict. Perry (2012) synthesized the current state of knowledge of determinants of subterranean climate conditions as it relates to bat hibernation, with the goal of supporting researchers and managers responding to the emergent threat of WNS. This synthesis explored the relationship between subterranean microclimate and surface climate, and how aspects of cave setting and configuration v www.esajournals.org may affect this relationship. In summary, subterranean temperatures are thought to stabilize at approximately the mean annual surface temperature (MAST) of a given location, once a sufficient distance from surface entrances is reached such that air exchange with the surface is limited. In this dark zone, where light does not penetrate (Field 1999), temperatures are often stable, whereas temperatures close to site entrances fluctuate with outside temperatures. However, temperatures in a given location within a cave or mine may vary in response to a variety of attributes of the subterranean space itself (e.g., cave size and configuration, number of entrances, and depth below the ground) as well as attributes of the landscape setting (e.g., tree or snow cover at the site entrance, and topographic aspect of the site entrance; Table 1). Many of these attributes are expected to promote air and water flow that drives subterranean temperatures below MAST. This expectation is supported by the fact that MAST at many known hibernacula in the eastern United States is well above the hibernation temperatures typically (but not always) selected by temperate vespertilionid bats in North America (2-10°C; reviewed in Perry 2012). On the one hand, MAST is often used as a best estimate of temperature conditions experienced by hibernating bats, but on the other hand, it clearly does not capture the colder conditions vespertilionid bats prefer. This lack of congruency indicates that a more nuanced estimate of the subterranean conditions available to bats over broad spatial scales is needed. To our knowledge, none of the hypothesized drivers of subterranean temperature reviewed by Perry (2012) have been quantitatively tested across diverse subterranean environments.
We assembled temperature data from a large number of caves and mines across the western United States and western Canada to quantitatively explore the hypotheses synthesized by Perry (2012). Our core objectives were to (1) model the relationship between MAST and subterranean temperature, and (2) quantify how subterranean temperature is affected by distance from surface entrances and other site attributes. Using the resulting model, we then (3) continuously mapped the range of subterranean temperatures that would be predicted to be available across western North America if a cave or mine of sufficient depth occurs there, and thus where temperatures suitable for bat hibernation may be most likely to occur. Subterranean conditions inevitably vary due to site configuration attributes that cannot be fully accounted for, and this study is not intended to explain all inter-and intra-site variability in subterranean temperatures. Rather, we suggest that our model offers a significant improvement compared to simply assuming that subterranean climates equilibrate to MAST. This has important implications for regional-level monitoring and managing subterranean fauna attuned to these climates, particularly bats, as they are exposed to pathogens such as P. destructans, climate change, and other stressors.

Subterranean microclimate data
We compiled temperature data collected from caves and mines (139 unique sites, and one published source (Oster et al. 2012). We only included datasets for which site location was provided with positional error <10 km. Data were collected using temperature loggers (iButton, Onset HOBO loggers), which recorded temperature at subdaily intervals that varied by site. Many sites contained multiple loggers, placed at a range of distances from the site entrance (479 total loggers, mean of 3.4 per site; Appendix S1: Fig. S1). Some datasets were accompanied by metadata describing key attributes that may influence temperature at the logger placement site (e.g., distance from site entrance, description of site size or shape, and site geology) or other logger placement details (e.g., wall, crevice, or ceiling). Although many loggers recorded humidity as well as temperature, logger saturation, a phenomenon in which the logger fails to accurately record humidity after reaching readings of 100% humidity, precluded the use of data from the majority of loggers. Therefore, we focused only on temperature here. We first screened all logger temperature data for deployment errors, particularly evident recording of surface temperature prior to and following placement of the logger at the site. We restricted our analysis to data recorded during the core winter months of December-February to restrict analysis to a consistent time frame that can be considered winter at all sites across the broad latitudinal range of our study. Some loggers were deployed over multiple winters. Logger-winters during which fewer than 14 d of readings were recorded were excluded from analysis (resulting in 131 sites, 427 loggers, and 699 logger-winters). We further restricted our dataset to loggers for which distance from the site entrance was provided or for which cave maps were provided that could be used to generate a logger placement distance (54 sites, 156 loggers, and 202 logger-winters). Recorded distances ranged from 0 to 324 m from the site entrance (median = 61 m).
Because the compiled and filtered raw logger dataset contained >1.2 million observations, which were highly variable in their distribution across loggers (Appendix S1: Fig. S1) as well as their frequency and regularity of recording, we conducted all analyses using a stratified random sample from the raw logger temperature dataset. We selected 250 records from each unique logger so as to retain all sites and loggers; for loggers with fewer than 250 records (19 loggers; 5%), all records were selected. Model results were insensitive to this choice of sample size. All logger data processing and analysis were conducted in R (version 3.4.1, R Foundation for Statistical Computing, Vienna, Austria).

Surface climate and landscape data
We derived MAST and predictors representing key landscape attributes using Google Earth Engine, a cloud-based computing platform supporting large-scale analysis on an extensive catalog of remotely sensed, climatological, and other geospatial datasets (Gorelick et al. 2017). We estimated MAST at each cave site using the DayMet version 3 dataset (Thornton et al. 2019), which provides gridded daily surface temperature at 1km resolution . Daily mean temperature was first calculated as the mean of daily minimum and maximum temperature. We then estimated MAST as the 20-yr (1998-2018) mean of daily mean temperatures. We also derived site-level landscape predictors that are believed to impact winter cave temperatures (Perry 2012). We extracted elevation of each site from the Shuttle Radar Topography Mission (SRTM) digital elevation model at 30-m resolution (Farr et al. 2007). Based on elevation, we derived a multiscale topographic position index (TPI), in which canyon and valley bottoms have low position and peaks and ridges have high position. Using a moving-window approach, TPI was calculated as the elevation of a focal raster cell minus the mean elevation within a given neighborhood surrounding the focal cell. We calculated TPI for three neighborhood sizes (500 m, 5 km, and 25 km squares) and then averaged these to produce a multiscale index. We extracted Continuous Heat-Insolation Load Index (CHILI), a surrogate for effects of solar insolation and topographic shading (e.g., due to slope aspect) on evapotranspiration, from the Global ALOS CHILI product at 90-m resolution (Theobald et al. 2015). Percent tree cover was extracted from the Terra MODIS Vegetation Continuous Fields product (NASA v www.esajournals.org Fig. 1. Sites at which subterranean temperatures were recorded. Data were sourced from both caves (hexagons) and mines (triangles). Sites at which no logger distances from site entrances were recorded are shown in gray and were not used in the analysis. 2019), which estimates percent tree cover within 250-m resolution pixels. We estimated the mean annual number of snow days at 500-m resolution from the MODIS Snow Cover dataset (V6; Hall et al. 2016) by averaging the number of days per year with at least 10% snow cover over the most recent 5-yr period available (2014-2018), as well as mean annual snowpack (snow water equivalent [SWE] on April 1) for the same time period at 1-km resolution from DayMet (V3; Thornton et al. 2019). We sampled groundwater table depth estimated from compiled point observations gapfilled with a mechanistic groundwater flow model (Fan et al. 2013) as a proxy for the likelihood of water flow within caves or mines, assuming that flow was more likely at sites with shallower (i.e., closer to the surface) groundwater.

Analysis of surface-subterranean temperature relationships
We quantified how well MAST approximates subterranean temperatures and the extent to which subterranean temperatures also depend on distance from the site entrance and other site attributes described above. Based on Perry (2012), we hypothesized that winter subterranean temperatures would increase with MAST, but that caves would tend to be colder than mines (due to greater airflow; Table 1). We further hypothesized that the relationship between subterranean temperature and MAST would depend strongly on distance from the site entrance. Specifically, we predicted that subterranean temperatures would increase with distance from the site entrance and that we would observe a significant interaction between site entrance distance and MAST. We hypothesized that other site attributes discussed by Perry (2012) would have negligible effects on subterranean temperatures. However, hypotheses driving our inclusion of these predictors are summarized in Table 1. We did not assess the influence of within-cave or within-mine attributes (e.g., volume, length, number of entrances, and logger position; i.e., wall, ceiling, and crevice) because these attributes were not consistently described in enough source datasets to draw reliable inferences.
We screened predictors for multicollinearity based on pairwise correlations and variance inflation factors, removing predictors that caused standard thresholds of 0.7 and 4.0, respectively, to be exceeded (Belsley 1991, Booth et al. 1994). All predictor values were centered and rescaled prior to model fitting. We then fit a linear mixedeffects model with fixed effects that included all remaining predictors and an interaction term between MAST and distance from entrances. We included random intercept terms for sites and loggers nested within sites to control for the nonindependence of repeated temperature observations from the same logger and for conditions unique to each site that could not be measured. We also accounted for nonconstant variance in temperature among loggers as a function of distance from site entrance using an exponential covariance structure, hypothesizing that recorded temperatures would be more variable closer to site entrances. We used ANOVA to assess whether inclusion of random effects and a variance component was warranted.
Using a model averaging approach and an information-theoretic framework (Burnham and Anderson 2002), we then fit all additive subsets of this global model to draw multi-model inference regarding predictive relationships and relative importance of predictors using the MuMIn package (Bart on 2018) for R. Multi-model inference produces parameter and error estimates that are not conditional on any one model, which is particularly advantageous when multiple models have similar weights of evidence or probability of best explaining the data (Burnham and Anderson 2002, 2004, Symonds and Moussalli 2011. All-subsets model averaging is recommended over selection of candidate model sets, provided there is strong logic for potential inclusion of each of the predictors considered (Symonds and Moussalli 2011), and has been repeatedly demonstrated to support more robust inference compared to selection of a single best model (Wasserman 2000, Burnham andAnderson 2004). The subset of models carrying 95% of the Akaike information criterion (AIC) weight was averaged based on their relative AIC weights, using a shrinkage estimation approach to produce unconditional parameter estimates that are not biased away from zero (Lukacs et al. 2010). Our final inferential model was derived from this 95% confidence set. We computed model-averaged regression coefficients, unconditional standard errors, 95% confidence intervals, v www.esajournals.org and cumulative AIC weights of evidence as a measure of variable importance Anderson 2002, 2004).

Prediction of available temperature ranges
We used the model-averaged result above to predict mean winter subterranean temperature conditions potentially available across the western United States and Canada, given the presence of a cave or mine of sufficient depth to create a stable dark zone. First, as an illustrative example, we mapped the range of mean winter subterranean temperatures predicted to be available in caves and in mines continuously across western North America, with ranges defined by distance from the site entrance of 10 and 100 m. These distances were chosen to represent an arbitrary minimum distance from a site entrance at which bats tolerant of temperature fluctuations might be expected to hibernate and a distance well into the stable dark zone, respectively, but the model can be applied to any distance or range of distances defined by a particular research question or management need.
To assess potential utility of the model relative to use of MAST alone as a proxy for subterranean temperature, we used a general estimate of the optimal temperature range for vespertilionid bat hibernation (2-10°C; reviewed in Perry 2012) to restrict our mapped results to areas where caves or mines, if present and of sufficient depth, would be most likely to offer suitable hibernation conditions, and compared the resulting window to that produced by MAST alone. We also summarized model-predicted temperatures and MAST at known hibernacula to the far north and south of our region of interest (Alberta, n = 34; and Texas, n = 21) to assess the model's capacity to estimate more plausible ambient hibernaculum temperatures (i.e., more frequently within the above optimal temperature range) than MAST alone. Finally, we summarized model-predicted temperatures and MAST at winter point locations of Myotis lucifugus (n = 90) collected from a variety of sources (Appendix S1: Table S1) and assessed their alignment with the distribution of ambient temperatures observed in hibernacula of M. lucifugus (n = 32) reported in the literature (Appendix S1: Table S2).

Surface-subterranean temperature relationships
The final inferential model is summarized in Table 2. Elevation and snowpack were dropped from the model due to high collinearity with MAST (Appendix S1: Table S3). MAST and distance from site entrance were the most important predictors of subterranean temperatures, followed by whether the site was a cave or a mine. As predicted, we observed a positive relationship between MAST and subterranean temperature, such that an increase of 1 standard deviation (2.82°C) in MAST corresponded to a 1.67°C increase in subterranean temperature (95% CI: 0.785-2.55). Similarly, subterranean temperatures v www.esajournals.org increased with distance from site entrances, such that an increase of 1 standard deviation (53.4 m) in distance from the site entrance corresponded to a 1.62°C increase in temperature (95% CI: 1.16-2.08). However, these positive effects of MAST and distance from site entrance on subterranean temperature were not strictly additive, as suggested by a negative interaction term (À0.32°C; 95% CI: À0.970 to 0.327), though the confidence interval on this estimate did include zero. After controlling for these influences, mines tended to be warmer than caves (by 3.08°C; 95% CI: 1.22-4.94). Other landscape predictors, such as elevation, topographic position, and solar insolation, had less discernible effects on subterranean cave temperatures as evidenced by lower variable importance, smaller standardized coefficients, and confidence intervals that spanned zero. ANOVA model comparisons indicated that inclusion of random intercept terms for loggers nested within sites, as well as distance from entrance as a covariance parameter, was warranted (Table 3). This parameter's coefficient, estimated as À0.614 in the top-ranked model, indicated that as predicted, variance in observed subterranean temperature decreased with increasing distance from site entrances. The topranked model was 56.8 AIC units better (i.e., lower) than a null model containing random intercept terms for loggers nested within sites and exponential covariance structure, but was 86,022.6 AIC units better than a null model with no random effects or covariance structure. These results suggest that the fixed effects of MAST, distance from site entrance, and other site characteristics approximated the data well; however, explicitly accounting for residual spatial structure in the data dramatically improved model performance.

Mapping predicted available temperature ranges
By applying the inferential model to predictor values at each 1-km raster cell across western North America, we mapped mean winter temperatures predicted to be available in caves and mines if such features exist at any given site, at distances of 10 and 100 m from site entrances as bracketing examples (Fig. 2); species tolerant of temperature fluctuations can be found as close as 10 m to site entrances (e.g., Corynorhinus townsendii; C. Lausen, personal communication), whereas 100 m is representative of the stable dark zone selected by many species for hibernation. We also mapped MAST across the same extent for comparison (Fig. 2e). These maps reiterate that predicted mine temperatures are warmer than predicted cave temperatures at a given distance from the site entrance and that for a given site type (cave or mine), temperatures were warmer with increasing distance from site entrance. We also show that modeled temperatures were more moderate than MAST; they were warmer than MAST in cold regions and colder than MAST in warm regions.
To relate these maps back to vespertilionid bat hibernation, we masked areas where mean temperatures of 2-10°C are not predicted to be found anywhere between 50 and 100 m from the site entrance (where bats other than C. townsendii typically hibernate; C. Lausen, personal communication.), again comparing these hibernation-suitable windows to those predicted based solely on MAST (Fig. 3). We show considerably broader windows of potentially suitable hibernation  v www.esajournals.org 9 October 2020 v Volume 11(10) v Article e03274 conditions using our model of subterranean temperature compared to that approximated by MAST. Notably, our model more frequently predicts suitable conditions at known hibernacula in Alberta and Texas that MAST alone would predict to be unsuitable (Fig. 4). Similarly, our model consistently predicts conditions within and centered on the range of hibernaculum temperatures tolerated by M. lucifugus as reported in the literature at sites where M. lucifugus has been observed in winter ( Fig. 5; Appendix S1: Tables S1, S2).

DISCUSSION
The central aim of this study was to quantitatively test the qualitative relationships between subterranean temperature, MAST, and site attributes presented by Perry (2012), using what is to our knowledge the largest set of available cave microclimate data compiled to date. We confirmed that although MAST is a strong determinant of subterranean temperature, MAST alone is a fairly poor direct proxy for ambient temperature of caves and mines. Caves, in particular, tended to be colder in winter than MAST would suggest, and temperatures deviated more from MAST as one approaches cave and mine entrances, becoming both colder and more variable. This has important implications when estimating the conditions likely to be experienced by hibernating bats in a given location. MAST is often used as a proxy for subterranean temperature in the absence of ground-based data (Grieneisen 2011, Hayman et al. 2016). Our results demonstrate that better estimates of ambient subterranean temperatures can be made by accounting for site type (cave vs. mine), distance from site entrance, and other landscape attributes that tend to drive subterranean temperatures below MAST. This holds despite inescapable variability around predicted means due to aspects of site configuration, such as number of entrances that would modulate airflow, that could not be accounted for here.
Many site attributes discussed by Perry (2012) had little or no detectable effect on cave and mine temperatures in our compiled dataset. Our failure to detect these relationships may have been due to any effects of surface-level landscape attributes (i.e., topography, vegetation, and snow cover) dampening with greater distances from site entrances. Our ability to detect these relationships may also have been obscured by noise in the data resulting from the coarse resolution of remotely sensed metrics (i.e., 30 m to 1 km) and/ or spatial error in reported site locations, although any such noise is likely minimal relative to real variation in subterranean temperatures within and among sites. This natural variability is likely driven by variability in site configuration attributes (e.g., cave size, number and orientation of entrances, and placement of loggers on walls vs. in crevices) that were generally not recorded or reported in the datasets available to us.
Our analysis focused on estimation of subterranean temperature, but temperature is only one determinant of suitable subterranean conditions. It has been widely hypothesized that hibernaculum humidity impacts the hibernation physiology of bats (Willis et al. 2011, McGuire et al. 2017, as well as habitat suitability for other subterranean species. Thus, a location that exhibits a suitable temperature regime may not be usable for some bat species that require high humidity to hibernate normally, which includes numerous species under threat from WNS (Langwig et al. 2012). However, it is likely that hibernating bats adaptively use high-humidity micro-refugia inside hibernacula (e.g., crevices and blind adits) to avoid otherwise excessively dry conditions, thereby increasing the suitability of a site through fine-scale selection of structure and substrate. Hibernacula that contain water, whether standing, flowing, or seeping, are also more likely to offer sufficiently humid conditions. Humidity is also a critical determinant of suitable conditions for growth of P. destructans (Marroquin et al. 2017) and the severity of its physiological effects on hibernating bats (Cryan et al. 2010). Consequently, it is important for future monitoring efforts to reliably record humidity and for future analyses to estimate how subterranean humidity may vary with surface conditions and site attributes as the present study has done for temperature, although this may prove an unavoidably difficult task as current loggers saturate and fail in the high-humidity conditions of most hibernacula.
Our model enables prediction of mean ambient winter temperature expected to be available at any given location and distance into a cave or mine, should such a feature exist. Together with knowledge about the locations and conditions within hibernacula that bats tend to select, we can approximate the mean conditions experienced during hibernation and where suitable conditions are likely to occur. However, we do not recommend application of the model at the individual site level. Given the limited nature of the available data, the model is intended for use in understanding general patterns and predicting over broad geographic extents. Although the model outperforms MAST alone for predicting ambient hibernaculum temperatures, it remains a coarse estimate that should be treated as a hypothesis to be further tested and refined as additional data become available. Cave microclimates are complex over space and time, as is selection of microclimate conditions by bats during hibernation. Cave and mine attributes, such as the number and orientation of entrances, amount of airflow, whether the site includes refugia from microclimatic variation (e.g., rock crevices or blind adits), as well as how hibernating animals dynamically use these features (i.e., by moving among different locations within a hibernaculum throughout the winter), will influence the suitability of a site on a finer scale than our model can resolve. Some bat species (e.g., M. lucifugus) require hibernacula with highly stable microclimates (Thomas and Cloutier 1992), whereas other species (e.g., Eptesicus fuscus) tolerate a wider range of or more variable microclimates (Kl€ ug-Baerwald and Brigham 2017). Intraspecific temperature preference and tolerance may vary geographically (Kl€ ug-Baerwald and Brigham 2017) and may also depend on whether hibernating bats are solitary or clustered (Langwig et al. 2012, Boraty nski et al. 2015. These factors may contribute to selection of hibernacula depending more on the specific morphology of the site that governs temperature variability than on average ambient temperature. Despite these considerations, we suggest that our approach is the best available given the data, offering a clear improvement over MAST as a simple proxy for understanding the relationship between surface and subterranean conditions at this scale. Furthermore, our model's ability to predict a range of mean available temperatures over a range of distances from site entrances may more reasonably capture the range of temperatures resulting from the complexity described above in habitat selection and hibernaculum use than MAST alone. MAST has been used widely as a proxy for subterranean temperature, for example, to estimate susceptibility of M. lucifugus to WNS (Grieneisen 2011) and to predict duration of winter survivorship given species-specific physiology and exposure to WNS (Hayman et al. 2016). Our model offers a refinement to this generalization that we suggest can improve inference and prediction in similar future studies. We also suggest that the more feasible hibernaculum temperatures predicted by our model can aid prediction of species occurrence in a given region should a sufficiently deep cave or mine exist. Our model may also be applied to smaller, shallower features used by some species as hibernacula (e.g., crevices and depressions) by assuming a minimal distance from site entrance, though results of this application should be interpreted with caution because the model was not informed by data from these types of features. Although predicted hibernation-suitable windows (Fig. 3) are likely to be geographically broad for many species, these could be further refined by, for example, karst, known mine locations, or other landscape attributes of known importance to a focal species to help target survey and monitoring efforts for potential hibernacula.
The model presented here is based on historic climate normals (20-yr MAST), but could be applied to projected future temperatures as a first step to assess the potential impacts of future climate change on subterranean environments and their continued suitability for hibernating bats and other subterranean species. For example, monitoring of bat populations could be targeted in areas where climate projections suggest that suitable hibernation conditions are likely to be lost (vulnerable hibernacula) or, conversely, where they are likely to be retained (stable hibernacula that may serve as refugia for bats). Similarly, the model could help managers prepare for possible range expansions into their jurisdictions under future climate scenarios. Changing climate conditions will also affect the spread and establishment of P. destructans in subterranean environments, as well as its impacts on infected bats (Campbell et al. 2020). We expect that efforts to predict and plan for these impacts can benefit from improved estimation of subterranean microclimate conditions over broad geographic extents (e.g., western North America).
We strongly recommend that future monitoring efforts and field studies of subterranean microclimate record not only temperature and humidity at multiple sites and at multiple positions within a given cave or mine, but also record site attributes (e.g., distance from entrance, number of entrances, presence of water, and airflow) in a standardized manner. This would enable future analyses to improve understanding of surface-subterranean relationships and to more fully account for the influence of site characteristics when estimating microclimate conditions. Spatial predictions of suitable hibernaculum conditions would also be improved by integrating both temperature and humidity and their relationships with surface conditions and site attributes. In turn, the combination of estimated temperature and humidity ranges for a given potential hibernaculum would better support prediction of P. destructans impacts on the bat species that may occur there. And although the present study has focused on estimating hibernaculum conditions for bats, we expect this information to have important applications to other cave fauna that are poorly understood and far understudied relative to bats but that may be highly dependent on particular subterranean conditions to persist (Mammola et al. 2019).

CONCLUSIONS
Our study offers an improved understanding of the determinants of subterranean climate and an improved ability to predict subterranean temperatures. We suggest that further refinement can be achieved by more standardized collection of cave and mine attribute data that may influence subterranean climate, and that similar assessments of patterns in humidity are needed to more fully understand and predict subterranean conditions for cave-dwelling organisms. This information has the potential to improve understanding of the distributions and unique habitat needs of not only hibernating bats, but also a variety of poorly studied cave-obligate species, many of which are rare and endemic. Furthermore, we suggest that a quantitative understanding of the link between surface and subterranean climate can help researchers and managers to better anticipate and plan for the impacts of future climate change on subterranean environments and the species they support, as well as how changes in climate may impact disease spread. Although we have focused on the importance of subterranean climate for bat hibernaculum suitability in western North v www.esajournals.org America, with implications for the severity of WNS impacts as the disease spreads through this region, our approach and findings are applicable to subterranean environments in other temperate regions and to the wide range of organisms that depend on them.