Subterranean movement inferred by temporary emigration in Barton Springs salamanders (Eurycea sosorum)

Introduction

Amphibians are declining in many parts of the world owing to numerous factors (Lannoo, 2005; Catenazzi, 2015; Grant et al., 2016). Studies linking movement patterns to underlying behavioral and ecological mechanisms (Nathan et al., 2008) are needed to better understand amphibian population declines (Pittman, Osbourn & Semlitsch, 2014). Tracking individual movement and fate has been suggested as a way to close this information gap by directly linking fitness consequences of movement to management strategies or environmental conditions (Bailey & Muths, 2019). Although appealing, this approach requires radio telemetry (Bailey & Muths, 2019), which is impractical for many amphibians because they are too small to be fitted with transmitters. Indirect assessments of movement behavior are often more practical for many amphibians (e.g., capture-recapture studies; Marsh et al., 2005; Lowe, Likens & Cosentino, 2006; Ringler, Ursprung & Hödl, 2009; Cayuela et al., 2020), including fossorial or cave-dwelling species (Balázs, Lewarne & Herczeg, 2020).

Many salamander species use subterranean habitat because it provides environmental conditions that are relatively stable compared to surface habitats. Some species are obligate residents of caves with highly specialized adaptations (troglobites), while others are facultative troglophiles that periodically venture underground (Gorički et al., 2019). Studies of these populations are challenging because individuals are difficult or impossible to observe directly. Because of access restrictions to cave habitats, high levels of endemicity (Culver et al., 2000; Trontelj et al., 2009), and apparent rarity (Krejca & Weckerly, 2007; Jugovic et al., 2015; DiStefano et al., 2020), answering basic ecological questions (e.g., Who lives there? How many are there?) remains a fundamental challenge for the study of cave-dwelling animals (Wynne et al., 2019). Cave-dwelling vertebrates are typically observed only in relatively small portions of their habitat, where humans can enter or remotely access their environment. Sampling generally occurs in human-accessible caves or through “windows” into their environment such as wells or spring outlets (Miller & Niemiller, 2008; Graening et al., 2010; Day, Gerken & Adams, 2016; Niemiller et al., 2016; Phillips et al., 2017; Krejca & Reddell, 2019). Moreover, many cave-adapted species are of conservation concern because of their high endemicity and the sensitivity of subterranean habitats to disturbance and pollution (Chippindale & Price, 2005; Miller & Niemiller, 2008; Fenolio et al., 2013; Devitt et al., 2019; Gorički et al., 2019), underscoring the need for basic research on the biology of subterranean fauna (Mammola et al., 2019).

Studies of cave organisms have relied mostly on capture-recapture data to make inferences about their movement behavior (Means & Johnson, 1995; Trajano, 1997; Lopes Ferreira et al., 2005; Băncilă et al., 2018; Balázs, Lewarne & Herczeg, 2020). Direct observations of recaptured individuals have been used to quantify movement patterns where habitat is accessible (Means & Johnson, 1995; Trajano, 1997; Balázs, Lewarne & Herczeg, 2020). Where animals retreat to inaccessible habitat, e.g., from a cave stream or spring outlet to the aquifer that feeds it, capture-recapture data may indicate immigration to or emigration from these areas (Means & Johnson, 1995; Day, Gerken & Adams, 2016). Explicitly modeling transitions between observable and unobservable states using capture-recapture data (so-called temporary emigration) allows for estimates of demographic parameters (Kendall, Nichols & Hines, 1997), breeding status (e.g., Kendall & Bjorkland, 2001), habitat use (Bailey, Simons & Pollock, 2004a; Cecala, Price & Dorcas, 2013), and movement (Băncilă et al., 2018). A temporary emigration model provides estimates of movement into and out of the sample area based on the transition rate between observable and unobservable states (Kendall, Nichols & Hines, 1997). In cases where observable states correspond to accessible areas (e.g., within a cave gallery or at a spring outlet), quantifying temporary emigration allows inferences about the characteristics of the population that use both the accessible and inaccessible areas rather than being limited to only the area directly sampled.

