Insectivorous bat occupancy is mediated by drought and agricultural land use in a highly modified ecoregion

California's Central Valley, one of the most productive agricultural regions worldwide, is home to a high number of at‐risk species due to habitat conversion. Amplifying the issue, the Central Valley faces severe droughts, creating water scarcity in surrounding natural areas. At least 14 insectivorous bat species live in this region, and prior studies show mixed results regarding the impact of agriculture and drought on bats. The aim of this study was to investigate how bats use agricultural areas during drought.


| INTRODUC TI ON
Habitat loss and climate change are two of the greatest threats to global biodiversity (Fahrig, 2003;Thomas et al., 2004). Conversion of natural ecosystems, primarily for agriculture, is widespread-in North America, this conversion is estimated at 50% (Hoekstra et al., 2005;Watson et al., 2016). Climate change is further modifying natural areas as extreme events including drought, flooding and heat waves become increasingly common (Luber & McGeehin, 2008;Trenberth et al., 2013). The combination of habitat loss and climate change drives species dispersal, phenology shifts and population declines (Morrison et al., 2012;Newbold, 2018;Parmesan, 2006).
Responses of individual species to these stressors vary, but specialists are often more likely to decline (Travis, 2003;Davies et al., 2004, but see Prugh et al., 2018).
Bats are a diverse taxonomic group with high interspecific variation in behavioural and physiological adaptations (Jones et al., 2009;Russo & Ancillotto, 2015). However, relatively "slow" life histories of bats (Barclay & Harder, 2003;Promislow & Harvey, 1990), characterized by low fecundity and high survival, make bats susceptible to stressors that cause adult mortality Jones et al., 2009;Kunz & Lumsden, 2003;Voigt & Kingston, 2016). Causes of adult mortality, such as disease, wind energy development and habitat loss, have caused dramatic reductions in North American bat populations Rodhouse et al., 2019;Weller et al., 2009). To recognize future change, we need to understand bat distributions under current climate and habitat conditions. Bats are impacted by agricultural intensification at broad spatial scales, but responses vary depending on crop type, surrounding habitat and cultivation strategy (e.g. organic versus conventional and magnitude of pesticide input; Williams-Guillen et al., 2016). Annual row crops and rice lack structural variation, decreasing abundance and activity of many bat species (Gehrt & Chelsvig, 2003;Toffoli & Rughetti, 2017;Williams-Guillen et al., 2016). Conversely, bat activity in orchards is often similar to natural areas, a relationship often attributed to habitat structure (Braun de Torrez, 2014;Williams-Guillen et al., 2016).
As duration and intensity of drought increases, the association between bats and agriculture may change. Highly modified agricultural areas could become suitable, providing a consistent source of surface water and crop pests, especially during drought (Elmqvist et al., 2016). During drought, irrigation diversions may deplete natural waterways, as less water is available (Faunt et al., 2016).
Water is diverted from natural areas using levees, which lack structural complexity, invertebrate diversity and shade (Hodkinson & Jackson, 2005;O'Brien et al., 2006). Although not a direct surrogate for natural areas, bats may track manmade and irrigation water sources to find resources in otherwise depleted landscapes (Amorim et al., 2018;Razgour et al., 2011). Drought can increase crop pest abundance by accelerating pest life cycles and driving early pest migrations (Leschin-Hoar, 2015; Trumble & Butler, 2009), a resource bats may exploit. Alternatively, the combination of drought and agriculture may cause bats to track resources to less disturbed areas.
Current understanding of the relationship between agricultural lands and bat biodiversity has largely focused on fine spatial or temporal scales. Our objective was to examine how bats use agricultural lands at a broad spatial scale during and after drought in one of the most productive agricultural regions in the world-California's Central Valley. We used single-species occupancy modelling to investigate our objective and hypothesized bat occupancy would be greater in agricultural areas during drought, especially orchards.
However, we also hypothesized overall occupancy would be lower for all species during drought. The Central Valley has faced fluctuating water conditions and intense drought in recent years (Faunt et al., 2016). The drought period from 2012 to 2016 represents the four driest consecutive years on record in California, a consequence of record high temperatures and record low cumulative precipitation (DWR, 2016). This extremely dry year was followed by atmospheric river storms in winter of 2016-2017. These storms are created by concentrated streams of water vapour that can result in powerful rainfall and flooding, which led to runoff conditions greater than 150% of average in 2017 (DWR, 2017). These contrasting extremes make the Central Valley a region where the effects of drought and agricultural land use on bats can be studied in tandem.

| Study area and site selection
California's Central Valley is 60-100 km wide and 720 km long (total area 47,000 km 2 ). Elevations range from 3 m to 91 m (WRCC, 2016). Climate varies along a north-south gradient, with a hot Mediterranean climate in the Sacramento Valley to the north and a Mediterranean steppe climate/low latitude desert in the San Joaquin Valley to the south. Summers are hot and dry and winters cool and damp (WRCC, 2016), with most precipitation occurring between October and May (TWC, 2017). This region has been converted from grasslands, riparian forests and wetlands to over 59% agriculture land (Boryan et al., 2011;Faunt et al., 2016;Sleeter, 2008;Soulard & Wilson, 2015). The dominant vegetation cover types are row crops, orchards, rice, grassland, shrubland, wetlands and wooded riparian corridors (Sleeter, 2008). The human population of the Central Valley was 6.8 million in 2010 and is projected to double by 2060 (Great Valley Center, 2014).
We chose sites using a spatially balanced random sample of hexagons from the United States Department of Agriculture (USDA) Forest Inventory and Analysis programme grid (hexagon radius ~ 2.6 km). We took equal-sized random samples from 10 to 20 roughly equal area zones of hexagons in the study area. From this sample, and opportunistic outreach, we contacted over 1,200 private and public landowners. Land access constraints resulted in a less than random sample, so we applied model-based inference to investigate our research objective (Gregoire, 1998). Specifically, we relied on covariate associations from regression to draw inferences and used those associations to predict occupancy across the study area (see Section 2.5). From positive landowner responses, we selected 1-2 survey locations within each hexagon and spaced sites at least 1,000 m apart (Rich et al., 2019). We chose this distance to maximize spatial replicates across the study area and target other study taxa (birds, mammals, herpetofauna;see Rich et al., 2019). All survey locations were treated as sites in analyses. Points were stratified by vegetative life-form: crop, orchard, grassland, alfalfa, rice, wetlands and riparian.

| Acoustic recording and analysis
We surveyed 274 sites between March and July of 2016 and 2017 ( Figure 1). Detectors were deployed for 7 nights, but equipment failure and logistics resulted in a range of 4-7 successful data collection nights per site. We sampled sites from March to July deploying approximately 30 detectors at a time. At each site, we deployed an SM3BAT detector with a wind-baffled ultrasonic microphone (SMM-U1; Wildlife Acoustics, Inc.). We elevated microphones using 3 m metal conduit and directed them towards open areas or flyways. We deployed recorders for 22 min 30 s per half hour from 30 min before sunset until 0,400. We activated the acoustic channel for 7 min 30 s to record audible nocturnal activity for a concurrent study.
We identified recorded calls to species using Kaleidoscope Pro V4.3.2 (Wildlife Acoustics). Noise filtering and auto-classification was completed using Kaleidoscope Pro's auto-classifier. Output included a maximum likelihood estimator p-value (p) per night, per site (Wildlife Acoustics, 2017). We vetted species-night-site events with a match ratio of p < .05 (high probability of positive species ID).

F I G U R E 1
Locations of 274 sites across California's Central Valley, plotted over USDA CropScape reclassified land cover. Sites were selected using a spatially balanced, partially random design in 2016 (black, n = 90) or 2017 (green, n = 184) A single trained observer confirmed or rejected each filtered event after review of at least 5 recording files per event.

| Occupancy and detection covariates
We expected land cover, distance to wooded wetland, site water, year and latitude would influence bat distribution (Smith, 2019).
We extracted land cover for 2016 and 2017 from USDA CropScape Cropland Data Layer (Boryan et al., 2011) and reclassified rasters into seven categories: row crops (all ground-based crops), orchard/ vineyard (all tree crops), semi-natural open habitat (grassland, shrubland, barren, and fallow), rice, herbaceous wetland/open water, wooded wetland and developed lands using "raster" in R version 5.3.2 (Hijmans, 2018; Figure 1). We extracted land cover percentages within buffered radii of 500 m, 1,000 m and 2000 m using each acoustic detector as centroid. These covariate buffer sizes represent expected sample unit area for bat species with varying home range sizes or movement ecology. If average use area was not available, we used influential radii for similar species. Western red bat (Lasiurus blossevillii) and canyon bat (Parastrellus hesperus) were modelled at a 500 m scale (Nicholls & Racey, 2006;Walters et al., 2007). Hoary bat (Lasiurus cinereus), California myotis (Myotis californicus), little brown bat (Myotis lucifugus) and Yuma myotis (Myotis yumanensis) were modelled at a 1,000 m scale (Evelyn et al., 2004;Henry & Thomas, 2002). Lastly, western mastiff bat (Eumops perotis) and Brazilian free-tailed bat (Tadarida brasiliensis) were modelled at a 2000 m scale (Olimpi, 2017;Vaughan, 1959). Euclidian distance to wooded wetland was used to describe distance to natural riparian edge habitat (Pierson et al., 2006). In addition to remotely sensed data, we reviewed site photographs to determine whether water was present in at least one of eight site photographs or absent at deployment to account for irrigation ponds, standing water and seasonal flooding. Year was a fixed effect with 2016 representing the drought year and 2017 representing post-drought. We included latitude to account for regional variation in occupancy from the northsouth regional climate and habitat differences between the northern Sacramento Valley and the southern San Joaquin Valley (Soulard & Wilson, 2015). We included interaction terms between year and the three crop types to test whether bat use in irrigated areas varied pre-and post-drought (e.g. year*row crop).
We expected Julian day, vegetative clutter, wind speed, call quality and serial correlation to affect landscape bat detectability (Smith, 2019). We included Julian day to model seasonal variation not described by measured covariates, such as bat and insect prey phenology and differences in detectability based on temperature and season (Hayes, 1997). Julian day was correlated with daily mean temperature (R = .71, p < .01) and was used as a proxy. We extracted vegetation clutter from 30m resolution NLCD 2011 canopy cover layers and calculated total tree canopy cover percentage within a 100m radius buffer (Homer et al., 2012;Stilz & Schnitzler, 2012), the visit level sample unit, to account for cover influence on sound transmission and call shape (Parsons & Szewczak, 2009;Patriquin et al., 2003). High wind reduces activity and recording quality (Parsons & Szewczak, 2009) so we extracted daily maximum wind velocity at each recorder from gridMET surface meteorological data (Abatzoglou, 2013), which estimates wind speeds at a 4 k resolution, but does not account for terrain. Since our vetting methods were biased towards high quality calls, it is possible we did not detect bats from low call quality recordings. The number of automated IDs obtained is a function of call quality and determined by many factors (Parsons & Szewczak, 2009). The number of calls available for vetting can increase species confirmation, so we included a natural logtransformed variable of acoustic calls auto-classified to species, at the ith site, on the jth occasion. This covariate serves as a proxy for unmeasured fine-scale features affecting detection probability (e.g. site-level precipitation; Banner et al., 2018). Lastly, we investigated serial autocorrelation in the data using a join count chi-squared test (Wright et al., 2016) which suggested correlation between subsequent detector nights, and violation of the assumption of independence of each occasion (Fischer et al., 2009;Hayes, 1997). We addressed this violation of model assumptions by including a first order Markov term (i.e. lag term) in the model, where detection was dependent on the prior night's detection (Hines et al., 2010).

| Single-season occupancy models
We used single-season occupancy models, estimated in a Bayesian framework using a Markov chain Monte Carlo (MCMC) algorithm, to investigate drought and agriculture influence for eight bat species (Appendix S1; Kery & Royle, 2016;MacKenzie et al., 2002). We chose not to use multispecies occupancy modelling (e.g. Iknayan et al., 2014) so we could model occupancy at different spatial scales.
We fit models using JAGS version 4.3.0 through package "R2jags" in program R version 5.3.2 (Plummer, 2003;R Core Team, 2018;Su & Yajima, 2015). We treated each survey location as a sample unit and each night (30 min before sunset until 0,400) as a separate survey occasion (n = 4-7). Because bats are a highly mobile species, they can cover many sample units during a given survey night. As such, "occupancy" hereafter should be interpreted as the probability a species used a given sampling unit during the surveyed period (MacKenzie, 2005;MacKenzie et al., 2018). We standardized all continuous covariates the same mean (0) and standard deviation (1) and tested for collinearity using Pearson's correlation coefficient. We used vague priors (Kery & Royle, 2016;Northrup & Gerber, 2018) and treated all terms as fixed effects. For intercept values, we used a Uniform (0, 1) prior on the probability scale. For all logit-scale parameters, we used a Normal prior with a mean of 0 and variance of 2 [Normal (0, 0.5 τ)]. We evaluated sensitivity to priors by considering a Uniform (−10, 10) for all logit-scale parameters (model coefficients) but found no difference in posterior mean estimates when covariates had high influence in the model (indicator variable above 0.5).
We ran 3 independent chains of 50,000 iterations, discarded 5,000 iterations as burn-in and retained every 10th iteration to obtain 13,500 posterior samples. For all models, we assessed convergence using the Brooks-Gelman-Rubin convergence diagnostic (r-hat) and assumed convergence at < 1.1. We also visually inspected chains for convergence.
Covariate influence was determined using indicator variables (Kuo & Mallick, 1998). Indicator variable selection enabled a single model run and provided a posterior mean indicating importance of each covariate (Hooten & Hobbs, 2015;Kuo & Mallick, 1998). This latent variable (w in below model) has a specified prior of Bernoulli (0.5). We investigated covariates when the indicator variable posterior mean was included in 50% or more iterations (greater than 0.5) (Allen et al., 2018;Whitlock et al., 2020). Covariates with an indicator variable value between 0.5 and 0.9 had substantial uncertainty, with 90% credible intervals overlapping zero.
All models we specified at the ith site on the jth occasion as follows: We report odds ratios and their credible intervals obtained using the antilog of posterior probability distributions from model coefficients. We tested model goodness-of-fit using a posterior predictive check on aggregated site detection history (Kery & Royle, 2016) and calculated a chi-squared test statistic (c-hat) comparing the ratio of the test statistic for observed and expected datasets (Darryl I. MacKenzie & Bailey, 2004). Values closer to one indicate a wellfitting model. We also used semivariance plots to assess potential spatial autocorrelation in our occupancy modelling results at two spatial scales (see Appendix S2; Bailey & Gatrell, 1995).

| Range-wide occupancy estimates
We used our models to extrapolate species predicted occupancy across the entire Central Valley ecoregion. We created grids by converting the buffered circle area to a square with the same area as the respective model (e.g. 1,000 m radius buffer = 1772.45 m sided square). Land cover covariate percentages per grid cell were extracted using the "raster" package in program R version 5.3.2 (Hijmans, 2018). Euclidian distance to riparian areas was averaged per grid cell. Grid cell latitude was the centroid. Site water, a survey-specific covariate, was not included in projections, and maps represent predicted occupancy at sites without water. We created predictive maps using posterior mean estimates for each covariate.
The posterior mean beta estimates were multiplied by each respective indicator variable to adjust mean estimates (i.e. low IV, low covariate influence). Model uncertainty ranges were mapped using the same methods for the upper and lower 90% credible intervals. We mapped cumulative species occupancy by aggregating species distribution predictive maps, for a cumulative occupancy ranging between 0 and 8. Bats modelled at finer spatial scales were coarsened to a 3.5 km by 3.5 km area by taking the mean predicted occupancy within intersecting cells.
We also estimated average occupancy for species across the entire Central Valley study area (see Furnas, 2020 for a similar approach for birds). For each MCMC iteration, we drew a random sample from the prediction grid equal to the sample size of our study (n = 274) and took the predicted occupancy mean for this sample.
The posterior distribution of this derived quantity represents our expected occupancy estimate for each bat species across the Central Valley. As noted above, we took this approach because our study locations were not a random sample of the Central Valley.

| RE SULTS
We surveyed 1,888 site-nights across 90 sites in 2016 (n = 657) and 184 sites in 2017 (n = 1,231). Land ownership varied, with 99 deployments on state lands, 53 on non-governmental organization lands, 77 on private, 44 on federal and two on local land. We identified 401,501 auto-classified recording files to 17 bat species and 21,326 species-night-site events. After filtering, 7,253 species-nightsite events were vetted manually. We confirmed less than 50% of these calls, for 3,301 total species-site-night events. We confirmed 15 with manual vetting (Figure 2). We completed single-species occupancy modelling for eight bat species. Townsend's big eared bat (Corynorhinus townsendii), pallid bat (Antrozous pallidus), spotted bat (Euderma maculatum) and long-eared myotis (Myotis evotis) were excluded from occupancy modelling because they were detected on less than 2% of survey nights. Big brown bat (Eptesicus fuscus) and silver-haired bat (Lasionycteris noctivagans) were excluded from analysis because these species have similar call structures, a common source of false positives. All models converged after initial burn-in iterations, with no r-hat values exceeding 1.1. Posterior predictive checks indicated our models fit relatively well, with c-hat values ranging from 0.32 to 1.3 for all species (with a c-hat value > 1 indicating possible over-dispersion and lack of fit to model assumptions; Appendix S3).

| Single-species occupancy models
Overall, mean detection probability ranged from 0.06 to 0.98 and average area-wide occupancy ranged from 0.28 to 0.98 (Figure 3).
Mapped occupancy was highly variable, with some areas of high uncertainty due to greater uncertainty in posterior means (Figure 4).
Additionally, inclusion of interaction terms resulted in greater uncertainty in 2017. We found limited evidence of spatial autocorrelation, which varied by species and spatial scale. Our assessments were particularly limited at finer spatial scales (<10 km) due to the logit Ψ i = + 1 * Year i * w 1 + 2 * Distance to Riparian i * w 2 + 3 * Row Crops i * w 3 + 4 * Row Crops i : Year i * w 4 + 5 * Orchard i * w 5 + 6 * Orchard i : Year i * w 6 + 7 * Open i * w 7 + 8 * Site Water i * w 8 + 9 * Latitude i * w 9 + 10 * Developed i * w 10 + 11 * Rice i * w 11 + 12 * Rice i : Year i * w 12 sparse amount of data for fitting pairs of locations at small distances (Appendix S2).

| Detection covariate response
As expected, auto-classified call quality positively impacted detection for all bat species (Smith, 2019). As auto-classified files increased, detection probability also increased. This effect was weak for western mastiff bat, the lowest frequency bat of our target   Figure 4).
Cumulative species occupancy ranged from 2.13 to 7.5 (out of 8; showed an increase over most of the study area ( Figure 5). This pattern was largely driven by western red bat, hoary bat and Brazilian free-tailed bat, which all had greater range-wide occupancy in 2017.

| Land cover and occupancy
No shared habitat features influenced occupancy universally for all species; but all models provided at least marginal support for one or more habitat covariates. An increase in orchard cover by 1 SD (21%) increased odds of use for hoary bat by a factor of 3.33 (90%

| D ISCUSS I ON
Our study examined broad-scale bat distributions in the California Central Valley. In this region, both agricultural land use and year (drought) predicted bat species occupancy. We found agricultural lands may serve as important habitat refugia during drought years and provide evidence that long-distance migrants may exhibit more plasticity in habitat use when facing drought stressors.

| Drought
Not all bat species expanded their range between 2016 and 2017, which coincided with the end of an extreme drought. Prior broadscale monitoring in the Pacific Northwest found evidence of higher species turnover for tree-roosting and migrating bats and lower turnover for less migratory bats that use permanent roost structures, which our study corroborates . In 2016, during drought, migratory species occupancy-western red bat, hoary bat and Brazilian free-tailed bat-was lower, with patchy areas of high predicted occupancy. Migration is an adaptation to avoid or exploit dynamic climate conditions or resource availability (Fleming & Eby, 2003;Pettit & O'Keefe, 2017;Popa-Lisseanu & Voigt, 2009), and long-distance migrants may track resources (Constantine, 1959).
Lower observed occupancy during 2016 may represent restrictions in species distribution or landscape abundance; however, without tracking individual movements, drivers of these patterns are unknown. Migratory species appeared to use their dispersal abilities to exploit resources after drought. An alternate explanation for lower occupancy lies in irrigation management during drought conditions.
Many farms use targeted irrigation during drought, which leaves orchards with less dense foliage. Further, riparian tree survival and biomass is decreased during drought (Garssen et al., 2014). Foliage density loss in the study area may have led to the observed differences in occupancy between years for tree-roosting obligates, western red bat and hoary bat. The increase we observed between years contrasts with an overall declining trend for hoary bats in the Pacific Northwest over an interval that included our study .
Myotis species, canyon bat and western mastiff bat occupancy did not vary between years. These bats hibernate (except western mastiff bat) and prefer permanent roost structures (tree and rock crevices, buildings, bridges, caves). Bats that roost in permanent structures have greater roost fidelity than foliage roosting bats (Lewis, 1995), an advantageous trait that minimizes potential risks of dispersal. However, roost site fidelity may limit the propensity to move to track productivity. Lack of change between the two years suggests fine-scale site fidelity may slow or limit species response to broad-scale changes in climatic conditions, such as drought cessation. Alternatively, these species may all be less impacted by drought, which seems less likely given their otherwise divergent life histories, or species abundance could have been affected during drought onset, which this study could not have detected.
Arid-adapted species, western mastiff bat and canyon bat, used cultivated areas more during drought; however, there was high un- however, production remained constant in orchards, requiring an increase in groundwater pumping. Groundwater depletion and subsidence decreases surface water flow and quality in riparian areas (Scanlon et al., 2012), which may have caused arid-adapted bats to use irrigated areas (Amorim et al., 2018 (Amorim et al., 2018;Braun de Torrez, 2014;Williams-Guillen et al., 2016). This observation may also have been caused by the wide diversity of behavioural and physiological adaptations represented by the 8 bat species.

| Land cover
Notably absent from this study were negative effects of row crops on bat species occupancy (Gehrt & Chelsvig, 2003;Toffoli & Rughetti, 2017;Williams-Guillen et al., 2016). Row crops did not clearly predict occupancy for any species, but a positive effect of crops on Brazilian free-tailed bat, which forages over open habitats and exploits crop pest populations was supported by indicator variable inclusion (McCracken et al., 2012). Although rice is often managed as wetland habitat for bird species (Sterling & Buttner, 2011), this crop did not impact occupancy for any bat species except little brown bat. This parallels a European study, which found one genus of bat foraged in conventional rice lands (Toffoli & Rughetti, 2017 (Katibah, 1984), making up 0.3% of total land use.
Comparatively, over 28% of the Central Valley is used for orchards and vineyards. As native hardwood riparian forests have been converted to cropland, it seems likely that foliage-dependent bats, such as western red bat and hoary bat, would also use orchard cover as a surrogate for native hardwood vegetation. Although this effect was not apparent for western red bat, the Central Valley region has consistently been characterized as a critical area for the species (Pierson et al., 2006). The positive association suggests bats may be using orchards as alternative habitat (Pierson et al., 2006). . The exception to this pattern was Brazilian free-tailed bat, which is generally more detectable due to echolocation frequency (Parsons & Szewczak, 2009). Big brown bat and silver-haired bat, two species that produce similar calls to Brazilian free-tailed bat, were omitted from analysis to reduce possible false-positive errors. We chose to use the same model for all bat species in this analysis to compare among species; however, one consequence of this decision appears to be possible overfitting for species with few detections.

| Broad-scale acoustic monitoring
This study reinforces the value of broad-scale acoustic monitoring efforts and their ability to detect change at broad spatial scales between years (Barlow et al., 2015;Loeb et al., 2015;Neece et al., 2018;Rodhouse et al., 2019). Although our study was at a large scale relative to other agriculture studies, we were confined to a highly modified agricultural region, and the effects of agriculture on bats need to be investigated at broader spatial scales, inclusive of a broader variety of habitat types. Moreover, this study emphasizes the need to monitor effects of landscape changes made during times of scarcity, especially if these decisions displace bats from their preferred habitats. Continuous monitoring would allow wildlife agencies to measure population trends and employ a dynamic occupancy modelling framework to determine extinction and colonization probabilities (Neece et al., 2018;Rodhouse et al., 2015).
The grid we deployed was substantially smaller than recent broadscale bat monitoring efforts (Barlow et al., 2015;Loeb et al., 2015;Neece et al., 2018;Rodhouse et al., 2012); however, we found limited evidence of spatial autocorrelation. Our finer grid allowed us to integrate our data collection with other sensor-based survey methods (e.g. acoustic recorders for birds, camera traps for terrestrial mammals), maximize coverage and provide stronger inference about stressor impacts on ecoregion biodiversity (Rich et al., 2019).
We recommend that researchers weigh their specific data collection needs during the study design phase (Reichert et al., 2021), but acknowledge that broad-scale deployments (e.g. NABat grid-based master sample) may reduce spatial autocorrelation caused by bat movement and allow for a more robust interpretation of occupancy (Loeb et al., 2015).
Understanding landscape bat occurrence in areas with high potential for human-wildlife interaction is important for bat conservation. These areas have a higher risk of zoonotic disease emergence , as more people interact with bats. In addition to zoonosis threats, humans can spread Pseudogymnoascus destructans, the fungus that causes white-nose syndrome, which was detected in California for the first time in 2018 (CDFW, 2019).
Perhaps most importantly, agricultural areas are an area where bats have the largest positive impact on human economies (Boyles et al., 2011;Kunz et al., 2011). Flight allows bats to move with relative ease across a matrix of private and public lands and our study highlights the importance of outreach for representative sampling. The most effective conservation strategies for bats will need continued public land maintenance as well as supporting potential bat habitat on private land. This latter component will require working with agricultural managers to reinforce the benefits of bats and ecosystem services. To achieve this goal, education about ecosystem services provided by bats, fine-scale research on how bats use different crop types, pesticide impacts and research on agriculture intensification and drought will be required to conserve bats in the Central Valley and other agricultural regions of the world.
Our study demonstrates that agricultural lands can be an important source of habitat refugia for many bat species in the Central Valley, especially during drought, and provides evidence that some bat species may exhibit more plasticity when facing drought stressors. Although native vegetation conversion to agriculture has repeatedly had negative impacts on bat activity and diversity (Williams-Guillen et al., 2016), agricultural lands, especially orchards, may have provided a habitat surrogate for unavailable natural areas in this highly modified ecoregion.

ACK N OWLED G EM ENTS
We thank Z. September 2016).

Pe e r Rev iew
The peer review history for this article is available at https://publo ns.com/publo n/10.1111/ddi.13264.

DATA AVA I L A B I L I T Y S TAT E M E N T
Occupancy and covariate data are available on the Dryad public repository at https://doi.org/10.5061/dryad.1vhhm gqs8