Prioritization of Vulnerable Species Under Scenarios of Anthropogenic-Driven Change in Georgia's Coastal Plain

Researchers expect future environmental change to have substantial and irreversible effects on many ecosystems around the world in the coming century. Predictions for global sea level rise (SLR) suggest an increase of up to 2 m by 2100, which would affect ecosystem factors such as tidal range, elevation, and proportions of total brackish, freshwater, and salt marsh (Nicholls et al. 1999; Craft et al. 2009; Kirwan and Megonigal 2013; Watson et al. 2017; Hansen et al. 2016; NOAA 2015). For higher-elevation habitats that SLR is less likely to directly inundate, estimates of future urban growth indicate that development could threaten ecosystem biodiversity and productivity through habitat loss or conversion (Swenson and Franklin 2000; Seto et al. 2012). Effective management of wildlife populations benefits from planning for the relative vulnerability of species to these anthropogenic factors, particularly for threatened or endangered species already at risk of potential extinction (Daniels et al. 1993; Benscoter et al. 2013; Walls et al. 2019). A better understanding of relative vulnerability will help managers to prioritize species and habitats at greatest risk from multiple stressors (Given and Norton 1993; Miller 2006; Gauthier et al. 2010; Barrett et al. 2014; Walls et al. 2019).

In defining criteria used to determine priorities, frameworks that include exposure to change, sensitivity to change, and resilience to change can help practitioners integrate multiple dimensions of vulnerability (Turner et al. 2003; Williams et al. 2008). However, estimating the measures used for evaluating some of these dimensions—such as sensitivity and resilience—requires information on life history, physiology, and adaptive capacity that is often lacking, particularly for rare species (Williams et al. 2008). Species distribution models (SDMs), which provide spatial predictions of species habitat suitability (and thus potential distribution) by relating species occurrence data to environmental variables across a landscape, provide a means to estimate exposure risk to change as one indicator of vulnerability, without requiring more detailed information on species traits (Elith et al. 2006; Elith and Leathwick 2009; Dawson 2011). Species are likely to vary in their exposure to different threats; for example, Brittain and Craft (2012) found that some coastal-dwelling avian species experiencing habitat loss from development were unlikely to experience similar habitat loss from SLR. Thus, effective prioritizations will include an understanding of the severity and form of exposure to potential habitat change that species may experience.

The coastal plain ecoregion in the southeastern United States is a global biodiversity hotspot, noted for high species endemism coupled with high vulnerability to declines in species richness (Myers et al. 2000; Noss et al. 2015). Since the early 20th century, widespread clearing of natural lands for agriculture and development has modified much of the region, resulting in habitat fragmentation and loss, and continued urbanization is likely as human populations in this area are predicted to increase by roughly 30% by 2030 (Turner and Ruscher 1988; Ross et al. 2006; Reece and Noss 2014). Researchers also expect SLR to influence the severity of habitat loss via heightened tidal inundation in the region (Dahl et al. 2017). In Georgia's coastal plain, evidence suggests that these factors are likely to negatively impact species' habitat throughout the coming century (Craft et al. 2009; Georgia Department of Natural Resources 2015; Hunter et al. 2015). A previous study characterized the habitat vulnerability of 24 Georgia coastal plain species to SLR but did not evaluate multiple types of change scenarios (Hunter et al. 2015); thus, researchers have done little work to quantify the vulnerability of listed species in this region to simultaneous threats.

In this study, we used SDMs for 15 species of conservation concern to map habitat suitability across the coastal plain in Georgia. We then assessed vulnerability to multiple threats by projecting the SDMs onto theoretical future landscapes after applying scenarios of urbanization and SLR. Finally, we ranked species vulnerability for each scenario and under multiple ranking schemes considering total exposure to potential habitat change. In addition to exposure, species prioritizations should consider the global and regional imperilment status of species, as well as the adequacy of existing protected areas and other conservation measures. Therefore, we included an evaluation of the fraction of predicted habitat in protected areas as well as the global conservation status of each species. Our objectives were to 1) build SDMs for 10 species and utilize previously published SDMs for 5 additional species to map predicted habitat suitability for a total of 15 species of concern in Georgia, 2) assess species vulnerability via habitat exposure risk using scenarios of potential habitat loss from SLR and urbanization, and 3) evaluate how species prioritization for regional management action may change under these different regimes, and contextualize these prioritizations around current available protected habitat and global vulnerability.

