Identifying climate‐resistant vernal pools: Hydrologic refugia for amphibian reproduction under droughts and climate change

Vernal pools of the northeastern United States provide important breeding habitat for amphibians but may be sensitive to droughts and climate change. These seasonal wetlands typically fill by early spring and dry by mid‐to‐late summer. Because climate change may produce earlier and stronger growing‐season evapotranspiration combined with increasing droughts and shifts in precipitation timing, management concerns include the possibility that some pools will increasingly become dry earlier in the year, potentially interfering with amphibian life‐cycle completion. In this context, a subset of pools that continues to provide wetland habitat later into the year under relatively dry conditions might function as ecohydrologic refugia, potentially supporting species persistence even as summer conditions become warmer and droughts more frequent. We used approximately 3000 field observations of inundation from 449 pools to train machine‐learning models that predict the likelihood of pool inundation based on pool size, day of the year, climate conditions, short‐term weather patterns, and soil, geologic and landcover attributes. Models were then used to generate predictions of pool wetness across five seasonal time points, three short‐term weather scenarios and four sets of downscaled climate projections. Model outputs are available through a website allowing users to choose the inundation thresholds, time points, weather scenarios and future climate projections most relevant to their management needs. Together with long‐term monitoring of individual pools at the site scale, this regional‐scale study can support amphibian conservation by helping to identify which pools may be most likely to function as ecohydrologic refugia from droughts and climate change.

