Introduction

Bats, a diverse and widespread order of mammals that provide important ecosystem services, have been experiencing significant regional declines due to the introduction and spread of disease, wind energy development, and habitat loss [1–6]. Understanding how these threats drive trends in bat populations at a landscape scale allows the development of effective management strategies [7–9]. In particular, human land use change (e.g., urbanization, loss of forest cover, fragmentation, increase in edge density) may influence the presence of bat species due to the importance of roosting and foraging sites for bat reproduction and survival [10–13]. For instance, the presence and activity of many bat species are positively associated with forest cover [14–16], and the loss and alteration of forest cover influences occupancy rates of different habitats [17,18]. Forest fragmentation can also negatively impact bat abundance, distribution, and use of space [19–23] although in some cases forest fragmentation can have a positive effect of bat activity [24]. Other landscape factors associated with urbanization and development may also influence habitat use by bats. Roads can serve as either conduits or barriers to movement, depending on the type of road and the surrounding habitat [25–31]. For example, rural unpaved roads surrounded by forests in Indiana are selected over paved roads surrounded by development and open areas [32]. Roads also act as barriers to movements as well as sources of mortality [33]. Other factors associated with urbanization and loss of forest habitat such as poor water quality and artificial lights may also affect bat habitat use and activity [34–36]. Additionally, climate change is expected to cause widespread changes in land cover and habitat for bats [2].

Shrinking habitat for terrestrial vertebrates has been predicted for the southeastern United States, especially in scenarios with expansion of urban and agricultural areas [37]. For instance, the percentage of developed land in the contiguous United States is predicted to double from 1997 to 2025, with the greatest increase in southern regions [27], leading to significant reductions in wildlife habitat and increased habitat fragmentation. Based on historical land use change, Wear et al. [38] predict urban land coverage to at least double in the southern United States from 1997 to 2060, and predict a 7–13% decline in forest coverage across the south, with the Piedmont region losing the greatest percentage of forest (up to 21%). Loss of forest cover and increasing urbanization are likely to greatly impact the abundance and distribution of bats across the southeast but few data are available to predict how various bat species will respond.

Effective conservation strategies that address various threats to bat populations in the southeastern United States require landscape scale monitoring and analyses [2,39]. Thus, our objective was to conduct an assessment of habitat factors influencing bat species habitat use at the landscape scale, which may provide valuable insight into how predicted land use and land cover change could affect bats. We used temporally dynamic (i.e., multi-season) analyses to estimate effects of habitat metrics on species occupancy rates and calculate annual colonization and extinction rates to account for year-to-year variation for a suite of bat species encountered across South Carolina, USA. We also generated species distribution maps based on our model predictions, which may be used in the future to reveal changes in species distributions in response to changing habitat [9]. Our species-specific findings can also be used to inform landscape management decisions that may affect bat populations, guide subsequent, finer-scale investigations in species-specific patterns of habitat use, and serve as a baseline for future comparative studies examining changes in bat habitat use and species distributions over time.

Methods

Sampling design

We used the North American Bat Monitoring Program (NABat) framework to acoustically survey bat species across South Carolina. The sampling frame for NABat consists of a continuous grid of 10 x 10 km cells across North America [8]. We identified priority survey cells within South Carolina based on the NABat master sample, which uses the generalized random tessellation stratified algorithm to assign priority numbers to cells to maintain a spatially balanced and randomly distributed sample (Fig 1 and S1 Table). We conducted stationary point surveys for four consecutive nights and mobile transect surveys on two of the four nights from mid-May through mid-July 2015 and 2016. In 2015, we surveyed 35 cells: 15 with mobile transects only, six with stationary point surveys only, and 14 with both survey methods. In 2016, we surveyed the 35 cells from 2015 and three additional cells (one with mobile transects only, and two with stationary points only) for a total of 38 cells surveyed: 13 with mobile transects only, eight with stationary point surveys only, and 17 with both survey methods. We surveyed three cells with stationary point and mobile transect surveys in 2016 which were surveyed with mobile transects only in 2015 (Fig 1). Each stationary point survey began 30 minutes prior to sunset and ended 30 minutes after sunrise, while each mobile transect survey began 45 minutes after sunset and was driven at 32 kph, with duration dependent on the length of the transect (25–48 km). During each year, we surveyed the same stationary point locations and mobile transect routes within each cell where possible.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 1. Distribution of NABat priority cells surveyed across South Carolina using mobile transects only, stationary point surveys only, or both survey methods, May-July 2015 and 2016.