As is common with rare or enigmatic species, sufficient data for mechanistic distribution models were not available; therefore, we constructed SDMs using occurrence records correlated with environmental variables (Elith et al. 2006; Williams et al. 2008; Elith and Leathwick 2009; Dawson et al. 2011). We followed the standards proposed by Araújo et al. (2019) such that most processes and outputs met a standard level of bronze or higher (Araújo et al. 2019). Where processes did not meet recommended community standards, we offer rationale for their omission (see Table S1, Supplemental Material). We consulted with wildlife experts from federal and state agencies, universities, and nonprofits to compile an initial list of approximately 50 target species considered of conservation interest in Georgia. From this list, we eliminated species with a Georgia Conservation status of S5 (currently stable), those species with fewer than 20 occurrence records (Stockwell and Peterson 2002), and species that lacked essential information on range or habitat preferences. Our final list included 15 avian, reptile, and amphibian species (Table 1). For 5 of the selected species (Table 1), SDMs using similar methods were already available from a complementary project (Crawford et al. 2020). We clipped the habitat suitability maps generated from the SDMs in Crawford et al. (2020) to our study region and used this output to examine habitat exposure vulnerability for these 5 species; for the remaining 10 species, we built SDMs and generated habit suitability maps following the methods described below.

Our study area comprised the combined known range of all 15 species within Georgia's coastal plain ecoregion. The ecoregion (Figure 1) is the largest geographical portion of the state, extending from the fall line in the north to the Atlantic Ocean. The region is marked by a wide variety of habitat types such as well-drained sandhills, xeric longleaf pine Pinus palustris and wiregrass Aristida stricta forests, salt marsh, bottomland hardwood swamps, and maritime forest. To better reflect the coastal range of several of our species, we divided the study extent into two parts (Figure 1). For species found exclusively in the lower coastal plain, we restricted the extent to coastally influenced areas and included the nearshore Atlantic Ocean as a land cover type. For all other species, we included the full extent of the ecoregion, as well as a 2-km spatial buffer extended above the fall line that encompassed our largest maximum biological window size (detailed below).

We completed all spatial analyses using ArcGIS version 10.6.1 (Esri 2019) and R version 3.2 (R Core Team 2019). For each species, we compiled a geospatial database of occurrence records denoting point locations where observers species had recorded species presence or noted them absent; the spatial extent for each species covered the extent of the species-specific range in our study area. We used three primary sources to collate our occurrence records: 1) the Georgia Department of Natural Resources Element Occurrence (EO) database, 2) the eBird citizen science data repository, and 3) various research and monitoring studies containing occurrence records for species in our study region.

We used the EO database as the primary data source for a majority of our 10 species. The Georgia Department of Natural Resources' EO database is in accordance with NatureServe guidelines for species data at the state level and consists of records from past and ongoing research studies that are maintained and updated annually (Georgia Department of Natural Resources 2013; NatureServe 2019). For each EO record, we used spatial polygon data to denote areas where a species or natural community is or was historically present within the appropriate seasonal range, originating in date from early 20th century (historical records) to present day. To avoid potential issues with older historical records and to ensure that records reflected recent land cover data, we eliminated all records of collections prior to the year 2000 as well as records with a precision of less than 500 m. We converted polygon records to point coordinates by taking the center point to represent the estimated coordinates at which a species was present.

We also supplemented our avian EO datasets with records from the eBird citizen science data repository (Sullivan et al. 2009). Regional filters help reduce the likelihood of misidentifications and double counts, but there is potential for spatial and temporal biases in the datasets. We followed the eBird recommended best practices for additional bias correction (Johnston et al. 2019). To capture observations that reflect recent updates to the eBird filtering protocol, we used only records from 2010 to 2019. For two of our species (the painted bunting Passerine ciris and the wood stork Mycteria americana), use of Georgia's coastal plain is typically associated with seasonal breeding behavior (Gaines et al. 1998; Springborn and Meyers 2006). To avoid using eBird sightings that may have occurred during migration (and thus do not necessarily reflect true habitat preferences within their species-specific range), we restricted the data to known breeding range and eliminated observations that occurred outside of breeding season for these species. All other avian species present in this study are year-round occupants of Georgia's coastal plain; thus, eBird records reflected occurrence records across all seasons. For all species, if two or more observations were located at a single coordinate set—which can occur due to multiple recorded sightings of single individuals—we randomly selected a single observation to eliminate duplicates. As an additional measure of filtering recommended as part of the eBird best practices, we used hexagonal grids with a 5-km spacing between grid centers d to randomly subsample the remaining points, so that the final output was one observation per 30 × 30 m cell.

For two of our species, we supplemented our EO data using individual research and monitoring studies. For the diamondback terrapin Malaclemys terrapin we used presence data from a multiyear (2008–2018) seining and drone survey in Georgia's tidal creeks and streams conducted during breeding seasons (Grosse et al. 2011; Maerz, unpublished data). We randomly selected a single point along creeks that had recorded diamondback terrapins (He and Gaston 2000). Data for the seaside sparrow Ammospiza maritima came from a multiyear survey of salt marsh bird distributions in coastal Georgia (Hunter et al. 2017).