Barton Springs salamanders (Eurycea sosorum) occur in both surface and subterranean habitat, but populations are typically only accessible for study near springs. When first described, the salamanders were assumed to be restricted to the immediate vicinity of Barton Springs in Austin, Texas (Chippindale, Price & Hillis, 1993), one of several large spring complexes in the Edwards Aquifer of central Texas that harbors endemic species (Bowles & Arsuffi, 1993; Krejca & Reddell, 2019). Because the subterranean conduits feeding Barton Springs are largely inaccessible, biologists questioned whether the species was primarily a surface dweller that occasionally enters subterranean habitat or an obligate subterranean species that is flushed to the surface accidentally (Sweet, 1978; Chippindale, Price & Hillis, 1993). Following its description, in 1997 this species was listed as endangered under the US Endangered Species Act of 1973 due to its exceptionally small range (at the time, known from only three springs in a city park) and presumed low abundances (Chippindale, Price & Hillis, 1993; City of Austin, 1997; US Fish and Wildlife Service, 1997). Although additional occurrences of this species have subsequently been recorded outside of the type locality (Bendik et al., 2013a; McDermid, Sprouse & Krejca, 2015; Devitt & Nissen, 2018), abundances at those sites are lower and surface habitat is more limited in extent compared to Barton Springs. Abundances can be high at the type locality (e.g., >500 individuals) but may fluctuate by orders of magnitude and occasionally reach zero, or nearly so (Bendik & Dries, 2018; Dries & Colucci, 2018). This suggests that either surface populations are occasionally extirpated and recolonized from subterranean populations or that emigration from the surface to the aquifer is occurring (Bendik & Dries, 2018). Complete extirpation without emigration into the aquifer would suggest that the surface habitat is a population sink (i.e., unable to be maintained without immigration; Pulliam, 1988) and that surface populations could be a dead-end in terms of species persistence. If so, changes to population management practices and recovery planning may be necessary, because a primary focus of conservation for Barton Springs salamanders includes management and restoration of the surface habitat (Dries et al., 2013; US Fish and Wildlife Service, 2016).

Here, we report the results of a four-year robust-design capture-recapture study to better understand movement patterns and subterranean habitat use of Barton Springs salamanders. We examined how extrinsic environmental factors (spring discharge, turbidity, dissolved oxygen, temperature) and intrinsic factors (body size) affect survival and temporary emigration. Population dynamics differ between juveniles and adults in Barton Springs salamanders (Bendik & Dries, 2018), and size is often correlated with the propensity for movement in amphibians (Cayuela et al., 2020). Similarly, spring flow is an important driver of environmental conditions (e.g., Mahler & Bourgeais, 2013) that are important to the ecology and physiology of groundwater Eurycea salamanders (Fries, 2002; Woods et al., 2010; Crow et al., 2016; Bendik & Dries, 2018) and may influence movement patterns as well.

Materials & Methods

Study site

We studied Barton Springs salamanders (Fig. 1) at Eliza Spring (Fig. 2), part of the Barton Springs complex in Austin, Texas, USA. Eliza Spring is adjacent to Barton Springs Pool, a popular spring-fed swimming pool. Barton Springs is the primary discharge point for the Barton Springs segment of the Edwards Aquifer, a sole-source karstic aquifer that responds rapidly to changes in precipitation via a system of sinking streams that occur over a 430 km² area (Hunt, Smith & Hauwert, 2019). Eliza Spring emerges into a concrete-bottomed pool (74 m²) with a water depth of approximately 0.3 m (Fig. 2A) and harbors the highest average densities of Barton Springs salamanders among known spring localities (4.32/m² ±  3.4 SD compared to ≤ 0.25/m² ± 0.3  at other sites within the Barton Springs complex; Dries & Colucci, 2018). Throughout most of our study, the spring pool discharged directly into a concrete pipe (0.6 m diameter) that flowed to Barton Creek, likely representing a point of no return for salamanders. During the last year of our study, the pipe was replaced by a 1 m wide overland stream that created 25 m² of additional habitat, but that still terminated at a concrete tunnel that conveys Barton Creek (Fig. 2B).

Results

Our capture-recapture data set included 2,046 observations of 1,578 Barton Springs salamanders, 333 of which were recaptured at least once between periods. Almost all salamander observations occurred within the spring pool (99%) and only 1% occurred within the newly constructed stream (the stream was not fully colonized yet; unpublished data). During the first stage of model fitting, the best model of detection was p = c with full-time variation across all survey events (w =1.00). Model-averaged values of $\hat{p}$ ranged from 0.25 (SE = 0.054) during the second survey of November 2019 to 0.82 (SE = 0.043) during the first survey of November 2017. The most general models with temporary emigration were ranked higher than the non-movement model (ΔAICc = 1065), with Markovian movement favored over the random movement model (w = 1.00; ΔAICc = 42.4). None of the models including environmental covariates were favored (Table 2; sum of w = 0). In fact, our results indicated the effects of environmental variables (DO, spring flow, turbidity) were far outweighed by the effect of body size as a determinant of surface availability and survival (Table 2). Models with time-varying survival and emigration parameters were also more parsimonious than those where temporal variation was explained by environmental covariates (Table 2). The best model overall was S (∼size), γ″ (∼time + size), γ′ (∼time), p = c (∼period:time) (w = 0.35; Table 2). Goodness-of-fit tests are not available for robust-design models, but we note that the order of the top models was invariant to adjustments of the variance inflation factor, $\hat{c}$.

