Global realized niche divergence in the African clawed frog Xenopus laevis

Abstract Although of crucial importance for invasion biology and impact assessments of climate change, it remains widely unknown how species cope with and adapt to environmental conditions beyond their currently realized climatic niches (i.e., those climatic conditions existing populations are exposed to). The African clawed frog Xenopus laevis, native to southern Africa, has established numerous invasive populations on multiple continents making it a pertinent model organism to study environmental niche dynamics. In this study, we assess whether the realized niches of the invasive populations in Europe, South, and North America represent subsets of the species’ realized niche in its native distributional range or if niche shifts are traceable. If shifts are traceable, we ask whether the realized niches of invasive populations still contain signatures of the niche of source populations what could indicate local adaptations. Univariate comparisons among bioclimatic conditions at native and invaded ranges revealed the invasive populations to be nested within the variable range of the native population. However, at the same time, invasive populations are well differentiated in multidimensional niche space as quantified via n‐dimensional hypervolumes. The most deviant invasive population are those from Europe. Our results suggest varying degrees of realized niche shifts, which are mainly driven by temperature related variables. The crosswise projection of the hypervolumes that were trained in invaded ranges revealed the south‐western Cape region as likely area of origin for all invasive populations, which is largely congruent with DNA sequence data and suggests a gradual exploration of novel climate space in invasive populations.