To minimize the potential for spatial bias that can arise from clustered records, we applied a filter over all occurrence records that randomly removed records occurring within a species-specific biological window of each other (Veloz 2009; Boria et al. 2014). Biological windows (neighborhoods) reflect the relationship between an organism and its surrounding landscape at a specified scale. We based the neighborhood size used to filter records on the average value of each species' core territory obtained from the literature (Table 2). Thus, the final set of records for each species had no more than one observation within the specified neighborhood radius of any other observation. Because of a lack of true absence data for the majority of our species (Table 2), we created sets of pseudoabsence points to compare the environment of known occurrences to background environments (Engler et al. 2004; VanDerWal et al. 2009). We randomly generated pseudoabsence points for each species so that all points 1) fell outside of the predefined neighborhood radius of presence points and 2) within the study area appropriate for each species. We generated points at a 1:4 ratio of presence : pseudoabsence points, so that after the removal of points within neighborhood radii, the final output consisted of presence points and pseudoabsence points at an approximate ratio of 1:3 (following Barbet-Massin et al. 2012). The final set of occurrence data for each individual species was a spatial layer of point locations of filtered presence and absence or pseudoabsence data from the individual or combined sources across the species-specific spatial extent.

For each species, we tested a suite of biotic and abiotic predictor variables hypothesized to influence species distributions, based on a literature search of each species' habitat preferences. We used 30-m raster datasets describing characteristics across the species-specific extent (Table 3). Our hypothesis-based set of variables captured each species' relationship to factors associated with 1) land cover, 2) vegetation characteristics, 3) topography and soil, 4) disturbance, and 5) climate (Table 4). Because species may use habitat differently at various spatial scales, it can be useful to investigate relationships between species and landscapes using multiple neighborhood sizes (Addicott et al. 1987; Johnson et al. 2005; Hagen-Zanker 2016). We selected a minimum and maximum neighborhood size for each species based on best available information regarding species habitat use at different scales and calculated the mean value for each variable at each scale.

We created landscape metrics for each variable using FRAGSTATS version 4.2 (McGarigal et al. 2012), the Spatial Analyst Toolbox in ArcMap, and the SpatialEco package in R (Evans 2020). For land cover, we created variables describing appropriate vegetation, barren or urban land, and wetlands or water factors related to individual species' habitat preferences using present-day habitat types included in the United States Geographical Service (USGS) 2016 National Land Cover dataset (NLCD 2016), the Sea Level Affecting Marshes Model (SLAMM Version 6.6 [Clough et al. 2010]), the Slope Land-use Excluded Urban Topology Hillshade Growth Model Urban Growth model (SLEUTH GA [Clarke et al. 1997]), the USGS National Wetland Inventory database, LANDFIRE Existing Vegetation Type datasets (LANDFIRE 2018, and USGS National Hydrography dataset. Each individual dataset was reclassified as 1, meant to define the variable of interest, or 0, encompassing all other land cover classifications. We then used these datasets to extract metrics such as percentage of habitat, mean habitat patch area, number of habitat patches, and habitat edge density within the specified neighborhoods for each species and scale. For certain variables, we defined unique components of the landscape relevant for individual species. For example, to capture the tendency of the Wilson's plover Charadrius wilsonia to choose nesting sites closer to marsh vegetation, we created a variable that captured the percentages of beach and dry land within 100 m of adjacent marsh (Derose-Wilson et al. 2013). We drew vegetation characteristics from the LANDFIRE Existing Vegetation Height dataset and the National Land Cover Database Tree Canopy Cover dataset. Existing vegetation height was used to describe average heights of vegetation relevant to species. We reclassified heights from 0 (shortest) to 1 (highest) in 0.25 increments, with each increment describing a range of heights. We also included a deciduous index, which captured the difference in winter and summer greenness. Topography and soil variables captured relevant characteristics of elevation and drainage. We used elevation data from the USGS 30-m Digital Elevation Model. For several species, we used a Topographic Position Index (TPI) to represent a location's elevation relative to its local surroundings, including components such as sandhills (positive TPI) and valleys (negative TPI). A soil drainage index from the gridded SSURGO raster from the Natural Resources Conservation Service (NRCS 2017) described poorly-drained to well-drained soils. For disturbance variables, we used a fire frequency dataset representing the percentage of years burned between 2001 and 2016 and a historical land disturbance dataset created from historical land cover data identifying places classified as developed or converted to agriculture between 1938 and 2001 (Crawford et al. 2020). Because air temperature in Georgia's coastal plain has limited spatial variability (and thus explanatory power), we did not include gridded temperature data as a covariate. Datasets representing the mean summer and winter precipitation were included for several reptile species. We extracted values from the gridded datasets for each covariate to the presence and pseudoabsence points for model fitting.