Despite models with environmental covariates not ranking highly, a comparison among models (n = 27) that only included environmental covariates as a function of survival and emigration parameters provides some insight as to the direction of different environmental effects (Table 3). An exhaustive model set of all possible covariate combinations (which is generally not recommended; Burnham & Anderson, 2002b) precludes use of model averaging of the beta estimates. Instead, we provide coefficients (on the logit scale) from the top five models (sum of AICc w = 0.88) within this subset to illustrate the direction of environmental effects (Table 4). Estimates for the top models were generally consistent and did not change sign (Table 4). Flow-temp-turbidity and flow-only models of survival outperformed turb-DO models, while turb-DO models for emigration were favored over others (Table 3). Temperature and flow were positively correlated with survival, which was consistent with our hypothesis (Table 4). Turbidity was positively associated with both emigration parameters, but only significantly so for γ″. DO was negatively associated with both emigration parameters, but again, this affect was negligible for γ′. The remaining coefficients had 95% confidence intervals overlapping zero (Table 4). Collectively, these results suggest that salamanders move below ground with increasing turbidity and decreasing DO, which was consistent with our expectation.

Growth was modeled from 901 measurements of 372 recaptured individuals. Mean parameter estimates with 95% credence intervals from the von Bertalanffy growth model were as follows: a = 31.57 (30.99–32.18), k = 7.81 (7.14–8.52), λ = 4.49 (3.34–6.01) and sdme = 1.25 (1.11–1.39). The bootstrap analysis using random draws from the posterior distributions of parameter estimates from the growth model showed little variation in AICc for the top model (range 10654.32–10655.75), indicating that different potential realizations of the individual growth trajectories (based on the population mean) have little impact on our inference. This is due to our study design sampling intervals being long enough and growth of salamanders fast enough (Fig. 3), that small shifts in growth rate parameters have a negligible effect on the predicted transition rates between size classes at the resolution of our data.

Body size had a positive effect on S and γ″. Estimates of the body size effect on the logit scale (beta) from the top model were 1.43 (SE = 0.17) and 1.45 (SE = 0.28) for S and γ″, respectively. Among both size classes, estimates of monthly apparent survival, $\hat{S}$, ranged between 0.75 (SE = 0.04) and 0.95 (SE = 0.02) (Fig. 4A), which corresponds to annual survival rates between 0.03 and 0.52. Temporary emigration varied among size classes ($\widehat{{\gamma }^{\prime \prime }}$ range: 0.16–0.97); the probability of moving from the pool to subterranean habitat was higher for larger individuals (mean $\widehat{{\gamma }^{\prime \prime }}=0.79$, SE = 0.02) compared to smaller individuals (mean $\widehat{{\gamma }^{\prime \prime }}=0.53$, SE = 0.06) (Fig. 4). Similarly, the probability of remaining in subterranean habitat ( $\widehat{{\gamma }^{\prime }}$ range: 0.59–0.96) was higher for larger individuals (mean $\widehat{{\gamma }^{\prime }}=0.86$, SE = 0.02) than smaller individuals (mean $\widehat{{\gamma }^{\prime }}=0.77$, SE = 0.17) (Fig. 4). Abundance ($\hat{N}$) at the surface during each primary period ranged from 62 (SE = 3) to 377 (SE = 9) (Fig. 5).

A typical recaptured individual was not observed for at least one survey period before being observed again, which is consistent with the high estimates of temporary emigration. The median time for between-period recaptures was 105 days, which is the average time interval between surveys. In some cases, individuals were not observed for several survey periods before being recaptured. For example, 20% of between-period recaptures (66 of 333) were not seen for more than one year at some point during the study, and 13 individuals had at least a two-year span between observations.

Discussion