Physiographic regions of South Carolina are also displayed.

https://doi.org/10.1371/journal.pone.0206857.g001

We used Anabat SD2 bat detectors with directional, stainless steel microphones (Titley Scientific, Columbia, MO, USA) and 2.5 m microphone cables for both survey methods. For stationary point surveys we mounted the microphone inside a water resistant PVC housing and attached it to the top of a 1.8 m high tripod. For mobile transect surveys we placed the microphone at the center of a vehicle’s roof with no waterproof housing.

Data analysis

We hypothesized that probability of occupancy would vary by ecoregion for species with limited ranges within our study area, but that it would not vary by ecoregion for species with statewide ranges (Table 1) [43]. Thus, we included a categorical covariate (Region) based on the primary U.S. Level III Ecoregion [40] within each cell.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Table 1. Predicted and observed effects of environmental variables on the probability of occupancy for each species surveyed in South Carolina, May-July 2015 and 2016.

https://doi.org/10.1371/journal.pone.0206857.t001

We hypothesized that the effects of land cover types and forest fragmentation on probability of occupancy varied among species, based on species summer roosting and foraging site preferences [11,16,26,45,46] (Table 1). We calculated percent land coverage within each cell from the National Land Cover Database (NLCD 2011) [47] and aggregated “Pasture/Hay” and “Cultivated Crops” as Ag, all classes of development as Dev, upland forest types as Forest, and “Woody Wetlands” for bottomland forest (F.Wet). We used our reclassified NLCD 2011 data as input in Fragstats version 4.2.1 [48] and calculated Contagion (a landscape measure which increases as land cover type interspersion decreases and dispersion increases) and measures of forest and forested wetland edge density (F.ED and F.Wet.ED) within each cell.

We hypothesized that stream length would have a positive effect on occupancy of many bat species because they prefer to forage near streams and riparian areas [16,34,35,49]. Additionally, streams can act as forest edges and may be important sources of water, so we also hypothesized positive effects of increasing stream length on probability of occupancy for species that commonly forage along edges, but may not be explicitly associated with riparian areas (Table 1). We calculated total stream length within each cell using ‘NHDFlowline’ data from the National Hydrology Dataset (NHD) [50].