We used logistic regression in a model-selection framework to create our presence and absence or pseudoabsence SDMs (Burnham and Anderson 2002), thus predicting relative suitability across each species' range. We completed all models and accompanying statistical analyses in R, version 3.5. For each species, we first grouped all variables by their type (e.g., beach or flat for the American oystercatcher Haematopus palliatus; Table 4) and performed model selection to compare neighborhood sizes and choose the scale most appropriate for each variable; we selected models using the Akaike Information Criterion (AIC_c) weights to select the model with the highest weight. We then tested all variables from the first stage of model selection for collinearity by using a Pearson correlation coefficient ≥ 0.7 or ≤ -0.7 to evaluate pairs of variables that were correlated, choosing the best-supported variable of the two. Finally, for each species, we tested all remaining variables at their best-supported scales, including quadratic and 1-way interaction terms for variables where the literature supported quadratic or interaction relationships, and chose the best-fitting model again using the highest AIC_c weight.

To examine the performance of the chosen models, we used fourfold cross-validation to conduct model evaluation using several methods (Araújo and Guisan 2006). We calculated the area under the curve of the receiver-operating characteristic plot, the classification accuracy of presence and pseudoabsence points (using the True Skill statistic [TSS] − sensitivity [proportion of true presences] + specificity [proportion of true negatives] − 1), and the point biserial correlation coefficient (Fielding and Bell 1997; Liu et al. 2011). We also produced variograms using model residuals to test for spatial autocorrelation and calculated the relative importance of each variable using hierarchical partitioning, which measures proportion of variance explained by each component of the model (Table S1, Supplemental Material; Chevan and Sutherland 1991). All best-fitting models from the initial model selection process met our evaluation criteria; thus, the final set of SDMs reflected the original suite of best-fitting models. We projected the best-fitting models (Table 5) out to the landscape to map continuous habitat suitability ranging from 0 (not suitable) to 1 (highly suitable). To evaluate fraction of potential exposure to habitat loss, it is necessary to select a threshold from this continuous range that denotes habitat as binary classes of suitable or unsuitable (Bean et al. 2012). We converted continuous suitability into binary maps of 1 or 0 (suitable or unsuitable) using the optimal threshold value based on the maximum TSS value; that is, the value of probability of species presence that—if used to classify the landscape into suitable and unsuitable patches—results in highest classification accuracy for presence and pseudoabsence (Allouche et al. 2006).

As model predictions always contain some error, it is desirable to remove any error for which the source can be determined. In our case, some models included cells containing either open water or high-intensity urban (or both) as potentially suitable habitat. This can occur because variables included in the best model for a species may not discriminate with respect to those characteristics (e.g., percentage of landscape in particular cover type or elevation). In these cases, we converted those cells to the unsuitable class. For coastal regions, we defined open water as all classes encompassing inland open water, riverine open water, estuarine open water, and open ocean using the SLAMM classes for Georgia's coast. In upland areas where SLAMM data is not available, we used the National Land Cover Database Open Water category. To define high-intensity urban areas, we used the SLEUTH urban growth model, which denotes current urbanization based on land cover, transportation, and topography. We created datasets describing currently protected habitat for each species by extracting the binary SDM outputs to the USGS Protected Areas Database (USGS 2019), and the Georgia Department of Natural Resources Conservation Lands Database (Georgia Department of Natural Resources 2019). These databases provide an inventory of lands both publicly and privately owned and managed that are registered in state or federal programs where conservation is a management goal (Crawford et al. 2020). We calculated present-day area (km²) and percentage of both total available habitat and habitat within protected lands.

We used the SLAMM predictions for coastal regions to define inundation from SLR and the SLEUTH urban growth model to represent urban growth in the future. We chose four future time points at which to examine the impacts of both SLR and urban growth: 2025, 2050, 2075, and 2100. For urbanization, SLEUTH defines raster classes denoting predicted probability of growth based on data describing slope, land use, exclusion, urbanization, transportation, and hill-shade (Clarke et al. 1997). Classes are numbered from 3 (0–2.5% probability of urban growth) to 16 (97.5–100% probability of urban growth). We included all classes greater than or equal to 10 (50–60% probability) and reclassified them to a single class to represent predicted urban growth for each future time point. To predict coastal change from SLR, SLAMM uses digital elevation data and National Wetland Inventory data to simulate processes involved in wetland conversion under different scenarios of SLR, resulting in a dataset representing changed land cover conditions. We chose a SLR scenario of 2 m of rise by 2100. We selected this representative value because 1) tests of 1-m and 1.5-m scenarios revealed no difference in rankings from the 2-m scenario and 2) recent studies suggest that higher scenarios of SLR may be more realistic given the current trajectory of global temperatures (Bamber et al. 2019; Kulp and Strauss 2019).