Our results contradict characterizations of Barton Springs salamanders as primarily surface-dwellers, or as subterranean inhabitants accidentally flushed to the surface (see, e.g., Sweet, 1978; Chippindale, Price & Hillis, 1993). Barton Springs salamanders at Eliza Spring can be found in abundance at the surface (Bendik & Dries, 2018; Dries & Colucci, 2018), but often move underground. High estimates of γ″ demonstrate that temporary emigration of individuals from the surface to the subsurface is frequent. A large proportion of individuals also remained underground between surveys, as shown by similarly high estimates of γ′. Some individuals were not recaptured until more than a year later, suggesting that salamanders may spend long periods of time below ground before returning to the surface. Movements between the surface and subsurface may be a large component of previously observed surface population dynamics, which are characterized by high variability in abundance (Bendik & Dries, 2018; Dries & Colucci, 2018).

High temporary emigration in Barton Springs salamanders is consistent with observations of other plethodontid salamander populations (Bailey, Simons & Pollock, 2004a; Bailey, Simons & Pollock, 2004b; Price, Browne & Dorcas, 2012; O’Donnell & Semlitsch, 2015; Bendik, 2017; Drukker et al., 2018). In general, surface activity of plethodontid salamanders tends to be low, with surface abundance being a small fraction of the superpopulation size (e.g., see Table 1 in O’Donnell & Semlitsch, 2015). Plethodontids migrate beneath the forest floor (Bailey, Simons & Pollock, 2004a; O’Donnell, Thompson & Semlitsch, 2014) or to aquatic subterranean refugia (Bonett & Chippindale, 2006; Bendik & Gluesenkamp, 2013; Steffen et al., 2014) in response to moisture conditions at the surface, which can influence apparent demographic patterns (Connette, Crawford & Peterman, 2015; O’Donnell & Semlitsch, 2015). Although temporary emigration in plethodontids is common, reasons for lower surface availability of Barton Springs salamanders are likely different than those terrestrial salamanders. Lungless salamanders must keep their skin moist for cutaneous respiration and to prevent water loss (Feder, 1983) which influences their surface activity. Neotenic Eurycea require suitable levels of DO (Woods et al., 2010) and flowing water to maintain activity at the surface (Bendik & Gluesenkamp, 2013), but otherwise should be unaffected by the challenges of variable moisture faced by their terrestrial counterparts.

Models with environmental effects performed poorly compared to models with individual effects (body size) or time-varying parameters. Despite this, environmental effects can still provide some insight into the temporal variation we observed in salamander movement. Turbidity was positively associated with movement below ground, suggesting that either conditions worsen on the surface during higher turbidity or they improve below ground relative to the surface (e.g., via increased organic matter availability in the subterranean ecosystem). This pattern is consistent with observations of lower surface abundance when fine sediment deposition is high (Bendik & Dries, 2018). DO values were negatively associated with salamanders moving below ground but not significantly associated with movement to the surface. Salamanders may move in response to declining DO in anticipation of adverse conditions, for example, as values below 4.4 mg/L start to have negative effects on salamander growth and survival (Woods et al., 2010). However, this is speculative because DO remained above 4.8 mg/L during our study (Table 1). Furthermore, DO does not vary substantially between the near-surface and surface habitats and is higher at the springs than deep underground (Mahler & Bourgeais, 2013), so salamander movement below ground in search of higher DO conditions is unlikely. However, subterranean habitat could be less metabolically taxing for individuals (e.g., less movement and lower stress due to lack of predation) and prompt migration underground as DO declines. DO is correlated with flow, which was positively associated with survival and consistent with our hypothesis. Other studies have indicated the importance of this variable for Barton Springs salamander population demographics. For example, lagged flow is positively correlated with reproduction and low flows are associated with low abundance (Bendik & Dries, 2018; Dries & Colucci, 2018). Spring flow may also be associated with food availability within the aquifer, although this hypothesis is untested. Temperatures during our study were within the range for optimal growth for Barton Springs salamanders (15–24 °C; Crow et al., 2016), so response to temperature may be an indicator for other ecological conditions that improve survival, such as food availability, rather than indicating a direct effect on mortality. DO, temperature and flow are intertwined in this system and can exhibit correlated, but complex non-linear relationships (Mahler & Bourgeais, 2013); ultimately it is hard to isolate these effects. While an exhaustive post-hoc analysis of environmental covariates was beyond the scope of this study, this approach might be useful in the future to explore a suite of more complex models.