Due to effects of roads on species presence [25,26,28–31,51], and because different road classes are often associated with different landscapes [52], we hypothesized effects of roads on bat species occupancy would vary based on road type (Table 1). We used National Transport Dataset (NTD) RoadSegment data [53] and U.S. Forest Service Roads [54], and classified roads into four categories, primarily based on Master Address File/Topologically Integrated Geographic Encoding and Referencing Feature Class Code Definitions (https://www.census.gov/geo/reference/mtfcc.html): primary roads included divided highways with access ramps (Pri), secondary roads included highways with intersections (Sec), tertiary roads included single lane rural and city roads (Ter), and quaternary roads included forest access roads (Qua). We then calculated the length of all four road classes within each cell.

Annual turnover rates (i.e., colonization and extinction) at the landscape scale could indicate changes in species range due to land cover changes or regional threats to bat populations. Given that we only monitored for two years, we predicted that it would be unlikely for land cover change to influence turnover rates. However, regionally, the arrival of white-nose syndrome (WNS) is ongoing in the southeastern US. Therefore, we predicted higher turnover rates for MYELUSE and PESU, because they are affected by WNS in the northwestern part of our study area [55]. In addition, the act of migration could influence turnover rates, and we accordingly predicted high turnover rates for LACI because they exhibit migratory behavior [42].

Due to limited access for stationary point surveys and constraints in establishing mobile transects [56], some stationary point locations were near cell edges and short segments of some mobile transects were up to 2.2 km outside cells. Additionally, we may have detected bats within cells that entered from areas outside cells. Thus, similar to other studies [11,17,18,30,35], we buffered cell boundaries by 2.2 km before we measured covariates to include relevant landscape effects.

We used a multi-season Bayesian occupancy modeling approach to evaluate the influence of hypothesized environmental factors on the probability of occupancy for each species and to calculate turnover rates. At the sample unit-level, this approach models the probability of occupancy within each sample unit during each sampling period (in our case, each year) as a temporally autoregressive function of intercept and covariate effects with parameterization for colonization and extinction. At the survey level, this approach models the probability of detection, which is dependent on the presence of the species, as a function of intercept and covariate effects [9]. We treated cells as our sample unit, considered each night at each point or transect as a separate survey occasion, and created presence/non-detection tables for each species on each survey occasion. We never detected CORA and we detected LABOLASE in every cell each year; therefore, we were unable to model occupancy for these species. Because we surveyed each cell within one week each year, we treated populations as closed within years, and open between years and calculated turnover rates for each species between the two years as the probability that an unoccupied cell became occupied (i.e., colonization), and an occupied cell became unoccupied (i.e., extinction). We used non-informative priors, treated all terms as fixed effects, and used JAGS version 4.1.0 (http://mcmc-jags.sourceforge.net/) through package ‘rjags’ [57] in program R version 3.3.3 (https://www.r-project.org/) to fit models. We ran three independent chains of 25,000 iterations, discarded an initial 5,000 iterations as burn-in, and retained every fourth iteration for a total of 18,750 iterations per model. We assumed model convergence when the potential scale reduction factor (i.e., the Brooks-Gelman-Rubin diagnostic; R-hat) of each parameter was < 1.1. We ranked models using the Widely Applicable Information Criterion (WAIC), which we calculated for each model using the package ‘loo’ version 1.1.0 [58], and considered models closely competing when they were within 2.0 WAIC from the top-ranked model [59]. We calculated ΔWAIC from the top-ranked model, then calculated each model’s relative likelihood, and finally calculated model weights to evaluate relative support of the models. We considered a covariate effect significant if the estimated 95% credible intervals did not include zero.

To account for imperfect detection, we modeled detection probabilities for each species independently and modeled occupancy as a function of intercept only. We used nine survey variables (Survey Type, Survey Duration, Vegetation Clutter [no, medium, or high amounts of vegetation within 5 m of the detector], Survey Date, Survey Issue [detector malfunction or incomplete mobile survey due to weather], Temperature, Relative Humidity, Wind Speed, and Rain) to test hypothesized effects of single-term models and additive effects models [56]. We included the covariates from the top-ranked detection model for each species in the occupancy modeling process (S2 Table).

Based on a priori hypotheses, we used single-term models and models with additive effects of some covariates to model occupancy by each species (Table 2). We also tested null and global models. Prior to model fitting, we standardized all continuous covariates to have a mean of 0 and standard deviation of 1. We used Pearson’s correlation to test for correlations among covariates and considered those with a Pearson’s |r| > 0.7 as correlated and did not include them in the same model. Pri, Sec, and Ter road classes were correlated with Dev, Sec was correlated with Ter, and F.Wet was correlated with F.Wet.ED (S3 Table). Therefore, Ter was not included in any models. Since DAIN and MYAU are associated with forested wetlands [16,60], we substituted Forest and F.ED with F.Wet and F.Wet.ED, respectively, and omitted F.Wet.ED from the global model for these species because it was significantly correlated with F.Wet.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Table 2. Reasoning for 11 a priori occupancy models that we tested for each bat species in South Carolina, USA, May-July 2015 and 2016.

https://doi.org/10.1371/journal.pone.0206857.t002

We used two approaches to evaluate model performance. First, we performed k-fold cross-validation to calculate area under the receiver-operating curve (AUC) for the top detection model, and each model with WAIC ≤ the null occupancy model. AUC values range from 0 to 1, where 0.5 indicates no predictive power and 1.0 indicates perfect predictive performance [61]. We randomly partitioned the data five times, with 66% as training datasets and the remainder of the data as testing datasets. We reviewed partitions to be sure each training dataset included at least one cell from each of the five ecoregions, and used the same partitions to evaluate models for each species. We used the package ‘ROCR’ version 1.0.7 [62] to calculate AUC for each model. We also evaluated model performance by comparing known species distributions, based on Menzel et al. [43], to predicted species distributions. To generate predicted range maps, we measured environmental covariates in all 893 10 x 10 km NABat cells (each buffered by 2.2 km) throughout South Carolina and calculated estimated occupancy rates for each cell based on the top-ranked occupancy model for each species. We then compared these predicted occupancy maps to historical range maps for each species based on Menzel et al. [43].

Results

We recorded 61,397 call files in 2015, 21,972 of which passed our noise filter, and 65,727 call files in 2016, 42,960 of which passed our noise filter. We classified 15,292 and 27,380 call files to species in 2015 and 2016, respectively.

Six of the eight species or species groups had occupancy models that ranked higher than the null model. The null model ranked highest for EPFULANO and MYAU, NYHU had two top-ranked models, and DAIN, LACI, MYLELUSE, and NYHU had at least three models that closely competed with the top-ranked model (i.e., within 2.0 WAIC; Table 3). Ecological models performed better than null models for four of the species, but the predictive performance of the top-ranked models varied, from 0.46 for LACI to 0.76 for NYHU, and was low for many of the models (Table 3). Additionally, DAIN, LACI, MYLELUSE, NYHU, and TABR had top-ranked models with equivalent or lesser predictive performance than competing models (Table 3). All covariates in all occupancy models reached convergence, except the Piedmont region in two models for EPFULANO.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Table 3. Occupancy probability model results for each species of bat surveyed using acoustic detectors across South Carolina, May-July 2015 and 2016.

https://doi.org/10.1371/journal.pone.0206857.t003

We observed support for the Region model for DAIN and MYLELUSE (Table 3). Mean estimated probability of occupancy for DAIN was significantly higher in the Southern Coastal Plain than all other regions except the Blue Ridge (Fig 2). The mean probability of occupancy for MYLELUSE was highest in the Blue Ridge region, but it only significantly differed from the estimate for the Southeastern Plains (Fig 2).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 2. Mean estimated probability of occupancy of northern yellow bats (DAIN) and small-footed bats, little brown bats, and northern long-eared bats combined (MYLELUSE) within each ecoregion of South Carolina, May-July 2015–2016.

Blue bars indicate 95% credible intervals. Within species, regions which share a letter above their intervals are not significantly different from one another.

https://doi.org/10.1371/journal.pone.0206857.g002

The Stream + F.ED model was the top-ranked model for LACI and NYHU (Table 3). Forest edge density significantly affected the probability of occupancy for LACI, but did not significantly affect the probability of occupancy for NYHU (Tables 1 and S3). From the lowest (0.50 m/ha) to the highest (101.43 m/ha) forest edge density, LACI probability of occupancy decreased from 98% to 5%, with a steep negative slope beginning at 50 m/ha (Fig 3). Stream length did not significantly affect the probability of occupancy for either species (Tables 1 and S4). In addition to the Stream + F.ED model, we also found equivalent support for the single-term F.ED model for NYHU (Table 3), where the effect of forest edge density was negative but not statistically significant (Tables 1 and S4).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 3. Estimated effect of forest edge density on hoary bat (Lasiurus cinereus) probability of occupancy across South Carolina, May-July 2015 and 2016.

Probability of occupancy is based on the top ranked model for hoary bats. Gray shading indicates the 95% credible interval.

https://doi.org/10.1371/journal.pone.0206857.g003

The top ranked model for PESU was Ag + Dev + Forest + Qua (Table 3). Probability of occupancy for PESU was negatively associated with increasing agricultural cover and positively associated with increasing developed land, forest cover, and quaternary road length (Table 1). However, none of these effects were statistically significant (Tables 1 and S4).

We found support for the Pri + Sec + Qua model for TABR (Table 3). Occupancy probability was positively associated with increasing lengths of all road classes (Table 1). Increasing lengths of secondary and quaternary road classes had stronger positive effects on occupancy than the primary road class. However, none of these effects were significant (Tables 1 and S4).

Based on the top ranked model(s), estimated mean probabilities of occupancy and turnover rates varied among species (Fig 4 and Table 4). MYLELUSE had the lowest estimated mean probability of occupancy (0.12) and NYHU had the highest estimated mean probability of occupancy (0.96; Fig 4). Both top ranked models for NYHU produced the same estimates of occupancy. Occupancy probabilities were lower for EPFULANO, LACI, MYAU, and NYHU in 2016 compared to 2015 and higher for DAIN, MYLELUSE, PESU, and TABR; however, estimates did not significantly differ between years for any species (Fig 4). Turnover rates ranged from 0.02 for NYHU (both top-ranked models) to 0.55 for MYLELUSE (Table 4) and 95% credible intervals were narrow for species with low turnover rates and wide for species with higher turnover rates (Table 4).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 4. Estimated mean probabilities of occupancy of each species across South Carolina, May-July 2015 and 2016.

Estimates are based on the top ranked occupancy model for each species. NYHU had two top-ranked models; NYHU 1 refers to the forest edge density (F.ED) model and NYHU 2 refers to the Stream + F.ED model. Blue bars indicate 95% credible intervals. Refer to Table 1 for species code definitions.

https://doi.org/10.1371/journal.pone.0206857.g004

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Table 4. Estimated turnover rates of species occupancy across South Carolina from May-July 2015 to May-July 2016.

https://doi.org/10.1371/journal.pone.0206857.t004

We found that predicted distribution maps based on top-ranked occupancy models differed among species, but closely matched 2003 known ranges (Fig 5). DAIN and MYLELUSE, species for which Region was the top-ranked model, each had high probabilities of occupancy in ecoregions that are completely within their 2003 known ranges, and neither species had a predicted occupancy greater than 30% outside these regions. In 2003, the known distribution of NYHU was statewide and both models predicted occupancy rates greater than 90% statewide. The known distribution of TABR was also statewide, and the model predicted occupancy rates greater than 90% in most areas, except in areas with fewer roads, but rates were consistently above 50%. Although PESU known distribution was statewide, predicted occupancy rates were lowest, down to 15%, in areas with proportionally high agricultural land cover, and high throughout the rest of the state. LACI predicted occupancy rates were highest in the Blue Ridge region, which fully encompasses their 2003 known range. However, LACI occupancy was also high in much of the Southern Coastal Plain region and areas of other regions where stream length was high and forest edge density was low.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 5. Predicted distribution maps for bat species across South Carolina.

Distributions are based on effect estimates in the top-ranked occupancy model for each species, if non-null, and measures of environmental covariates in each cell. Black-outlined squares indicate cells where species were detected in 2015, 2016, or both years. Known summer ranges are based on Menzel et al. [43]. Refer to Table 1 for species code definitions.

https://doi.org/10.1371/journal.pone.0206857.g005

Discussion

Similar to the diversity of threats facing bat populations, our findings highlight the diverse, and often species-specific, environmental factors that predict bat occupancy. In the southeast, where there is ongoing rapid urban and agricultural expansion [27], our results highlight the need to understand which factors are important to the occurrence of each species in designing appropriate conservation and management plans. In addition, WNS was only known to occur in a small northeastern portion of our study area during our research, but has expanded farther across the state and the southeastern US (https://www.whitenosesyndrome.org/). Our data can provide a critical baseline for future analyses of how land use changes influence bat distribution and habitat use over time [9].

Our results suggest that for some bat species, existing U.S. Level III Ecoregion delineations [40] incorporate many of the features we predicted to affect occupancy. Despite the Blue Ridge region making up only a small part of our study area and having only one priority cell (which led to a very wide 95% credible interval for the estimate in this region, Fig 2), we observed strong region-specific patterns in species occupancy and distribution for some species. For example, as we expected based on previous historical range maps [43], we found the highest probability of occupancy for DAIN in the Southern Coastal Plain region. By contrast, for MYLELUSE, we found a significantly higher mean estimated probability of occupancy in the Blue Ridge region than the Southeastern Plains region, as we expected based on their known occurrence in the Blue Ridge region [43]. However, we detected MYLELUSE in three of the five regions (although they were not detected in many cells throughout each region), including one detection in the Coastal Plain. Although this detection was far outside their known historical range [56], these findings are consistent with a similar study in North Carolina [63]. The low detections of MYLELUSE and detections across the state likely led to the relatively low predictive performance (AUC = 0.62) of our model and suggests further research is needed on these species. In contrast, we did not find support for a predicted effect of Region on LACI and MYAU occupancy. LACI appear to have more widespread summer distributions within our study area than expected (Fig 5) as has also been found in North Carolina [63]. MYAU were not detected in many cells within the regions they are thought to occupy, possibly due to their preference for forested wetlands [16,60] that occur in small isolated patches and therefore were not heavily surveyed. MYAU is more of a habitat specialist and typically roosts in tree cavities in wetland habitats and forages in riparian areas [16,60]. Due to its specific habitat requirements, measures of forested wetlands and stream lengths at the cell level may not have been appropriate for modeling MYAU occupancy, and a finer scale analysis may be more suitable. Thus, ecoregion may be a useful tool for determining management activities for some bat species, but not others, especially if species occur in regions in which they were not previously found or have specialized habitat requirements within the ecoregions in which they occur.

While our approach was sensitive for large-scale differences among ecoregions, it likely was not sensitive to site-specific habitat conditions that drive fine-scale habitat selection [35]. For example, apart from region, the only habitat covariate retained in our top candidate model set and for which parameter estimate confidence intervals did not overlap zero was a negative association between forest edge density and LACI occupancy. This selection pattern contradicts other studies that found LACI prefer to forage near forest edge and openings [14,64], but others have observed LACI select for forest interior [65] or exhibit no strong preference between forest edges and opening interiors [66]. Our observation must be considered within the context that LACI exhibit migratory behavior and the majority of individuals may move north, out of our study area during summer [42]. Further, LACI are high flyers [67] and do not always echolocate when they are commuting [68,69]. This may explain why we found a very low mean probability of detection (4.0%) [56], relatively high turnover rate (0.42; Table 4), and a low predictive performance for the top-ranked occupancy model. Although LACI occupancy estimates were relatively high (averaged about 49%), the higher turnover rate may be a further indication of transient individuals opportunistically using habitat as well as low detectability. Further investigation of LACI summer habitat use (e.g., radio tracking of individuals) may be needed to explain the potential effects of forest edge density on their habitat use. The null model likely ranked highest for EPFULANO because it occurs statewide [16], so landscape scale environmental covariates were not able to explain slight differences in occupancy. Although environmental models ranked higher than the null model for NYHU, PESU, and TABR, none of the effects were significant. These species, along with EPFULANO, occupied most of our study area and used a variety of habitats, which could explain why we failed to find significant effects of any landscape-scale environmental factors for these species. In a similar occupancy study in Missouri, Starbuck et al. [70] found greater effects of habitat conditions at the landscape scale than at more localized scales. However, they analyzed occupancy and assessed effects of habitat conditions at each survey point, while we analyzed occupancy and assessed habitat conditions at the scale of 100 km² cells, and landscape factors that affect bat occupancy at points near cell edges may not be represented by cell-level metrics. Thus, conducting multi-scale analyses, with landscape occupancy at the cell level and site use at the level of stationary points and mobile transects, accompanied by more localized measures of habitat conditions, may produce results that better reflect bat species roosting and foraging preferences [35].

While we were not always able to gain insight into individual habitat conditions driving bat occupancy, our models were still useful and insightful for assessing distributions and will be valuable for assessing long-term patterns in landscape-level changes in bat populations. For instance, models for DAIN and MYLELUSE predicted high occupancy rates in regions where they were known to occur. Additionally, NYHU and TABR generally had high predicted occupancy rates statewide and are typically found throughout the state. PESU are positively associated with forest cover [18] and their predicted distribution was lower in areas with low forest coverage and high agricultural coverage. Many of the species which had models with low predictive performance also had low probabilities of detection [56] (Table 3), which could be an indication of biased AUC estimates due to false negatives (i.e., non-detection of species where they were actually present) [71]. Additionally, our data may not have been suitable for assessment with k-fold cross validation and AUC because we had a small sample size (i.e., 38 cells), which was divided into subsets of 25 cells for training data and 13 cells for testing data. Our top-ranked models for most species may therefore be sufficient for predicting landscape occupancy and could be used to guide future mist netting efforts and updating species range maps, even in cases where we determined covariate effects were not significant and found models had low predictive performance.

With longer-term monitoring, data from our two-year study can be utilized in future studies or serve as a baseline for additional analyses to assess landscape-level patterns and change in bat communities in the southeastern US. Overall, populations of some species may be declining in our study area [55] but it was not evident in our study at the landscape scale, likely due to only having two years’ of data. However, other studies at similar spatial scales detected changes in bat populations over time, but these studies were at longer temporal scales (eight to 15 years) [7,9,72]. Thus, if the monitoring we initiated is continued, it is likely that managers will be able to better detect impacts of WNS and other threats to bat populations in our study area [73]. Estimated turnover rates appear to be related to species detection probabilities, where species with high detection probabilities had low turnover rates, averaging about 5%, and species with lower detection probabilities had higher turnover rates, averaging about 40% [56] (Table 4). For example, high turnover rates and low probabilities of detection for DAIN and MYLELUSE could indicate false negatives (i.e., non-detection where species were actually present) each year. The higher turnover rate for LACI may be related to our potential detection of transient individuals, due to their migratory behavior [42]. However, species with higher turnover rates also had wide credible intervals, so these results should be interpreted with caution, and further study is needed.

Overall, by establishing standardized, statewide acoustic monitoring of bats, we found ecoregion and landscape-scale effects on bat occupancy, and provided baseline data which may be used to further evaluate bat species distributions and occupancy at a landscape scale in the face of WNS and rapid land use change. Future studies could utilize our predicted distribution maps and incorporate our data to conduct analyses at various spatial and temporal scales to potentially reveal additional effects of land cover variables on bat species occupancy, and changes in bat populations over time if monitoring is continued [8,73]. Results of our study and further analyses of our data (e.g., with finer scale habitat measures, and more years of data) can therefore increase ecological knowledge of bats and be used to inform conservation priorities [7–9], which is critical to the sustainability of bat populations due to the numerous threats they currently face.

Supporting information

Acknowledgments

We thank Mary Bunch for her collaborative efforts, which helped initialize and support this research. We thank the South Carolina Department of Natural Resources, U.S. Forest Service, U.S. Fish and Wildlife Service, U.S. Army Corps of Engineers, South Carolina Forestry Commission, South Carolina Department of Parks, Recreation, and Tourism, Folk Land Management, Naturaland Trust, Greenville Water, and other agencies, organizations, and volunteers for allowing access and otherwise helping to collect field data. We particularly thank Johnny Stowe, Rob Harrison, Donald Cockman, Daniel Belken, Daniel Jones, and Ashley Watford for devoting their time to collecting data. We also thank Yoichiro Kanno and Catherine Bodinof Jachowski for advice and comments on an earlier draft and M. Kalcounis-Rueppell and an anonymous reviewer for helpful comments. Use of trade, firm, or product names is for descriptive purposes only and does not imply endorsement by the U.S. Government.