We reclassified SLAMM classes denoting riverine open water, estuarine open water, and open ocean to a single class for inundation under the 2-m SLR scenario. To assess the impact of future SLR and urbanization on coastal plain species, we used several measures of change in total available habitat and protected area habitat to define potential habitat loss. We first overlaid the binary SDM outputs for the entire extent and within protected areas for each threshold with the binary datasets conveying urban growth and inundation in R. We converted each SDM cell that overlapped cells denoting urban growth or inundation to “unsuitable,” so that the output was a raster with a changed sum of total habitat cells, indicating species exposure to potential habitat loss at each future time point. We used these outputs to calculate percentage of change, percentage of total, and area of range-wide and protected habitat for each scenario. Because scenarios of urban growth are unlikely to impact currently protected habitat, we only evaluated future habitat loss from inundation for protected areas.

We used several approaches to rank species for our final set of prioritization schemes. We first ranked species by their percentage of exposure to potential habitat change as a result of SLR and urbanization, where ranks of 1 represented highest percentage of potential habitat loss or change. We also ranked species by area of total available habitat using 1 to represent least amount of area (and thus top priority). To assess the amount of habitat falling within conservation areas and thus offering some potential protection via state or federal jurisdiction, we also ranked species by the total percentage of available habitat within protected lands, where 1 represented the lowest percentage of protected habitat. We compared these ranks for urbanization to SLR, choosing a cutoff rank of 8 to represent species within distinct categories of vulnerability for each type of change. Finally, we used Global Rank, a metric used for conservation status assessment created by national conservation nonprofit NatureServe (NatureServe 2019), to convey species' relative global vulnerability to potential extinction. We used a cutoff rank of 3 (“vulnerable”) to represent species above or below high global vulnerability.

Our filtering methods resulted in ≥ 49 total presence points for each species, consistent with recommended sample sizes for SDMs (Table 2). Best-fitting models exhibited adequate to good model performance (Table 6). Area under the curve values of 0.7 or higher are considered reasonably good at distinguishing presence from pseudoabsence; our values ranged from 0.77 (adequate) to 0.95 (excellent). Biserial correlation coefficients were above 0.30 for all models (Elith et al. 2006), with a minimum biserial correlation value of 0.44 and a maximum of 0.84. Accuracy (the maximum TSS) was above 65% for all models. Residual variograms revealed no concerning patterns of residual spatial autocorrelation. The red-cockaded woodpecker Dryobates borealis, and Bachman's sparrow Peucaea aestivalis displayed some evidence of spatial structure, but the range of spatial correlation was at a fine enough scale that no bias was expected. Details on variable contribution and coefficients for each model are in Tables S2 and S3, Supplemental Material. Detailed maps of predicted habitat suitability are in Figure S1, Supplemental Material. Total area of predicted habitat ranged from 343.3 to 35,798.2 km², using the optimal cutoff threshold (Table S4, Supplemental Material). Habitat falling within protected areas ranged from 107.5 to 3,649.7 km², or 7.3–31.3 % of total habitat.

We present prioritization schemes for species vulnerability to future exposure to potential habitat loss or change, using habitat classified as available under the optimal cutoff threshold. Species restricted to the lower coastal plain ranked highest for exposure to potential habitat loss from SLR by 2100, ranging from 7 to 72% (Figure 2). In contrast, species occupying broader ranges experienced little exposure to potential habitat loss from SLR (0–6%). The highest percentage of potential habitat lost by a species via urbanization was 22%. Species ranked for high vulnerability to exposure from SLR generally ranked low for vulnerability to change from urbanization, and vice versa (Figure 3A). Two exceptions were the painted bunting and wood stork, which ranked moderately high for habitat loss to SLR, and high or highest for habitat loss due to urbanization (Figure 3B). Several species ranked moderately low for vulnerability to both components (Figure 3C). Species ranking highest on the basis of percentage of habitat loss from SLR were also ranked highest based on total predicted area, but tended to rank lower based on the amount of habitat currently protected (Table 7). The top-ranking species under scenarios of urbanization ranked lowest based on total predicted area but tended to have less protected habitat (Table 7). Species ranked lowest (“apparently secure” to “secure”) for global vulnerability (Figure 4A [top panel]), typically had the highest regional vulnerability to SLR, while those ranked for higher global vulnerability (“vulnerable” to “imperiled”) had the lowest (Figure 4D [top panel]). Several species ranked high for global vulnerability also ranked reasonably high for regional vulnerability to urbanization (Figure 4B [bottom panel]). Species ranking low for urbanization vulnerability also ranked low for global vulnerability (Figure 4C [bottom panel]).