Smaller salamanders (<24 mm body length) had lower survival and a lower probability of moving underground, which was consistent with our prediction. Juvenile Eurycea have difficulty swimming against strong currents (Barrett et al., 2010) and are weaker swimmers than larger individuals (Azizi & Landberg, 2002). The spring outlets at Eliza Spring can reach high velocities (mean = 0.16 m/s, SD = 0.16, max = 0.65; City of Austin, 2020) in a range that has been shown to flush juvenile Eurycea cirrigera from substrates without gravel (Barrett et al., 2010). If smaller salamanders have more difficulty seeking subterranean refuge at the spring surface, they will be less able to respond to increases in competition for space or predation pressure compared to larger individuals, which more readily move into subterranean refugia. This may, in part, explain why survival is also lower for juvenile salamanders, and could result in a pattern of negative density dependence, as observed in a time-series analysis from two-decades of monthly Barton Springs salamander counts (Bendik & Dries, 2018). If this were the case, we should expect survival to vary with both size and year, given the differences in abundance we observed over the course of our study. In our analysis, however, time-dependent models of survival carried less weight (sum of w = 0.44) than those with size alone (sum of w = 0.66). Furthermore, a post-hoc comparison of S with an interaction between year and size did not improve upon our top model, suggesting that if there is a survival-density relationship on smaller salamanders, it is not apparent from these data. Negative density-dependence at the surface may also manifest itself through patterns of movement. There is a positive correlation between estimates of $\hat{N}$ at the beginning of the time interval and $\widehat{{\gamma }^{\prime \prime }}$ for that interval, and this relationship is stronger for the smaller size class (Pearson’s product-moment correlation, ρ = 0.60, P = 0.02) compared to the larger size class ( ρ = 0.49, P = 0.08). Therefore, salamanders may be moving underground in response to increased density at the surface, and there is evidence of density-dependent movement between habitats in other amphibians (Grayson, Bailey & Wilbur, 2011). Alternative statistical approaches to explore the relationship between density and demographic rates may be a useful avenue of future research (e.g., Kissel, Tenan & Muths, 2020).

Barton Springs salamanders may move into subterranean refugia for reproduction, or to avoid predation. Egg-laying occurs below ground and courtship may as well, and reproduction is not strictly seasonal as is the case for many other amphibians. Most individuals in the small size class were reproductively immature (85% were <22.5 mm in body length; Chippindale, Price & Hillis, 1993), which may partly explain why temporary emigration was lower for this group. Both predators and prey of Barton Springs salamanders are more abundant in surface habitat, but we did not measure variation in either of these factors for this study. Potential predators include crayfish (Procambarus clarkii), mosquitofish (Gambusia affinis) and larval damselflies (e.g., Argia spp.), while their prey includes a variety of invertebrates (DeSantis, Davis & Gabor, 2013; Gillespie, 2013; Owen et al., 2016; Davis, De Santis & Gabor, 2017). Below ground, the picture is quite different, where Barton Springs salamanders are likely the top predator (along with the sympatric E. waterlooensis), but invertebrate density is expected be much lower owing to a dearth of primary production in groundwater ecosystems (Gibert, Danielopol & Stanford, 1994). Barton Springs salamanders may therefore face a tradeoff between increased food availability at the surface vs. little to no interspecific predation in subterranean habitat. Growth is size-dependent and occurs rapidly during the first year of life (Fig. 3), so salamanders may take advantage of the resource-rich environment at the surface and therefore spend more time above ground when young. However, more research is required to understand ecological consequences of surface and subterranean habitat for these salamanders. Unfortunately, the subterranean ecology of the Barton Springs segment of the Edwards Aquifer is largely undocumented; even basic research on the diversity and distribution of its fauna is scarce (but see Hutchins, 2018; Nissen et al., 2018), despite otherwise extensive study of its hydrogeology and geochemistry (Hunt, Smith & Hauwert, 2019 and references therein).

Conclusions

Since Barton Springs salamanders were described from a single spring locality in 1993, knowledge of their population ecology has been gleaned from surface count surveys, limiting any inference about what we now know is a fraction of a larger superpopulation (City of Austin, 1997; US Fish and Wildlife Service, 2016; Bendik & Dries, 2018; Dries & Colucci, 2018). Using an indirect approach to detect movement, we quantified temporary emigration of Barton Springs salamanders into subterranean habitat from capture-recapture data. Temporary emigration was frequent, indicating that apparent abundance at the surface is not solely a function of immigration, recruitment and/or survival. Furthermore, these findings suggest that periodic declines in surface abundance of Barton Springs salamanders at Barton Springs are not wholly caused by mortalities, but may include emigration to adjacent, unobservable habitat. Our results demonstrate the utility of using estimates of temporary emigration to understand movement and habitat use in species where direct observation is not always possible.

Supplemental Information