Edwards, 2010; Petitpierre et al., 2012;Strubbe, Broennimann, Chiron, & Matthysen, 2013), these shifts are assumed to correspond to the release from dispersal barriers or biotic constraints rather than representing novel physiological adaptations (fundamental niche shifts; Tingley et al., 2014). A recent study suggests that environmental niches might be conserved within phylogenetic lineages (Schulte et al., 2012). However, few studies exploring the genetic foundation of the rate of adaptation and the role of intraspecific spatial variation in environmental niche structure exist. Furthermore, those which are available cover only few taxonomic groups (e.g., Krehenwinkel, Rödder, & Tautz, 2015;Lavergne & Molofsky, 2007;Schulte et al., 2012;Sillero et al., 2014;Vandepitte et al., 2014). Shifts in fundamental niches are only traceable in experimental setups and are difficult to disentangle using correlative models (Peterson, Papes, & Soberón, 2015). Therefore, most studies exclusively focus on shifts or stasis in realized niches.
Anthropogenic introduction pathways facilitate the spread of invasive species across biogeographical barriers and well beyond the species natural dispersal capacities (Wilson et al. 2009). This increased potential for large distance dispersal makes novel parts of the global environmental space colonizable, which may not be available in the native range. As unanticipated experiments, biological invasions represent model systems to study niche conservatism and dynamics (Jimenez-Valverde et al., 2011;Tingley et al., 2014). From a management perspective, the identification of variables restricting an invasive species' range is pivotal for successful management of established populations and to prevent further spread (e.g., Jimenez-Valverde et al., 2011). The comparison of realized niche spaces of native and invasive populations facilitates the identification of range limiting constraints of environmental variables that often correspond to physiological limits (Rödder, Solé, & Böhme, 2008;. Species distribution modeling (SDM) approaches have become promising tools for both basic and applied research (Kozak, Graham, & Wiens, 2008;Peterson et al., 2015) and have been frequently applied to predict the establishment success prior to issuing import permits for species of biosecurity concern (Bomford, Barry, & Lawrence, 2010;Kopecký, Patoka, & Kalous, 2016;Kumschick & Richardson, 2013). By comparing the environmental conditions of a species' native range with the conditions at the area of planned introduction, climate matching quantifies the risk of establishment, for example, for amphibian species in the European Union (Kopecký et al., 2016), plant taxa in Australia (Kumschick & Richardson, 2013), or freshwater fishes worldwide (Bomford et al., 2010). As such climate matching approaches assume climate niches to be conserved during the invasion process, they exclude the possibility that the native range does not correspond to the entire set of conditions a species can live in with its standing genetic background or that evolutionary processes expand the fundamental niche. The niche conservatism assumption has been frequently challenged in the past (Broennimann et al., 2007;Early & Sax, 2014;Rödder & Lötters, 2010;Urban, Phillips, Skelly, & Shine, 2008) partly because it neglects nonclimatic range limitations that also shape species spatial distributions. In addition, these biotic and abiotic factors might differ between the native and the invaded range (Colwell & Rangel, 2009;Early & Sax, 2014). While such limitations are undeniably extremely important for some species (e.g., regarding species interactions), these effects have been demonstrated to be highly scale dependent, for example, in the case of facultative predator-prey systems and rarely affect distributions at large extent and resolution if they depend on abundances (Soberón & Nakamura, 2009). Predator and prey may coexist if the abundance of the predator is below the carrying capacity of the local habitat. On the other hand, biotic interactions may be limiting irrespective of scale obligatory systems, for example, in specific butterflies where larvae require a specific host plant.
Furthermore, widespread, ecologically dominant species may be less affected by fine scale bionomic effects (Early & Sax, 2014).
The African clawed frog Xenopus laevis (Daudin, 1802) is a predominantly aquatic amphibian species native to the Mediterranean climate zone in southern Africa (Measey et al., 2012). The species was widely used for human pregnancy testing and distributed panglobally as a laboratory species for scientific research (Measey et al., 2012;van Sittert & Measey, 2016) where escapees and voluntarily released individuals established numerous invasive populations outside its native range (Measey et al., 2012). As X. laevis is abundant and widely distributed in its original range, has a high genetic variability, rapid growth, early sexual maturity, a broad diet, performs overland dispersal, is tolerant to various environmental conditions, and readily accepts heavily modified anthropogenic habitats, the species possesses great invasion potential (Measey et al., 2012. Among invasive amphibians, X. laevis has one of the highest recorded impacts on native fauna . These comprise competition for resources, direct predation of native amphibians (Amaral & Rebelo, 2012;McCoid & Fritts, 1980;Measey et al., 2015), a negative impact to the reproduction of native amphibians (Lillo et al., 2011) and has been shown to be an asymptomatic carrier of the chytridiomycosis fungus Batrachochytrium dendrobatidis (Fisher & Garner, 2007).
Being a model organism in a broad range of scientific disciplines, a variety of information on environmental constraints on X. laevis is available, including data on climatic factors that determine diurnal and annual activity patterns, reproduction and thermal tolerances (Balinsky, 1969;Casterlin & Reynolds, 1980;Eggert & Fouquet, 2006;McCoid & Fritts, 1989;Miller, 1982;Nelson, Mild, & Lovtrup, 1982;Wilson, James, & Johnston, 2000;Measey, 2016). Hence, comparisons of experimentally quantified properties of its fundamental niche with its realized niches are possible by comparing environmental conditions as observed at species occurrences with physiological and behavioral information obtained from the literature. This information may provide important insights into niche dynamics and the predictive ability-or shortcomings-of correlative SDMs.
Most specimens of X. laevis were exported through a single organization and represent a genetic lineage originating from the southwestern Cape region of South Africa (De Busschere et al., 2016;Lillo, Dufresnes, Faraone, Lo Valvo, & Stock, 2013;Lobos, Mendez, Cattan, & Jaksic, 2014;van Sittert & Measey, 2016). However, other lineages were detected in France suggesting that this invasive population is composed of animals from multiple source populations (De Busschere et al., 2016). Although X. laevis predominantly colonizes areas with equivalent climatic conditions, populations have also been established in temperate regions characterized by an oceanic climate (Fouquet & Measey, 2006;Measey & Tinsley, 1998;Rubel & Kottek, 2010). These conditions likely represent the edge of the species' physiological capacity, and possibly extend beyond environmental conditions within the native range.
Being invasive on multiple continents, including Europe, South, and North America, X. laevis represents a pertinent model organism to study realized niche dynamics during invasion processes. In this study, we assess whether the realized niches of invasive populations represent subsets of the realized niche in the species' native distributional range or if niche shifts are detectable. If so, human-mediated dispersal may have facilitated the exploration of fundamental niche space beyond the species' realized niche in its native range. Alternatively, niche shifts could indicate novel physiological adaptations. We expect potential niche shifts to be restricted to a subset of environmental variables as selective pressures may operate differently depending on the climatic setup of the area of introduction. This expectation appears to be reasonable given that environmental variables may interact very differently with the species' physiology: There may be hard limits such as critical thermal minimum or maximum being lethal or variables representing rather soft limits, for example, affecting the species behavior such as triggering reproduction. From a geographic point of view, the combination of environmental conditions may strongly vary across space creating potentially very different selective landscapes. Therefore, we ask whether the realized niches of the invasive populations still harbor information of the niche of the ancestor, that is, if it is possible to reconstruct the origin of the source populations in the native range.

| Data acquisition and preparation for environmental niche analyses
As environmental layers, we obtained remotely sensed bioclimatic variables (Beaumont, Hughes, & Poulsen, 2005;Busby, 1991) with a spatial resolution of 3 arc min (Deblauwe et al., 2016), representing minima, maxima, and average values of monthly, quarterly, and annual ambient surface temperature as well as precipitation (Table 1) Ihlow et al. (2016). We follow the taxonomic interpretation by Furman et al. (2015) confining the native range of X. laevis sensu stricto to southern Africa (South Africa, Lesotho, Swaziland, Namibia, Zimbabwe, and parts of Botswana, Mozambique and Malawi).

| Quantification of potential niche shifts
In order to quantify and assess potential niche shifts in the African clawed frog, we performed both univariate analyses using density profiles and multivariate hypervolume analyses (Blonder, 2016;Blonder, Lamanna, Violle, & Enquist, 2014) for the species' native distributional range in southern Africa and all known invasive populations in Europe, South, and North America. For univariate comparisons, we computed density profiles using the relevant functions of the sm package for Cran R (Bowman & Azzalini, 2014). The multivariate hypervolume analysis is designed to capture the environmental niche of the target species (or its populations) following Hutchinson's original idea of a Grinnellian niche space as an n-dimensional hypervolume (Blonder et al., 2014;Hutchinson, 1957). By computing the geometry and topology of the native and all invasive populations, hypervolumes can be quantified and compared in terms of shape, total volumes, niche positions, intersections, and unique parts (Blonder et al., 2014;Guisan et al., 2014). As environmental background an area defined by a circular buffer of 200 km surrounding all records from the species' native range was used in order to capture the available climate space and hence the potential for pre-adaptation.
The distribution of species occurrence records in orthogonal niche space is generalized based on multivariate kernel density estimations across principal components to remove effects of spatial autocorrelation and different availabilities of specific environmental combinations in geographic space (for details, see Blonder et al., 2014;Blonder, 2015Blonder, , 2016. Based on the kernel density estimations a new set of random records with a homogenous density across the environmental space captured by the original occurrence records was created. The total hypervolume is derived from this set of random records. This procedure allows the determination of the total niche volume of a species regardless of the geometry of the niche space (or within its respective native or invasive range) and is comparatively insensitive to low sample size. This niche volume can be projected back into geographic space. The resulting maps indicate areas exhibiting environmental conditions which are part of the species' niche volume (realized niche/realized distribution). The overlap between hypervolumes was determined using the Soerensen index [i.e., for hypervolumes A and B: S = 2*|A int B|/(|A| + |B|)]. We used the hypervolume package (Blonder, 2015) to delineate hypervolumes with two different approaches; (1) using a bandwidth of 0.5 enclosing all occurrences in environmental space (termed bdw herein) for delimitation of the multidimensional kernel) or (2) using a multivariate minimum convex polytope (mcp). The bandwidth approach assumes that all environmental conditions within the volume enclosed by a multidimensional buffer of 0.5 PC units around the species records is suitable.
In order to identify those environmental conditions of the native range which are also realized in invasive ranges, and vice versa, we projected the hypervolumes computed for the native and invasive ranges crosswise into geographic space. Assuming some degree of niche conservatism these areas are the most likely areas of origin. At the same time, our approach allows to track any shifts in realized climate niches and to quantify their direction and magnitude.
These test statistics were computed for the native and introduced populations using a set of 1,000 random points situated in a circular buffer of 200 km enclosing the occurrence records using the dismo package  for Cran R.

| RESULTS
The PCA revealed four PCs with eigenvalues >1 (Table 1) accounting for 90.2% of the total variation. The first PC explained 38.4% of the total variance and was mainly correlated with the maximum temperature of the warmest month and quarter, annual temperature range,   annual precipitation, as well as precipitation of the wettest month and quarter. PC 2 (24.1%) was mainly correlated with precipitation of the driest month and quarter, precipitation seasonality and isothermality, whereas PC 3 (17.5%) was mainly correlated with the minimum temperature of the coldest month and quarter. The fourth PC summarizes 10.3% of the total variance and was mainly correlated with the mean temperature of the wettest quarter and precipitation seasonality. For further details, see Table 1.
Univariate comparisons among bioclimatic conditions at native and invaded ranges revealed that the invasive populations are well nested in the variable range of the native population in the mean annual temperature range, temperature seasonality, the maximum temperature of the warmest month, the minimum temperature of the coldest month, the temperature annual range, the mean temperature of the warmest quarter, the annual precipitation, and the precipitation of the wettest month and quarter. However, bioclimatic conditions in at least one invaded region exceed those in the native range in the annual mean temperature, isothermality, mean temperature of the wettest, driest, warmest and coldest quarter, precipitation of the driest month and quarter, precipitation seasonality, and precipitation of the warmest and coldest quarter (Figure 1). The most deviant invasive population were those from Europe. Bimodal peaks in some density plots suggested that European populations can be split into two distinct bioclimatic groups (Figure 1). While the left peak in temperature related variables (BIO 1,5,6,10,11,and 15, Figure 1) corresponded to the invasive populations from France and Great Britain which are both characterized by oceanic climate, the right peak referred to invasive populations from Portugal and Sicily, areas characterized by a Mediterranean hot climate (Rubel & Kottek, 2010). In precipitation related variables (BIO 14,17,and 18) Test statistics revealed the predictive performance of the fourdimensional hypervolume models to differ among regions, wherein F I G U R E 2 Density plots for the four principal components with Eigenvalues >1. Background conditions were extracted within a 200-km buffer enclosing the native species records. See Table 1 for details  (Table 2).
The total volume of the realized niche based on the complete occurrence data set was 395.34 for the bdw and 1403.70 for the mcp approach. Assessing populations from the native and the three invaded regions separately, models revealed that the native population possesses the largest realized niche volume, followed by the invasive populations in South America. The two approaches revealed different results for Europe and North America (Table 3). Intersections between hypervolumes were largest between the native and the invasive populations from North America, and between populations from North and South America. In the other comparisons intersections were low, yielding low niche overlaps in terms of Soerensen indices (Table 3).
The largest centroid distance was observed between populations from North America and Europe and the smallest between populations from North and South America (Table 4). PC 1 had the highest influence on total niche volumes in all cases wherein the contributions of the remaining PCs varied between 0.49 and 0.92 (Table 5).

| DISCUSSION
Our results suggest some degree of realized niche shift between the native and all invaded regions, which is restricted to 10 of 19 bioclimatic variables and well visible in the density plots showing PC 2 and PC 4 ( Figure 2). This differentiation is also reflected in the results of the niche overlap analyses. Crosswise projection of the hypervolume models trained in invaded ranges onto the native range revealed the south-western Cape region as the likely area of origin for all invasive populations.

| Niche shifts
Physiological performance of ectothermic taxa, such as amphibians strongly depends on environmental conditions (Feder & Burggren, 1992). Affecting physiological functions such as digestive and locomotive performance, temperature is among the most important extrinsic factors (Angilletta, Niewiarowski, & Navas, 2002;Bennet, 1990;Huey & Kingsolver, 1989). In X. laevis and other pipid frogs, temperature was found to affect traits related to fitness like sprint performance (Miller, 1982), swimming velocity (Herrel & Bonneaud, 2012;Wilson et al., 2000), endurance capacity (Herrel & Bonneaud, 2012), reproduction, and larval development rates (Balinsky, 1969). Previous research suggested X. laevis to prefer an ambient temperature of 22 °C (Miller, 1982). Thermal tolerances were found to range from 14 °C to 32 °C (Casterlin & Reynolds, 1980). However, other authors suggested that X. laevis withstands more extreme thermal conditions (Nelson et al., 1982), and critical thermal limits of 2 °C and 39 °C, respectively, were obtained for the invasive population in France (J. Courant and A. Herrel, pers. comm.). Laboratory studies indicate that life stages might not be equally affected by extreme thermal conditions (Balinsky, 1969;Wu & Gerhart, 1991) and ontogenetic shifts during embryonic development seem to alter thermal tolerances (Nelson et al., 1982). However, the variability of the critical thermal limits among the native populations of X. laevis is not known and available information is limited to adult specimens originating from populations from the Western Cape.
A comparison of the experimentally established preferred temperature range and critical limits with our results reveals the annual mean temperatures in colonized areas to be lower than the species' thermal preferences except for Chile. However, ambient temperatures generally remain within the critical temperature range of the species.
During the warmest quarter temperatures correspond to the species' preferences in the native range and North America, while remaining lower in Europe and exceeding this range in South America (Figure 1). were also reported to actively avoid exposure to critical conditions, for example, by moving into cooler parts of the water body when exposed to high temperatures (47 °C) in Arizona (G. J. Measey, unpublished data) or excavate pits into the soft bottom mud where water temperatures remain lower (McCoid & Fritts, 1980). Albeit freezing may be a major mortality factor (Eggert & Fouquet, 2006), frogs were recorded moving at 5 °C (Wilson et al., 2000) and throughout winter in France (J. Courant and A. Herrel, unpublished data). While adults may still be active and foraging, even under ice-bound water, development of tadpoles under these conditions is hampered (G. J. Measey, unpublished data). However, the alleged extinction of two invasive populations that persisted for decades in Great Britain was linked to extreme weather events (cold and drought) of successive winters suggesting that climatic conditions exceeded the buffering capacity of the microhabitat inhabited by X. laevis (Tinsley, Stott, Viney, Mable, & Tinsley, 2015).

Maximum monthly temperatures
The observed deviation in minimum annual temperature between the native and invaded ranges suggests niche shifts of the realized niches approximating the critical minimum temperatures for X. laevis.
Induction of reproductive activity is correlated with both temperature and precipitation. Balinsky (1969) suggested that spawning in the native range in Gauteng at approximately 1,000 m a.s.l. can be triggered by a period of rainfall, whereas Green (2002)  bodies and hence availability of water is frequently not limiting reproductive activity. Temperature represents another potential constraint for reproduction. With water temperatures regularly reaching 20 °C, climatic conditions in California are less extreme than in the species' native range (McCoid & Fritts, 1980) leading to an almost year-round reproductive activity with continuous ovulation (McCoid & Fritts, 1989). It was suggested that low ambient temperatures might hamper reproduction (Nelson et al., 1982). However, the species was found to reproduce during most of the year in France with a strong peak during the summer months (April to June; Courant et al., in press). Hence, low winter temperatures may reduce recruitment but do not necessarily prevent the establishment of invasive populations. (c)

| Origin of invasive populations
1. Nonadaptive process exposition on historic climate cycles may allow a pre-adaptation beyond the current realized niche of a species leading to the possibility of acclimatization or phenotypic plasticity, as was suggested for, for example, the invasive freshwater crayfish Procambarus clarkii (Chucholl, 2011).
From experimental studies, Wilson et al. (2000) observed significant acclimation effects on swimming velocity, suggesting a pre-adaptation capacity to cooler environments, which may partly explain the invasion success of X. laevis in Europe.
However, SDM projections of a broad range of different South African amphibian species on paleoclimate scenarios revealed only marginal shifts in potential distributions under last glacial maximum conditions (Mokhatla, Rödder, & Measey, 2015;Schreiner, Rödder, & Measey, 2013) suggesting an overall limited probability of a historic pre-adaption to cooler environments.

2.
A first mechanism involving evolutionary changes would be an establishment of the initial invasive populations in climates similar to the environmental conditions within the area of origin and subsequent adaptation and niche expansion. This scenario is supported by the finding that the invasive populations in, for example, Chile originated from a single introduction and show a signature of a recent bottleneck and are genetically very similar indicating a rapid spread from a small source population (Lobos et al., 2014).
Although, multiple introductions occurred in the other invaded areas, in most cases all individuals originated from a single source region. Based on simulations, adaptation to novel environments may be facilitated in initially small population sizes (Holt, Barfield, & Gomulkiewicz, 2005). Selection for a certain beneficial allele might lead to faster fixation in a small population than in a large population. However, before their escape or release, frogs have been kept for several generations under artificial conditions in laboratories or captive breeding facilities, which might also influence the predisposition of the individuals.

3.
A second mechanism for rapid adaptation increasing thermal tolerances was suggested by Krehenwinkel et al. (2015). These authors

| CONCLUSION
In the present study, we demonstrate that invasive populations of X. laevis are established well beyond the species' multivariate realized niche in southern Africa. This finding has important implications for both macroecological niche theory and practical aspects of risk assessments using climate matching.
By suggesting that climate niches remain stable during the invasion process, the observed dynamics in the realized niches of X. laevis challenge previous findings. Although we cannot disentangle shifts in the fundamental niche from shifts in the realized niche alone, comparisons with natural history data including environmental tolerances and triggers for reproduction suggest that most niche shifts observed can be explained by realized niche shifts. One exception where this explanation is not evident is the invasive population in France. Here, further studies are required to test whether hybridization of different lineages has enabled a shift in the species' fundamental niche.
Given the magnitude of the detected niche shifts, the usefulness of climate matching approaches to assess invasion risk is challenged, as it might frequently underestimate the true potential distribution of a species when a geographic subset of the species' realized distribution is used for model training. As suggested by Broennimann and Guisan (2008) inclusion of all available distribution information may improve predictions but underestimations are still possible.
Furthermore, as previous authors have noted a careful analysis of the available environmental space within the training area of the models and a quantification of areas requiring extrapolation beyond this training range may improve the reliability of assessments by exclusion of those areas with high uncertainty (Elith, Kearney, & Phillips, 2010;Measey et al., 2012). Unfortunately, this procedure will restrict the application of climate matching to areas characterized by environmental conditions most similar to the training range and for the majority of other areas predictions which require extrapolation remains unreliable.
From a management point of view, it is possible to assess the reliability of correlative SDMs by comparing the species' realized environmental niche with what is known from its fundamental niche in terms of physiological performance and physiological limits. As demonstrated herein for X. laevis, it can be expected that its true invasion potential is larger than its estimated potential distribution based on its currently realized niche. In this case, the climate matching approach is well able to identify areas with high risk of further invasion but is less reliable in identifying actually unsuitable habitats.
Further conservation measures need to prevent additional introductions in suitable habitats which may arise due to anthropogenic climate change (Ihlow et al., 2016), not only within the already invaded regions but also in those areas which are highlighted by our global hypervolume model.

ACKNOWLEDGMENTS
This study was funded by BiodivERsA (project title: "Invasive biology of Xenopus laevis in Europe: ecology, impact and predictive models") and the Deutsche Forschungsgemeinschaft (DFG RO 4520/3-1). The publication of this article was funded by the Open Access fund of the Leibniz Association.