Our results indicate that all species evaluated are likely to experience habitat loss or change from either SLR or urbanization by 2100, but that few species will experience relatively high exposure to loss or change from both stressors. The severity of potential habitat loss experienced by species was highly dependent on their range and amount of habitat predicted to be currently suitable. Coastally restricted species had less predicted habitat area initially and experienced more severe declines in predicted habitat in the future than those species occupying some portion of the larger coastal plain. Given the restricted range of both salt marsh and coastline habitats in Georgia, this is unsurprising (U.S. Fish and Wildlife Service 2007). We found that species utilizing salt marsh, beach, and other coastal habitat types will be most vulnerable to potential habitat change from SLR, which is consistent with previous models of salt-marsh- and beach-reliant species (Galbraith et al. 2002; Brittain and Craft 2012; Hunter et al. 2015; Valdes et al. 2016). We also felt that focusing on results due to change by the end of the century was appropriate, as top-ranking species for early (2025) and midcentury (2050) habitat loss deviated little from results by 2100. A majority of species occupying some or all of the entire coastal plain experienced little habitat loss from SLR, averaging losses below 0.5% even for 2-m SLR scenarios. Instead, inland coastal plain species were predicted to experience relatively higher habitat loss (3–20%) from urbanization. These results agree with other studies suggesting development is a primary threat to Georgia's inland terrestrial species in the region (Gibbon et al. 2000; Plentovich et al. 2007; Breininger et al. 2012; Leonard et al. 2017).

Because ranks based on the fraction of protected habitat were a function of total available habitat area, assessing prioritizations based on protected lands was difficult. Several of the lowest ranking species for priority on the basis of protected habitat ranked low for vulnerability to habitat loss, yet still ranked moderately high for action based on total available habitat area. For example, the red-cockaded woodpecker, Bachman's sparrow, and gopher frog Lithobates capito all had roughly one-third to half of their total habitat falling within existing protected areas yet ranked within the top 10 for priority based on low total available habitat. Researchers have historically documented Georgia populations of these species in longleaf pine ecosystems (Plentovich et al. 2007; Maerz and Terrell 2016). Due to the widespread decline of these ecosystems and the recognition of their importance for multiple regional species, a large portion of remnant or restored longleaf pine habitat in Georgia tends to fall within some protected land. This indicates that while species utilizing these landscapes may have much of their available habitat currently secure, they still rank relatively high for action based on total available area due to specific habitat requirements. This pattern was also true for species that use salt marsh and beach habitats. Those species had one-fourth to one-third of their habitat protected along Georgia's coast, yet have little total predicted habitat to begin with (e.g., seaside sparrow, American oystercatcher). Further, protected or not, these habitat types will be vulnerable to the effects of SLR, whereas inland protected areas are unlikely to be vulnerable to development. This suggests that while managers may be able to use the amount of habitat falling within protected area to offset the impacts of future change, giving lower prioritization to species solely because they rank high on the basis of proportion of total habitat protected may be an ineffective strategy, and managers could adapt more holistic prioritization schemes by considering the amount of both protected and total available habitat for a species.

All species ranking highest for exposure to habitat change from SLR were classified as “apparently secure” or “secure” globally, indicating these species reportedly have large to medium populations currently showing no extreme declines throughout their range (Clay et al. 2014; NatureServe 2019). The challenge in weighing global status against regional vulnerability is that the magnitude of species' vulnerability to potential habitat loss within a region may outweigh priorities for global conservation, particularly when global priority schemes may be data deficient. For many of our focal species, researchers have not reevaluated population status since 1996 (NatureServe 2019). This means that global listings may not be accurate reflections of range and population resiliency, a problem documented in several large-scale global priority systems (Ramesh et al. 2017). This does not necessarily mean that global conservation status should be discounted as a tool for regional management decision, but rather that managers may need to use global conservation status as a secondary measure for final priority ranking if the magnitude of projected regional habitat loss for a species is extremely high (e.g., 35–100% of habitat projected to be lost). Although species ranking high for both global vulnerability and regional vulnerability to habitat loss (e.g., the Southern hognose snake Heterodon simus and the striped newt Notophthalmus perstriatus) will still likely rank high for priority for conservation action, managers will need to weigh the magnitude of regional vulnerability to habitat loss as a primary factor for those species ranking lower globally.