be among the greatest threats to many amphibian populations over the coming decades (Miller et al., 2018;Ray et al., 2016;Scheele et al., 2012). Climate projections in many parts of the world indicate hotter and longer growing seasons coupled with intensification of both heavy precipitation events and droughts (Ahmadalipour et al., 2017;IPCC, 2013;Trenberth et al., 2013;USGCRP, 2018).
These changes have potential to degrade amphibian habitat suitability and persistence, especially for species whose recruitment depends on climate-sensitive wetlands (Walls et al., 2013).
Seasonally inundated ponds provide critically important breeding habitat to a wide variety of amphibian species (Semlitsch & Skelly, 2008), including threatened and endangered species in North America (e.g., Ambystoma cingulatum, Ambystoma bishopi, Rana sevosa, and Pseudacris streckeri illinoensis; Chandler et al., 2017;Palis et al., 2006;Richter et al., 2003;Trauth et al., 2006) and globally (e.g., Pseudophryne pengilleyi; Scheele et al., 2012). In the northeastern United States, vernal pools-sensu Zedler (2003), also known as seasonal forest ponds (Brooks, 2000) or ephemeral forest pools (Brooks, 2009)-are small, geographically isolated wetlands that fill from fall through early spring and generally become dry by mid-to-late summer (Leibowitz & Brooks, 2007). The amount of time a vernal pool is available for amphibian breeding, and larval development can range from several weeks to most months of the year; however, annual drying in most years prevents establishment of predatory fish populations (Calhoun et al., 2014;Leibowitz & Brooks, 2007).
The reproductive success of pool-breeding amphibians is closely linked to pool hydroperiod (annual duration of inundation) and to the rate and seasonal timing of pool drying (Chandler et al., 2017;Green & Bailey, 2015). Although amphibian species differ somewhat in their wetland habitat requirements and reproductive phenology, all species have certain minimum hydroperiods required for successful reproduction and larval development (Brooks, 2004;Paton & Crouch, 2002). As a result, pool-breeding amphibians are especially vulnerable to climate-driven hydrologic changes that affect seasonal availability of inundated breeding habitat (Chandler et al., 2016;Green et al., 2013;Ray et al., 2016;Walls et al., 2013).
Droughts can shorten pool hydroperiods, causing pools to fill incompletely and/or to dry more rapidly than normal (Brooks, 2004;Werner et al., 2009). These drought effects can reduce amphibian reproductive output and, in some cases, lead to recruitment failure (Chandler et al., 2017;Palis et al., 2006;Scheele et al., 2012;Seigel et al., 2006). Many pool-breeding amphibians are adapted to interannual variability in hydroperiod, and populations may persist through the years of catastrophic reproductive failure (Kinkead & Otis, 2007;Taylor et al., 2006). However, climate-change-induced increases in the number of years with reproductive failure could adversely impact populations over time (Chandler et al., 2016).
Climate change can alter the hydrology of vernal pools through processes that govern both water inputs and losses. Because pool filling and maintenance of water levels are often strongly coupled with the timing and amount of precipitation, vernal-pool availability for amphibian breeding may be sensitive to changes in precipitation magnitude and timing-including more intense storm events separated by longer dry periods-and by snow-to-rain transitions in winter precipitation (Brooks, 2004(Brooks, , 2009Walls et al., 2013). Water loss from vernal pools occurs primarily through surface evaporation, soil drainage (i.e., shallow groundwater recharge), and lateral water transfers to adjacent uplands to replace soil water lost to transpiration (Brooks, 2009;Leibowitz & Brooks, 2007). The relative importance for pool drying of these water-loss processes depends in part on the geometry and substrate characteristics of pool basins, along with local groundwater characteristics and landcover and topography of the landscapes in which pools are embedded (Brooks & Hayashi, 2002;Calhoun et al., 2014;Jones et al., 2018). Because evapotranspiration is strongly temperature dependent, overall growing-season warming and earlier seasonal temperature increases (i.e., earlier spring increases in evaporative demand, along with earlier leaf-out leading to earlier ramp-up of transpiration) are expected to increase pool-drying rates and have the potential to shift pool-drying dates earlier in the season (Brooks, 2005(Brooks, , 2009Leibowitz & Brooks, 2007;Rodenhouse et al., 2009). Combined climate impacts on both water inputs to, and losses from, vernal pools have raised concerns about hydroperiod shortening, potentially interfering with amphibian reproduction (Brooks, 2004(Brooks, , 2009Chandler et al., 2016;Palis et al., 2006;Rodenhouse et al., 2009). However, vernal pools-even those in close geographic proximitymay vary substantially in their hydrologic responses to climate drivers, owing to differences in pool characteristics (e.g., size, morphology and groundwater interactions), as well as in features of their immediate surroundings and the larger landscape (Brooks & Hayashi, 2002;Chandler et al., 2017;Jones et al., 2018;Leibowitz & Brooks, 2007).
Variability in the underlying geomorphology and landscape setting of vernal pools suggests that some pools may be less hydrologically sensitive to climate change than others, for example, less responsive in terms of drying date and hydroperiod. In this context, pools that retain water later in the growing season under adverse climatic conditions (e.g., droughts, longer, and/or hotter growing seasons) may provide hydrologic refugia for pool-breeding amphibians. Hydrologic refugia (Cartwright, Dwire, et al., 2020;McLaughlin et al., 2017) are an important subset of climate refugia more broadly, which are locations that may help buffer natural communities from the effects of climate change . Vernal pools that exhibit slower or more gradual hydrologic shifts in response to climate change thus could support species persistence in stochastically variable climate regimes and potentially buy time for more comprehensive management responses to climate impacts . While identification of hydrologic refugia is an emerging component of wetland conservation in arid and semi-arid landscapes (Cartwright, Dwire, et al., 2020;Davis et al., 2013;Russell et al., 2020), few studies have sought to identify hydrologic refugia among wetlands in humid temperate regions, such as vernal pools of the northeastern United States.
We examined the drivers of pool-inundation patterns over the seasonal window associated with pool drying (i.e., late spring through mid-summer) in the northeastern United States. Using approximately 3000 inundation observations from 449 vernal pools located on protected areas across eight northeastern states from West Virginia to Maine, we trained machine-learning models to predict the likelihood of pool inundation based on seasonal timing (Julian day), precipitation inputs and drought indicators, pool characteristics (size and morphology metrics), soil and landscape characteristics, and long-term climate variables. These predictors represented hypothesized relationships drawn from the existing literature on drivers of vernal-pool inundation (Table 1). We then used the models to generate predictions for pool-inundation probability under various scenarios defined by combinations of seasonal time points, weather conditions, and future climate projections. The resulting inundation-probability predictions can be used to identify subsets of pools that might function as hydrologic refugia under a variety of seasonal, weather, and climate contexts, using various thresholds to define vernal pools as refugia (e.g., depending on species of interest). These model predictions can support vernal-pool managers in planning for climate adaptation by suggesting which pools may be most vulnerable to hydrologic effects from climate changeperhaps warranting additional hydrologic monitoring-and, conversely, which pools may warrant further study as potential hydrologic refugia for pool-breeding amphibians. Units contained from 1 to 89 pools (mean = 28 pools per unit). These pools were previously identified using probabilistic sampling methods, described in detail by Van Meter et al. (2008) and summarized briefly here. First, a 500-m grid of points was generated for each unit, and a random sample of points was used to define 50-m by 50-m plots that were visited by field crews in March through June, 2004 and For each plot in which a pool was found, neighbouring plots were established and searched; this process continued until no additional pools were discovered. Additionally at one unit (Patuxent Research Refuge), leaf-off colour infrared digital photographs were used to identify potential pools, a subset of which was field verified in May through June 2004.
Pools in this study were represented as point locations with horizontal positional accuracy ≤10 m. Pool elevation ranged from 1.4 to 1036 m above sea level. Mean annual temperature and mean annual precipitation ranged across units from approximately 4.1 C to 13 C and from 92 to 139 cm, respectively.

| Pool inundation observations
Pool inundation observations were conducted as part of the US Geological Survey's Amphibian Research and Monitoring Initiative (ARMI) wetland-breeding amphibians surveys. Pools were generally visited between March and July, from 2004 through 2016, though the specific dates varied by latitude according to amphibian breeding phenology (i.e., seasonal timing). At each visit, the inundated width and length of each pool were recorded, from which we calculated where A is the elliptical area as a proxy for pool inundated area, and L and W are pool inundated length and width, respectively (Brooks, 2005). Inundated depth was also recorded at each visit at the deepest point in each pool. Because pools were rarely observed to be dry in March or April and we were primarily interested in dry-down patterns later in the season, we used observations from May through July in inundation modelling (n = 3004 observations from 449 pools).
Our dataset did not include geomorphic or bathymetric surveys from which to calculate pool surface area or volume. Therefore, as a rough proxy for pool surface area we calculated mean April inundated area for pools with non-zero April inundation observations. This metric correlated strongly with maximum observed pool inundated area and with mean observed March inundated area (Spearman's ρ = 0.94 and 0.92, respectively), suggesting its usefulness as an overall indicator of pool surface area. We calculated pool area-to-depth ratio by dividing average April inundated area by average April inundation depth. Larger area-to-depth ratios generally represent pools with flatter, shallower bathymetry whereas smaller area-to-depth ratios represent deeper, narrower pools.

| Climate and landscape data
To evaluate the effect of climate-change projections on the hydrology of vernal pools in this study, we compiled 1-km gridded climate variables relevant to climatic water balance for a baseline historical period   (Daly et al., 1994) 17.2 (18.5) Positive Brooks (2004Brooks ( , 2005Brooks ( , 2009); Greenberg et al. (2015); Leibowitz and Brooks (2007) 6-month standardized precipitation evapotranspiration index Hydrogeomorphic setting (e.g., pool location on slopes, flats, depressions and other landforms) influences pool hydroperiod through surface-runoff and groundwater dynamics. g Soil Group A indicates well-drained soils with lower runoff potential; Group D indicates poorly drained soils with higher runoff potential.
h Inundation in pools on more permeable geologic substrates tends to be maintained in part by groundwater rather than solely by precipitation.
i Greater forest cover surrounding pools is expected to increase water losses by transpiration, however, greater canopy shading could reduce pool water losses from surface evaporation. j The National Land Cover Database (NLCD) 'cultivated crops' and 'hay/pasture' categories were combined into a single category 'agriculture'. The 'shrub/scrub' and 'herbaceous' categories were combined into a single category 'grass/shrub'. Only those landcover categories representing at least five vernal pools were retained for modelling.
precipitation inputs (on the scale of days) in addition to longer-term weather indicators (on the scale of months or years) because variability in short-term precipitation accumulation has been shown to be an important driver of pool inundation patterns (Brooks, 2004). For short-term water inputs, cumulative antecedent precipitation was calculated as the sum of daily precipitation from the Parameter-elevation Regressions on Independent Slopes Model (PRISM; Daly et al., 1994) for a specified number of days (3, 5, 7 and 10) leading up to the day of each inundation observation. After performing a sensitivity analysis (Table S2), we chose 5-day cumulative antecedent precipitation to represent short-term water inputs. To represent medium-term waterbalance conditions, Standardized Precipitation Evapotranspiration Index (SPEI) data were obtained from the WestWide Drought Tracker (Abatzoglou et al., 2017). SPEI is conceptually similar to the Standardized Precipitation Index (SPI) but also incorporates temperaturerelated effects on potential evapotranspiration (Vicente-Serrano et al., 2010). SPEI accounts for the difference between precipitation and potential evapotranspiration over varying time frames, normalized to local conditions. For example, for a June observation of pool inundation, 3-month, 6-month, and 12-month SPEI represent precipitation and potential evapotranspiration from April-June, January-June, and May-June, respectively. We performed a sensitivity analysis using these SPEI time frames (Table S2) and selected 6-month SPEI to represent medium-term water-balance conditions. Finally, we repre-  Table S1). We anticipated that pool inundation would be positively related to 5-day cumulative antecedent precipitation and SPEI, and negatively related to indices of aridity such as annual heat-moisture index (Table 1).
We also compiled landscape attributes (e.g., landcover and soil characteristics) that are known to influence local hydrology and pool-inundation patterns (Table 1;  Inset map in lower right shows vernal-pool locations at an example refuge, Lake Umbagog NWR therein). All dataset values were extracted directly to pool point locations, except for landcover type and forest aboveground biomass, which represent the upland areas surrounding pools. To assign values to pool locations for these datasets, we performed a sensitivity analysis using three possible buffer distances around pool point locations, with radii of 20, 50, and 100 m (Table S3). All data compilation was conducted in the R statistical language (R Core Team, 2017) using the packages sp, raster, rgdal, mclust and rgeos.
Data and R processing scripts used in this study are available from Cartwright, Morelli, and Grant (2020).

| Inundation models
We represented vernal-pool inundation using four binary inundation We investigated potential drivers of pool inundation by constructing separate boosted regression tree (BRT) models for each of the four inundation metrics. BRT models are a class of machinelearning algorithm in which regression trees are fitted iteratively, with the primary trees explaining the largest deviation in the response variable (in this case, binary inundation metrics), with subsequent trees seeking to maximally explain the remaining error (i.e., to explain distribution in the residuals). Advantages of BRT models over other classification frameworks include their ability to capture complex interactions and non-linear relationships and their capacity to handle missing values and multiple predictor types, that is, continuous and categorical variables (De'ath, 2007;Elith et al., 2008;Hastie et al., 2001). All BRT models were constructed with Bernoulli error distribution, tree complexity (i.e., number of splits) = 5, bag fraction = 0.5, and 10-fold cross-validation. We used the gbm.step procedure in the R package dismo to systematically select the optimal learning rate and number of trees for each model so as to maximize the deviance explained and minimize overfitting, in each case ensuring >1000 trees for each model (Elith et al., 2008;Hijmans et al., 2016).
We evaluated BRT model performance using three statistics: (1) the percentage of deviance in the response variable (inundation metric) explained by the model, (2) the cross-validated correlation coefficient and (3) the mean area under the curve of the receiveroperator characteristic (ROC-AUC). ROC-AUC values can range from 0 to 1, representing the probability that a randomly selected positive inundation observation (pool was observed to be wet) has a higher predicted wetness probability than a randomly selected negative observation (pool was observed dry). Thus, an ROC-AUC value of 0.5 would represent model performance no better than chance, and values increasing towards 1 indicate improved model ability to discriminate wet from dry inundation observations.
We constructed preliminary BRT models for the H1 through H4 inundation metrics using all candidate predictors in Table 1 ('full models'), a total of 22 predictors describing pool attributes, landscape characteristics, and weather and climate variables relevant to pool inundation patterns (see references in Table 1). We then conducted a three-step model simplification procedure to arrive at more parsimonious final models (Tables S4 through S6). The goal of model simplification was to remove predictors that were redundant and/or that could be removed without adversely affecting model predictive performance (Elith et al., 2008). Ten predictors were removed to produce simplified final models with 12 predictors. Model simplification did not substantially affect model performance (Table S6).

| Model predictions
BRT models for inundation metrics H1 through H4 were used to generate predictions of pool wetness under a variety of conditions.
Predicted wetness probability was expressed as continuous probability values ranging from 0 to 1. To compare observed and predicted inundation for model accuracy calculations, a threshold of 0.6 was applied to these continuous probabilities to differentiate predictions of 'likely dry' from 'likely wet' pool inundation status, because this threshold maximized model prediction accuracy while approximately balancing type 1 and 2 error rates ( Figure S1).
To test model prediction accuracy, we conducted leave-one-out cross-validation (LOOCV) by systematically withholding all observations from one pool at a time (as a test dataset) and generating predictions using all observations from all other pools (Tables S7 and S8).
Thus, LOOCV accuracy statistics for each pool represent the effectiveness of the modelling approach in predicting that pool's inundation patterns using a model trained on all other pools in the dataset.
Because pools varied in their LOOCV accuracy and in the numbers of observations per pool (from one to 17), we assigned pools to five categories of confidence in model predictions based on combined criteria for LOOCV accuracy and numbers of inundation observations. These confidence categories allow users of the prediction datasets to screen pools based on relative confidence in the accuracy of model predictions.
We created plots of predicted wetness probability on a daily time the three weather scenarios assumed the corresponding SPEI value on that day and the corresponding cumulative precipitation for the five days leading up to that day (e.g., SPEI = 1 and a total of 26.0 mm of precipitation from June 1 to June 5 for a prediction on June 5 under the wet-conditions scenario). For each scenario, prediction curves for each pool were used to calculate inferred drying dates, that is, the Julian days at which the predicted wetness probability curves fell below the 0.6 threshold separating 'likely wet' from 'likely dry' status.
We then generated 300 sets of modelled wetness probability predictions for the pools in the dataset, using the four inundation definitions (H1 through H4), three weather scenarios (dry, average and wet), five climate scenarios (historical, and 2050s and 2080s under suggest potential for increased summer drying of some vernal pools.  (Table 2).

| Drivers of pool inundation
Climate and weather patterns also influenced pool inundation probability. Both 5-day cumulative antecedent precipitation and 6-month SPEI had generally high relative influence values across models (Table 2) Figure 3g) and elevation had generally low relative influence across models ( Table 2), suggesting that elevation was not a strong driver of pool inundation after accounting for weather and climate conditions.
Overall, the strongest predictor of inundation probability was average April inundated area (Table 2), with larger-area pools more likely to be inundated throughout the season (Figure 3d). Although average April area and April area-to-depth ratio were strongly positively correlated across pools (Spearman's ρ = 0.89), they showed opposite patterns in partial-dependence plots (Figure 3d,e). This finding suggests that, for pools of a given size at a given seasonal time point, pools with deeper, narrower geometry were more likely to be inundated than pools with broader, shallower geometry, possibly related to the role of surface evaporation in pool drying.
Several landscape characteristics also helped explain pool inundation patterns, although landscape characteristics were generally less influential than pool surface area, seasonal timing (i.e., Julian day), and recent precipitation (Table 2). Predicted wetness probability showed moderate declines with increasing forest aboveground biomass in a 50-m radius around pools (Figure 3h). Wetness probability was generally lowest for pools in predominantly agricultural areas and highest for pools surrounded by herbaceous wetlands (Figure 3j). location-showed such minimal influence in preliminary models that they were dropped by the model simplification procedure (Table S4).
F I G U R E 3 (a-l) Partial-dependence plots from simplified boosted regression tree models, relating likelihood of vernal-pool inundation (vertical axes) with timing, climate, pool and landscape characteristics (Table 1). Coloured lines represent loess-smoothed averages across 10 iterations per model. To reduce the influence of outliers, interpretation of relationships is focused on the region of each plot between the dashed vertical lines, representing 5th and 95th percentiles of the predictors. H1 through H4 inundation metrics are binary representations of pool inundation defined by the following thresholds: H1, inundated depth and area >0; H2, inundated depth ≥5 cm and area ≥5 m 2 ; H3, inundated depth ≥10 cm and area ≥15 m 2 ; H4, inundated depth ≥15 cm and area ≥25 m 2 . 6-month SPEI is standardized precipitation evapotranspiration index calculated using 6-month antecedent conditions. In plot j, landcover categories Evergreen, W. wetlands, developed, deciduous, and H. wetlands represent evergreen forest, woody wetlands, developed open space, deciduous forest, and herbaceous wetlands, respectively

| Predictions of pool wetness
Pool-specific model prediction accuracy based on LOOCV (i.e., the percentage of a given pool's inundation observations that were modelled correctly using a model trained on all other pools) ranged from 0% to 100%, with mean prediction accuracy ranging from approximately 80% to 87% and median prediction accuracy ranging from 90% to 100% across the H1 through H4 models (Tables S7 and   S8). Using combined criteria of LOOCV prediction accuracy and number of observations per pool, approximately 54% of pools were assigned confidence of 'moderate' or better, 42% had confidence of 'moderately high' or better, and 13% had 'very high' confidence in model predictions (Table 3).
Pools demonstrated a variety of modelled seasonal dry-down patterns in response to short-term weather and future climate scenarios.  Note: Values are presented as means (standard deviations in parentheses) across 10 iterations for each model. H1 through H4 inundation metrics are binary representations of inundation defined by the following thresholds: H1, inundated depth and area >0; H2, inundated depth ≥5 cm and area ≥5 m 2 ; H3, inundated depth ≥10 cm and area ≥15 m 2 ; H4, inundated depth ≥15 cm and area ≥25 m 2 . For each inundation metric, the top 3 predictors with greatest mean relative influence are bolded. 6-month SPEI is standardized precipitation evapotranspiration index calculated using 6-month antecedent conditions. ROC-AUC is the mean area under the curve of the receiver-operator characteristic. All models were constructed with a learning rate of 0.015. suggest that relatively few pools across the northeastern United States will be inundated (according to either the H2 or H4 definitions) through early-to-mid July in the dry weather scenario, and that this scarcity of late-season pool inundation will be exacerbated under climate change (i.e., by the 2080s under the RCP 8.5). However, some pools show relatively high late-season predicted wetness probability even under dry weather and future climate scenarios (example in Figure 4C), suggesting their potential as hydrologic refugia from droughts and climate change.
If we define vernal-pool refugia as pools with ≥60% probability of holding any amount of water (H1 definition) by July 1 under dry weather and historical climate conditions (  (Paton & Crouch, 2002).
Average April inundated area was a strong predictor of observed inundation patterns (Table 2), with larger-area pools having more Numbers of observations in the stated criteria refer to the number of inundation observations for each pool (which ranged from 1 to 17); accuracy refers to pool-specific accuracy calculated based on leave-one-out cross-validation (supporting information sect. 5), with accuracy having been averaged for each pool across the H1 through H4 inundation metrics. b Numbers of pools and percentages of pools for each confidence category represent the pools meeting the stated criteria minus the pools having a higher category of confidence. c For each confidence category, the cumulative percentage of pools represents the percentage of pools having that level of confidence or greater.
frequent inundation than smaller-area pools (Figure 3). This relationship was reflected in predicted wetness probability values, which were positively correlated with pool area across inundation metrics, prediction dates, and weather and climate scenarios (Tables S9 through   S12). Thus, larger pools tended to have greater predicted wetness probability and might be more likely to serve as hydrologic refugia than smaller pools under a range of weather and climate conditions.

| DISCUSSION
In this study, pool-inundation models successfully represented inundation patterns for several hundred vernal pools at a regional scale, using information on a variety of landscape and pool-specific attributes in conjunction with seasonal timing and climate variables. While relationships of some of these predictors to pool inundation were well understood from previous investigations (e.g., antecedent precipitation) others have been less well studied (e.g., forest biomass surrounding pools; Table 1). A machine-learning approach enabled us to evaluate >20 such potential influences on inundation, including categorical variables, while ultimately drawing conclusions from more parsimonious, simplified models (Table 1). In contrast to previous studies that examined hydrologic dynamics in only one or a handful of pools at the site scale (e.g., Brooks, 2000Brooks, , 2004Dodd, 1994;Greenberg et al., 2015)  In vernal pools and related seasonal wetland ecosystems, hydrologic regimes are of primary importance in determining habitat quality for amphibians. Hydroperiod influences wetland geochemistry, water quality, litter decomposition, carbon cycling, invertebrate community composition and food webs (Babbitt et al., 2003;Boven et al., 2008;Brooks, 2000;Cartwright & Wolfe, 2016;Chandler et al., 2016;Golladay et al., 1997). Potential exists for earlier and/or more rapid pool drying-for example, from drought intensification and increased evapotranspiration under climate change-to interfere with reproductive success and long-term population viability for a number of amphibian species of conservation concern (Brooks, 2009;Chandler et al., 2016;Daszak et al., 2005;Greenberg et al., 2015;Ray et al., 2016;Walls et al., 2013). These species generally function as habitat specialists, and many have larval growth and development periods that are adapted to historical ranges of timing for pool drydown (Babbitt et al., 2003;Baldwin et al., 2006;Brooks, 2004;Seigel et al., 2006). While timing of metamorphosis is a plastic trait for many temporary-wetland-associated species, minimum hydroperiods are necessary for successful recruitment and faster metamorphosis may induce a lifetime fitness cost (Amburgey et al., 2016;Cabrera-Guzmán et al., 2013;Semlitsch et al., 1988). As a result, shifts in seasonal timing of pool filling and drying of even a few weeks can affect populations of temporary wetland-breeding amphibian species (Baldwin et al., 2006;Calhoun et al., 2014;Seigel et al., 2006).
While there are a number of putative causes of amphibian population declines, droughts and associated early pool-drying events have been linked to catastrophic reproductive failure, population crashes and local extirpations of pool-breeding amphibians (Palis et al., 2006;Richter et al., 2003;Scheele et al., 2012;Seigel et al., 2006;Semlitsch, 1987;Taylor et al., 2006). Many poolbreeding amphibian species exhibit high levels of interannual variability in reproductive success and have evolved to cope with F I G U R E 5 Distributions of predicted wetness probability of vernal pools in the modelling dataset for (a and b) the H2 inundation metric and (c and d) the H4 inundation metric, using model outputs with (a and c) historical climate variables and (b and d) climate projections for the 2080s under RCP 8.5. Dry, average and wet weather conditions are represented in orange, grey and blue, respectively. Because H4 is a more stringent metric than H2 for classifying pools as inundated, distributions of predicted wetness probability for H4 (c and d) are generally lower than for H2 (a and b). Distributions represent only the subset of pools with 'moderately high' or 'very high' prediction confidence stochastic years of reproductive failure (Brooks, 2004;Semlitsch et al., 1996;Taylor et al., 2006). Although periodic reproductive failure may not result in metapopulation-level declines, increases in the frequency of years with reproductive failure could harm longterm population viability (Brooks, 2004;Daszak et al., 2005;Taylor et al., 2006).
In this context, effective management of vernal pools and the amphibian populations they support requires the ability to anticipate potential hydrologic changes, which in turn requires understanding of the interacting factors that shape vernal-pool hydrologic regimes.
Findings from this study confirm that vernal-pool hydrology is influenced by a complex array of drivers, ranging in scale from regional Predictions are displayed only for pools categorized as having 'moderately high' or 'very high' prediction confidence based on cross-validation accuracy and numbers of observations; all other pools are displayed as hollow triangles. In (d), blue triangles (i.e., pools with ≥60% predicted wetness probability) could indicate hydrologic refugia to site-specific (Brooks, 2004(Brooks, , 2005Brooks & Hayashi, 2002;Calhoun et al., 2014;Grant, 2005;Leibowitz & Brooks, 2007). Poolinundation models in this study successfully represented known seasonal dry-down patterns (Brooks, 2004(Brooks, , 2005, with inundation likelihood decreasing from early May to late July (Figures 3 and 4).
Seasonality was represented by Julian day in models, which serves as a proxy for accumulated water losses from evapotranspiration. Julian day was not strongly correlated with other model predictors (Table S5) including recent precipitation, in agreement with the assessment by Brooks (2004Brooks ( , 2005) that vernal-pool inundation seasonality is largely driven by seasonality of evapotranspiration rather than precipitation. As expected, dry-down curves varied substantially among pools (e.g., Figure 4). Relative-influence values (Table 2) suggested strong roles of precipitation and evapotranspiration over relatively short (daily-to-monthly) time scales in shaping pool hydrology, in agreement with previous studies (Brooks, 2004;Davis et al., 2019;Greenberg et al., 2015). This importance of weather conditions was notable in the modelled dry-down curves for many pools under dry conditions relative to wet or average conditions (examples in Figure 4) and in distributions of inundation probability under various weather scenarios ( Figure 5). Overall, this study lends credence to concerns that vernal-pool hydrology may be sensitive to droughts and, consequently, that habitat for pool-breeding amphibians could be vulnerable to drought intensification and stronger evaporative demand under climate change (Brooks, 2004(Brooks, , 2009Greenberg et al., 2015;Leibowitz & Brooks, 2007;Ray et al., 2016;Rodenhouse et al., 2009;Scheele et al., 2012).
Differences among vernal pools in their hydrologic regimes and their responses to weather and climate drivers can be linked to a variety of pool-specific attributes. We found that larger and deeper pools may be more likely to remain inundated later into the summer, consistent with our hypothesis (Table 1) and with previous studies and expectations based on the physical dynamics of surface evaporation (Brooks, 2005;Brooks & Hayashi, 2002;Leibowitz & Brooks, 2007;Chandler et al., 2017; but see Snodgrass et al., 2000). Future modelling efforts could benefit from more detailed data on pool-basin bathymetry, such as from Lidar or field surveys (Brooks & Hayashi, 2002;Faccio et al., 2016), to enable more precise evaluation of relationships between pool hydrologic dynamics and pool-basin area, depth, and morphology. Site-scale data on soil and stratigraphy T A B L E 4 Numbers and percentages of pools across units (protected areas) meeting one possible definition of refugia (≥60% predicted probability of holding any amount of water on July 1 under dry weather conditions) Note: (1.) In the 'number of pools column', the first number denotes all pools; the number in parentheses is the number of pools having prediction confidence of 'moderate' or better. (2.) All results in the 'refugia' columns were calculated according to one possible definition of refugia, that is, ≥60% predicted probability of holding any amount of water (H1) on July 1 under dry weather conditions. Numbers of pools meeting the definition are listed, followed by percentages of pools (in parentheses). These results show only one possible approach to defining refugia, whereas other definitions using other inundation metrics, seasonal time points, and weather conditions could also be proposed. (3.) Results in the 'refugia (confidence filter)' columns represent only those pools for each unit having prediction confidence of 'moderate' or better; confidence categories are defined in Table 3 would also be very useful, including measurements of physical characteristics of pool-basin substrates (e.g., particle size, organic matter and drainage characteristics). Lacking such site-scale data at the regional scale of this study, we relied on soil-survey information (Soil Survey Staff, 2016). Our finding of greater predicted wetness probability for pools on poorly-drained than well-drained soil types reflects the general understanding of relationships between pool inundation dynamics and soil permeability and hydraulic conductivity (Grant, 2005;Leibowitz & Brooks, 2007;Rheinhardt & Hollands, 2008). For regional-scale studies in which substrate sampling is impractical across large numbers of vernal pools, soil-survey information may allow for valuable inferences. However, we acknowledge that fine-scale heterogeneity in pool substrate characteristics may be hydrologically important yet poorly represented by soil-survey information (Calhoun et al., 2014;Rheinhardt & Hollands, 2008), which may help explain the generally low relative-influence values for the soil-survey attributes representing water-table depth and soil hydrologic group (Table 2). Similarly, further investigation of potential groundwater interactions in vernal pools would benefit from site-level hydrogeologic data collection (e.g., installation of monitoring wells and piezometers; Barton et al., 2008;Chandler et al., 2017;O'Driscoll & Parizek, 2003). In some cases, groundwater influences can be inferred from water chemistry (i.e., specific conductance) and pool substrate characteristics (Brooks, 2005;Calhoun et al., 2014).
Pool substrate characteristics and degree and type of groundwater interactions are commonly a product of a pool's topographic position and geomorphic setting in the larger landscape (Rheinhardt & Hollands, 2008). In preliminary models, we represented landscape position using 30-m landform categories (i.e., pool was located on a ridge, upper slope, lower slope, or valley), but this attribute was removed from final models through the model simplification procedure (Table S4). Instead, landscape setting may have been partially captured by the surrounding landcover attribute, for example, greater inundation likelihood in pools surrounded by herbaceous wetlands relative to more upland landcover types ( Figure 3). A primary purpose of our model of vernal-pool inundation patterns was to identify the characteristics of pools that might be most likely to serve as hydrologic refugia under droughts and changing climatic conditions. This process identified the most likely refugia from the set of monitored sites (Table 4). However, the defining features of a refugium may vary depending on location and management focus (e.g., species of interest), so we provide a website for exploring multiple possible definitions of refugia depending on geography, seasonal timing, and management context and objectives (Cartwright, 2020).
For example, the greater prevalence of refugia among more northerly pools using predicted wetness probability on July 1 (Table 4) likely reflects climatic influence on timing of snowmelt and pool drying more than any special characteristics of these pools per se. Using the website (Cartwright, 2020), managers may wish to adjust dates of interest depending on pool geography to reflect such climatic gradients in the timing of pool drying.
Ultimately, we view these efforts as only a first step in the larger process of identifying and conserving pools that could function as hydrologic refugia for amphibians of conservation concern. Although many questions remain concerning the potential function of hydrologic refugia in supporting amphibian species persistence under climate change, existing knowledge can help identify hypotheses and research needs. Even without identifying specific refugia from a known set of surveyed sites (as we do here), a landscape with a diversity of pool sizes and hydroperiods is expected to maximize the probability that some pools within a protected area could serve as refugia during droughts or other adverse climatic conditions (Dodd, 1994;Shoo et al., 2011), lending support to calls to conserve diverse types and sizes of wetlands (Babbitt et al., 2003;Calhoun et al., 2017;Snodgrass et al., 2000). Our findings support the notion that larger pools and/or those with generally longer hydroperiods may be more likely to function as hydrologic refugia under droughts or climate change (Baldwin et al., 2006;Brooks, 2009). However, we stress that smaller 'dry-end' wetlands also serve important ecological functions and should not be considered expendable (Calhoun et al., 2003;Snodgrass et al., 2000). Beyond hydrologic suitability, a number of other pool characteristics may be important in determining wetland refugial capacity, including water quality, presence of disease agents, surrounding land use, and legal protections or lack thereof (Babbitt et al., 2009;Calhoun et al., 2003;Daszak et al., 2005;Evans et al., 2017;McLaughlin & Cohen, 2013 Our use of three relatively simplistic weather scenarios and climate variables from ensemble averages across 15 climate models implies that pool-wetness predictions do not capture the full range of possible future hydrologic changes. For example, if future climatic drying exceeds that represented by ensemble averages for the region, then our projected changes in pool inundation could be conservative.
Furthermore, if both droughts and storms were to intensify, pool hydrology-and thus availability of amphibian breeding habitat-could respond in complex ways not captured by either the 'dry' or 'wet' weather scenarios (Walls et al., 2013).
A further set of caveats stems from the fact that models of pool hydrology presented here are statistical, ultimately representing correlations in the dataset rather than proof of causation. While other statistical approaches have also produced successful predictive models of depressional wetland hydrology (e.g., Barton et al., 2008;Chandler et al., 2016;Greenberg et al., 2015), process-based hydrologic models offer an alternative approach using a water-budget framework and hydrologic accounting on a fixed (e.g., daily) time step (e.g., Garmendia & Pedrola-Monfort, 2010;Qi et al., 2019;Sun et al., 2006). The mechanistic nature of process-based models enables temporal examination of states and fluxes of water throughout the atmosphere-upland-wetlandgroundwater system, which could be useful to test hypotheses and scrutinize findings from this and other studies concerning the drivers of vernal pool-hydrology. However, the process of structuring, parameterizing and calibrating process-based hydrologic models for large numbers of vernal pools at a regional scale could be complicated by issues of process complexity and spatial resolution and complexity (Clark et al., 2017), for example, given the small size of vernal pools and the substantial heterogeneity of their hydrologic responses to climate patterns. Perhaps for these reasons, process-based models have typically been applied to only one or a handful of wetlands at a time (Qi et al., 2019;Sun et al., 2006;Wolfe et al., 2004). Further examination of results from this study, then, could be performed using processbased models at a subsample of vernal pools, ideally spanning the gradient from short-to-long hydroperiods.
The utility of refugia to amphibian metapopulations depends strongly on population responses to shifting mosaics of temporally dynamic breeding habitat. Increased pool drying reduces the connectivity of the network of remaining breeding habitats (Fortuna et al., 2006). The extent to which this affects metapopulation persistence depends on species' capacities for local dispersal and migration among increasingly disconnected 'islands' of breeding habitat Gamble et al., 2007;Werner et al., 2009). For example, Green and Bailey (2015) estimate that at least 50 breeding sites are necessary for metapopulation persistence, but increased isolation will likely increase this minimum number.
As amphibian species respond to climate change through shifts in phenology and geographic ranges (Gibbs & Breisch, 2001;Lawler et al., 2010), research is needed on the relationships between these changes and shifting availability (in both time and space) of breeding habitat. Because amphibians face a variety of threats related to habitat destruction, pathogens and climate change (Rodenhouse et al., 2009;Ryan et al., 2014;Semlitsch & Skelly, 2008;Wake & Vredenburg, 2008), the identification and protection of hydrologic refugia is only one part of a much larger conservation picture.

| CONCLUSIONS AND MANAGEMENT IMPLICATIONS
We modelled vernal-pool hydrology across the northeastern United States from a large sample of habitats within Federally protected areas, using thousands of inundation observations collected over more than a decade. Variability in hydroperiod across pools implies potential for some to function as hydrologic refugia from droughts and climate change. Regional-scale models are helpful for identifying pools that could provide refugia but could be strengthened by more detailed data collection at the site scale for hydrologically relevant pool characteristics such as pool substrate and stratigraphy, canopy cover and groundwater fluxes. Better understanding is also needed on how traits of different amphibian species may affect the value of refugia in practice, including dispersal and migration, breeding-site fidelity and phenological shifts.
This and other studies have suggested that earlier pool drying has the potential to affect the availability of breeding habitats across large landscapes, which has important management implications.
Conservation strategies that focus on maximizing current-day connectivity among pools may not adequately account for connectivity changes in drying landscapes. In addition to wetland restoration efforts, artificial wetland creation and hydroperiod manipulation may be useful strategies to help support persistence of breeding habitats (Calhoun et al., 2014;Green et al., 2013;Seigel et al., 2006;Shoo et al., 2011). Identifying refugia can be one part of a comprehensive management programme to increase the probability of amphibian persistence.