While all salt marsh- and beach-dwelling species ranked high based on potential habitat loss, we projected the seaside sparrow to experience particularly high potential habitat loss. Wildlife managers consider seaside sparrows to be habitat specialists, relying primarily on salt marsh for both nesting and foraging activities (Hunter et al. 2015). Unlike the other species that ranked top for habitat loss to SLR, all of which utilize beach habitats for nesting, foraging, or both types of activities in conjunction with salt and brackish marsh, the seaside sparrow is limited by its dependence on this singular habitat type for both breeding and foraging activities. Under the SLAMM model, salt marsh habitats are degraded substantially, with a conservative estimate of 6% loss by 2100 under a 1-m SLR. In contrast, the model predicts increases in tidal flats and estuarine beaches, meaning that beach habitats will be surprisingly persistent despite inundation; this likely accounts for the gap between habitat loss values between the seaside sparrow and species ranked directly below it. For upland species, the contribution of urbanization to the conversion of alternative inland habitat types could also inadvertently buffer the impact of habitat loss. This is most applicable in the case of the painted bunting and wood stork, both of which rank moderately high for vulnerability to habitat change from both SLR and urban growth. Both species concentrate around coastal wetlands or swamps (freshwater and brackish or salt marsh) and shrub-scrub and maritime forest habitats for nesting, but are also commonly found near inland agricultural, pasture, or even low or moderately developed areas in and outside of breeding season, consistent with habitat variables used in our top models (Gaines et al. 1998; Kopachena and Crist 2006; Brittain and Craft 2012 . While potentially more susceptible to both forms of stressors due to their use of habitat along the coastal–inland gradient, increases in urban-adjacent landscapes may help to offset the loss of natural landscapes for these species (Kopachena and Crist 2006; Lee and Carroll 2014).

We note two particular challenges with models simulating future SLR and urbanization in this region. First, we chose to present results using a 2-m SLR scenario. Predictions of SLR are variable, but recent studies suggest that a 2-m scenario is more realistic than a 1-m scenario (Kopp et al. 2017; Kulp and Strauss 2019), and our rankings were robust to the choice of scenario. Additionally, the relative rank of species was consistent in less severe scenarios, and we expected that these patterns would persist under even greater scenarios of SLR. Using this justification, we felt our choice of a 2-m scenario was appropriate. However, the magnitude of exposure does vary between scenarios, and agencies may choose to consider more conservative scenarios, potentially resulting in slightly different interpretations of results. Secondly, models for development often fail to account for some human responses to future stressors such as SLR. The SLEUTH model relies on information from current distributions of development and does not presently include information about potential human responses to SLR that may inadvertently impact coastal wildlife populations, therefore potentially underestimating the true impact of future development. For example, the building of seawalls is rapidly becoming a common urban planning technique to address SLR. Researchers have shown seawalls to have negative impacts on diamondback terrapin nesting habitat and to be correlated with reduced occupancy (Isdell et al. 2015; Winters et al. 2015). While efforts to construct seawalls are presently limited in Georgia, it is reasonable to assume that armoring of shorelines will be a mitigation tactic employed in the future, meaning that SLEUTH is presently unable to capture these relationships. Thus, our estimates of exposure to habitat change from development may be conservative for some coastal species that will likely be prone to vulnerability from development-related mitigation efforts.

While the nuances of impacts to species from future anthropogenic threats are difficult to completely assess, we present an initial attempt to help managers consider prioritization schemes under multiple types of change. Our results suggest that managers may need to prioritize species (or their habitats) based on the regional total amount of available habitat, and the projected magnitude of habitat loss. Using a scenario of a 2-m SLR, we projected that species restricted to the lower coastal plain would lose up to 40–70% of their habitat due to SLR, roughly twice the amount of habitat loss projected for top-ranking species under urbanization. Although these species currently rank low for global vulnerability, this result has implications for populations outside of Georgia, as SLR may have similar effects on populations along the rest of the Atlantic coast, potentially changing the nature of species' global status (Hayes 1994). It may also be valuable for managers seeking to develop short-term conservation plans to use priority schemes for timescales closer to midcentury, as several inland species ranked higher for immediate habitat loss due to urbanization by 2025 than species experiencing habitat loss from SLR by 2025. Additionally, efforts to prevent population declines or extinction of rare species may benefit more from long-term planning than short-term triage, meaning that vulnerability metrics for longer time horizons are optimal (Wilson et al. 2011). Differing conservation timelines and goals will mean that managers need multiple lines of evidence in order to make informed decisions and appropriately allocate resources. We offer a multiscenario, broad-ranging set of results that may help to contextualize potential management actions and provide a first step for addressing species vulnerability in coastal Georgia.

Please note: The Journal of Fish and Wildlife Management is not responsible for the content or functionality of any supplemental material. Queries should be directed to the corresponding author for the article.

Figure S1. Maps describing species-distribution model (SDM) spatial projections of predicted available habitat across relevant study areas in Georgia for 10 species, as of 2019. Species-specific results are given for (a) American oystercatcher Haematopus palliatus, (b) Bachman's sparrow Peucaea aestivalis, (c) diamondback terrapin Malaclemys terrapin, (d) eastern diamondback rattlesnake Crotalus adamanteus, (e) eastern indigo snake Drymarchon couperi, (f) painted bunting Passerina ciris, (g) red-cockaded woodpecker Leuconotopicus borealis, (h) seaside sparrow Ammospiza maritima, (i) Wilson's plover Charadrius wilsonia, and (j) wood stork Mycteria americana.

Available: https://doi.org/10.3996/JFWM-20-089.S1 (5.75 MB DOCX)

Table S1. Scores for species-distribution models (SDMs), following guidelines for scoring models outlined in Araújo et al. (2019). Our SDMs fall under the “prediction” category of model use; that is, modeled species-environment relationships mapped to potential distributions in the same time period and geographical region. For all four categories (response variables, predictor variables, model building, and model evaluation), reported scores represent the minimum score for the metric in question for all models created. In some instances, models may have scored higher for certain metrics (e.g., for response variable: spatial accuracy, identification of taxa, etc., may be gold–silver). We offer explanations for any scores that are deemed deficient.

Available: https://doi.org/10.3996/JFWM-20-089.S2 (23 KB DOCX)

Table S2. Proportion of variance (% contribution) explained by individual predictors in top models for each species in coastal Georgia as of 2019, using a hierarchical partitioning scheme.

Available: https://doi.org/10.3996/JFWM-20-089.S3 (24 KB DOCX)

Table S3. 2019 coastal Georgia species-distribution model (SDM) estimates of predictor effects, for all species (with confidence intervals).

Available: https://doi.org/10.3996/JFWM-20-089.S4 (30 KB DOCX)

Table S4. Metrics summarizing present-day characteristics of potential habitat, coastal Georgia, 2019. Total represents the percentage and area (km²) of potential habitat across Georgia's coastal plain. Total PA represents the total percentage and area of potential habitat within protected areas. Bold lettering indicates species with a study extent restricted to the lower coastal plain.

Available: https://doi.org/10.3996/JFWM-20-089.S5 (22 KB DOCX)

Reference S1. Araújo MB, Anderson RP, Barbosa AM, Beale CM, Dormann CF, Early R, Garcia RA, Guisan A, Maiorano L, Naimi B, O'Hara RB, Zimmermann NE, Rahbek C. 2019. Standards for distribution models in biodiversity assessments. Science Advances 5(1):48–58.

Available: https://doi.org/10.3996/JFWM-20-089.S6 (887 KB PDF) and https://doi.org/10.1126/sciadv.aat4858

Reference S2.Georgia Department of Natural Resources. 2015. Georgia State Wildlife Action Plan. Social Circle, Georgia: Georgia Department of Natural Resources.

Available: https://doi.org/10.3996/JFWM-20-089.S7 (16.49 MB PDF)

Reference S3. Maerz J, Terrell V. 2016. Status assessment of gopher frog (Lithobates capito) breeding wetlands and upland habitats in Georgia. Georgia Department of Natural Resources.

Available: https://doi.org/10.3996/JFWM-20-089.S8 (718 KB PDF)

Reference S4.[NOAA] National Oceanic and Atmospheric Administration. 2015. State of Georgia: coastal management program and program document. NOAA Office of Ocean and Coastal Resource Management and Georgia Department of natural Resources Coastal Resources Division.

Available: https://doi.org/10.3996/JFWM-20-089.S9 (1.43 MB PDF)

Reference S5. Ross CL, Barringer J, Benjamin S, Foard A, Leone de Nie K, Pierce D, Wilkins J. 2006. Georgia coast 2030: population projections for the 10-county coastal region. Georgia Institute of Technology.

Available: https://doi.org/10.3996/JFWM-20-089.S10 (438 KB PDF)

Reference S6. Watson A, Reece J, Tirpak BE, Edwards CK, Geselbracht L, Woodrey M, La Peyre MK, Dalyander PS. 2017. The Gulf Coast Vulnerability Assessment: mangrove, tidal emergent marsh, barrier islands, and oyster reef. Forest and Wildlife Research Center, Research Bulletin WFA421, Mississippi State University.

Available: https://doi.org/10.3996/JFWM-20-089.S11 (2.09 MB PDF)

Funding for this research was provided by the Georgia Department of Natural Resources' Coastal Resources Division, via the Coastal Incentive Grant (CIG) Program, as well as the Warnell School of Forestry and Natural Resources. The CIG Program is made possible by a grant through the National Oceanic and Atmospheric Association, pursuant to congressional funding from the Coastal Zone Management Act. Our thanks to Dr. Joshua Reece from California State University, Fresno, and all affiliated lab members who assisted with this work for their dedication to this project, and their role in building a significant portion of the groundwork for this research. We thank all individuals who contributed data and insight that facilitated the success of this project, including Dr. Elizabeth Hunter, Dr. Clint Moore, Greg Krakow, and Brett Albanese. We also thank all agencies and entities involved in data collation, including the Georgia Department of Natural Resources, the Cornell Lab of Ornithology (eBird), NatureServe, University of Georgia Savannah River Ecology Laboratory, Joseph W. Jones Ecological Research Center, U.S. Army (Fort Stewart), U.S. Air Force (Eglin Air Force Base), Longleaf Alliance, The Nature Conservancy, Partners in Amphibian and Reptile Conservation, and HerpMapper. Finally, we thank journal editors and reviewers for their time and effort spent in refining this manuscript.

Any use of trade, product, website, or firm names in this publication is for descriptive purposes only and does not imply endorsement by the U.S. Government.