Ecological niche models display nonlinear relationships with abundance and demographic performance across the latitudinal distribution of Astragalus utahensis (Fabaceae)

Abstract The potential for ecological niche models (ENMs) to accurately predict species' abundance and demographic performance throughout their geographic distributions remains a topic of substantial debate in ecology and biogeography. Few studies simultaneously examine the relationship between ENM predictions of environmental suitability and both a species' abundance and its demographic performance, particularly across its entire geographic distribution. Yet, studies of this type are essential for understanding the extent to which ENMs are a viable tool for identifying areas that may promote high abundance or performance of a species or how species might respond to future climate conditions. In this study, we used an ensemble ecological niche model to predict climatic suitability for the perennial forb Astragalus utahensis across its geographic distribution. We then examined relationships between projected climatic suitability and field‐based measures of abundance, demographic performance, and forecasted stochastic population growth (λs). Predicted climatic suitability showed a J‐shaped relationship with A. utahensis abundance, where low‐abundance populations were associated with low‐to‐intermediate suitability scores and abundance increased sharply in areas of high predicted climatic suitability. A similar relationship existed between climatic suitability and λs from the center to the northern edge of the latitudinal distribution. Patterns such as these, where density or demographic performance only increases appreciably beyond some threshold of climatic suitability, support the contention that ENM‐predicted climatic suitability does not necessarily represent a reliable predictor of abundance or performance across large geographic regions.


| INTRODUC TI ON
A core aim of ecology and biogeography is to understand how environmental conditions influence the abundance or demographic performance of species across space. Improved insight regarding conditions that support higher abundance or demographic performance of species can in turn increase understanding of the factors constraining range limits and inform conservation efforts aimed at identifying sites capable of sustaining species of conservation concern (Araújo & Williams, 2000). Long-standing biogeographical theory predicts that in the absence of density dependence, abundance and demographic performance of a species should peak at the center of its geographic distribution and decrease monotonically with proximity to the range boundary. This is theorized to be due to linear degradation of environmental suitability from the center to the edge of the range (the Abundant Center Hypothesis-hereafter ACH; Hengeveld & Haeck, 1982;Brown, 1984;Sexton, McIntyre, Angert, & Rice, 2009). However, empirical support for the ACH is weak (Abeli, Gentili, Mondoni, Orsenigo, & Rossi, 2014;Sagarin & Gaines, 2002;Sagarin, Gaines, & Gaylord, 2006;Sexton et al., 2009), suggesting that variation in environmental suitability for the occurrence and performance of a species is often more nuanced than can be explained by geographic position alone (Holt, 1997;Weber, Stevens, Diniz-Filho, & Grelle, 2017). Rather than the generally inaccurate practice of using range position as a proxy for abundance and demographic performance, metrics of habitat suitability incorporating more realistic representations of geographical variation in environmental conditions could be used to more accurately predict density and demographic performance across the landscape (e.g. Tôrres et al., 2012;Van Couwenberghe, Collet, Pierrat, Verheyen, & Gégout, 2013;Thuiller et al., 2014;Matthiopoulos, Field, & MacLeod, 2019). This concept lies at the heart of the niche centroid hypothesis (Martinez-Meyer, Diaz-Porras, Peterson, & Yanez-Arenas, 2013), which predicts that centrality in the sense of proximity to ideal niche conditions rather than geographic centrality is predictive of the abundance of a species across space. The niche centroid hypothesis posits that cases in which the abundance of a species is correlated with its proximity to the distribution's geographic center are a coincidental outcome of geographic centrality aligning with niche centrality.
One way in which niche centrality is estimated is through the use of ecological niche models (ENMs), which describe the environmental niche of a species in order to predict where suitable environmental conditions for its presence occur across the landscape. ENMs make predictions based on correlations between spatially explicit records of a species' presence (or presence and absence) and the environmental conditions at the location of the record (Elith et al., 2006;Peterson & Soberón, 2012). In other words, these models determine the environmental conditions that correspond to the presence of a species and use this information to predict suitability for its occurrence across the landscape. While newer approaches to distribution modeling may directly or indirectly incorporate predictors describing the presence of interacting species, dispersal, or indicators of nonclimatic environmental conditions such as disturbance (Bucklin et al., 2015;Engler et al., 2009;Linder et al., 2012;Randin, Vuissoz, Liston, Vittoz, & Guisan, 2009;Wisz et al., 2013), the majority of ENMs solely use large-scale abiotic variables such as climate and topography as predictors of habitat suitability (Austin, 2007;Bucklin et al., 2015;Pearson & Dawson, 2003).
Species' distributions have often been described as a spatial manifestation of their niche (sensu Hutchinson, 1957), with the assumption that the processes which define and govern presence within niche space also determine the abundance and performance of the focal species within it (Holt, 2003). Patterns of occurrence, abundance, and demographic performance are often theorized to be linked, such that higher occupancy of a species in an area indicates higher habitat quality and thus higher abundance and demographic performance of the species (Holt, 1997;Gaston et al., 2000;Holt, Gaston, & He, 2002; but see Dallas & Hastings, 2018). The practice of utilizing ENM-derived predictions of environmental suitability as an indicator of a species' abundance or performance assumes a relationship between projected suitability for presence and its local abundance (suitability-abundance relationship) or demographic performance (suitability-demography relationship) driven by the same environmental conditions associated with its presence. This assumption hinges upon the expectation that the environmental conditions which determine a species' presence similarly influence its abundance and demographic performance.
Although it has been claimed that the relationship between landscape-scale occupancy and local abundance is among the most strongly supported in ecology (Holt et al., 2002), the ubiquity of suitability-abundance relationships remains a topic of substantial debate. On one hand, Weber et al. (2017) found support for a general pattern of positive abundance-suitability relationships in meta-analysis that examined 450 species spanning 30 studies. On the other hand, a recent analysis of the relationship between ENM-predicted suitability and abundance of 396 mammal and tree species found no such relationship (Dallas & Hastings, 2018). The smaller number of studies that have evaluated the relationship between demographic performance and projected suitability yield similarly equivocal results: positive suitability-demography relationships have been observed in some cases (Brambilla & Ficetola, 2012;McLane & Aitken, 2012;Monnet, Hardouin, Robert, Hingrat, & Jiguet, 2015;Searcy & Shaffer, 2016;Sheppard, Burns, & Stanley, 2014), but appear to be the exception rather than the rule (Bacon et al., 2017;Bayly & Angert, 2019;Chardon, Pironon, Peterson, & Doak, 2020;Csergő et al., 2017;Oliver et al., 2012). The lack of consensus regarding the ubiquity of suitability-abundance and suitability-demography relationships may be due to local-scale constraints not captured in the large-scale predictors often used in ENMs (Davis, Jenkinson, Lawton, Shorrocks, & Wood, 1998;Guisan & Thuiller, 2005;Lembrechts, Nijs, & Lenoir, 2019;Pearson & Dawson, 2003;Varner & Dearing, 2014;Wisz et al., 2013), such that suitability predicts a maximum abundance or performance at a site that is then altered by attributes of the local environment (VanDerWal, Shoo, Johnson, & Williams, 2009).
Studies that examine the relationship between projected environmental suitability and a species' abundance and its demographic performance across its geographic distribution are rare. Those which do so tend to find contrasting patterns in suitability-abundance and suitability-demography relationships, and in that suitability is often positively correlated with abundance but not with population performance or estimated persistence (Oliver et al., 2012;Bean et al., 2014;Thuiller et al., 2014; but see Pironon, Villellas, Morris, Doak, & García, 2015). More studies of this kind are needed to further our understanding of whether ENM predictions of suitability might prove broadly useful for predicting species' abundance and performance across their distributions.
Here, we explore suitability-abundance and suitability-demography relationships across the latitudinal distribution of the perennial forb, Astragalus utahensis. Specifically, we ask whether predictions of climatic suitability derived from correlative ENMs accurately predict variation in the abundance and demographic performance of this species across its latitudinal distribution and what environmental predictors contribute most strongly to estimates of climatic suitability across the distribution.

| Study system
Astragalus utahensis occurs in spatially isolated populations on sparsely vegetated hillsides throughout its latitudinal distribution, which extends from southern Utah, USA (roughly 38°N) to southeastern Idaho, USA (roughly 43°N; Figures 1 and 2). It grows as a spreading basal rosette of villous compound leaves; the primary growing season is April through early July. Reproduction is solely by seed: plants produce racemes of 3-10 large magenta flowers which develop into hirsute seed pods containing medium-large gravity-dispersed seeds that fall in late June-early July (Green, 1976). Although the distribution of A. utahensis is relatively small, abiotic conditions vary substantially across the latitudinal range. Climatic data (from 1970 to 2000) show higher mean annual temperatures in the center of the latitudinal range compared with its northern and southern edges and a gradient of high to low annual precipitation from the northern edge of the range to its southern edge (Worldclim v. 2.0; Fick & Hijmans, 2017; see Appendix S1).
To characterize the current distribution of A. utahensis, we ob- The "southern edge" study populations were not the furthest south within our query of herbarium records, as access to the remote areas represented by records at the extreme southern edge was deemed too difficult for logistical feasibility. Nonetheless, the four southern edge study populations included in our study were within 120 km of the furthest south recorded populations for this species (Figure 2; Table 1; populations 1-4).

| Ecological niche model construction
We removed duplicate records from herbarium records and field observations of A. utahensis, as well as records that were misidentified and those without at least two decimal degrees precision for latitude and longitude coordinates. This resulted in a final dataset of 106 presence records. Our sampling frame for the construction of ENMs was a polygonal convex hull that buffered presence points by 100 km. We used a sampling frame larger than the distribution of presence records because a previous study found that habitat suitable for the presence of A. utahensis exists more than 200 km beyond the northern edge of its current distribution (Baer & Maron, 2019).
An average nearest neighbor analysis revealed significant spatial autocorrelation between presence records within the sampling frame; we corrected for this by randomly thinning presence records to a minimum distance of 13 km between presence records (spatstat package; Baddeley & Turner, 2005). This thinning distance maximized the number of records retained while ensuring that average nearest neighbor distance between presence records was not significantly different from that of an equally sized sample randomly drawn from within the sampling frame (Monte Carlo simulation one-sided p = .16). Of the 106 original presence records, 83 were retained post-thinning. We built models using a dataset consisting of the 83 thinned presence records and 830 pseudoabsences (10 times the number of presences; 10n) selected randomly from fishnet cells within the study area that had not contained a presence record removed during the thinning process. Although Lobo and Tognelli (2011) espouse the use of 100n pseudoabsences, using 10n pseudoabsences yielded models with higher mean accuracy and is in keeping with methods advocated by Chefaoui and Lobo (2008).
Data describing climatic conditions from 1970 to 2000 associated with each presence and pseudoabsence record were queried from the Worldclim 2.0 dataset at a 30 arc-second (~1 km) resolution and used as explanatory variables in ENMs (Fick & Hijmans, 2017).
We selected a subset of climatic predictors for inclusion in the ENMs to eliminate multicollinearity among variables. To do this, we first performed univariate generalized linear models for A. utahensis presence using each climatic predictor, retained the variable with the greatest explanatory power, and discarded all variables correlated with it at |r| > 0.7. We repeated this process until no additional variables remained, leaving us with the following list of climatic predictors, all correlated at |r| ≤ 0.7: precipitation of the driest month (bio14), temperature seasonality and precipitation seasonality (bio4 and bio15, respectively), isothermality (bio3), and mean temperature of the coldest, wettest, and driest quarters (bio 11, bio 8, bio 9, respectively).
ENMs for A. utahensis were built using the biomod2 package (version 3.3-7.1; Thuiller, 2019) in R (version 3.5.3; R Core Development Team, 2019). To avoid bias associated with the use of a single model algorithm and leverage the strengths of multiple algorithms (Araújo & New, 2007), we utilized an ensemble modeling approach. Ensemble models comprised a total of six algorithms: three regression algorithms and three machine learning methods.
We performed 10 runs of each algorithm; each model run was trained on a random selection of 90% of presences and pseudoabsences and internally cross-validated with the remaining 10% of presence and pseudoabsence data. We evaluated the accuracy of F I G U R E 2 Map showing locations of thinned A. utahensis presence records included in the construction of ENMs (from herbarium collection records and field surveys conducted for this study, gray points), study populations (squaressouthern edge populations, circles-central populations, and triangles-northern edge populations), and predicted suitability from the ensemble ENM. Gray indicates areas of predicted absence and green indicates areas of predicted presence, with darker gray representing lower suitability and darker green higher suitability each model run using the area under the curve of the receiver operating characteristic (AUC; Hanley & McNeil, 1982) and the true skill statistic (TSS; Allouche, Tsoar, & Kadmon, 2006). The value of AUC ranges from 0 to 1, with AUC = 1 representing a perfect prediction and AUC ≥ 0.7 a useful prediction (Swets, 1988). TSS is similar to Cohen's Kappa (Cohen, 1960) but is less sensitive to prevalence (Allouche et al., 2006). TSS ranges from −1 to 1, where TSS < 0 indicates a prediction no better than random, TSS = 1 indicates perfect prediction and thresholds for evaluating model performance mirror those for Cohen's kappa, with a TSS value of ≥0.4 representing a useful model (Allouche et al., 2006).

| Density and demographic monitoring
We

| Statistical analyses
All statistical analyses were performed using R Statistical Software (version 3.5.3; R Core Development Team, 2019). We used linear regressions to test for linear and exponential relationships between ENM-predicted climatic suitability (hereafter, suitability) and (a)

| Ecological niche model
The

| Suitability, density, and demography
The log-transformed summed density of all life stages in 2016 and mean summed reproductive and juvenile plant density showed a significant J-shaped relationship with suitability ( Figure 3; Table 2).
Log-transformed seedling density in 2016 exhibited a marginally significant J-shaped relationship with suitability (Table 2)

| D ISCUSS I ON
Across the entire latitudinal range, climatic suitability was not a strong predictor of A. utahensis density or demographic performance. and performance are also only influenced by climatic suitability beyond some threshold.

| Suitability versus abundance
The lack of a linear suitability-abundance relationship across A.
However, suitability-abundance relationships are generally weaker in plants than animals, and evidence suggesting that the strength of these relationships is unaffected by the proportion of the range examined is derived solely from studies of mammals and birds (Weber et al., 2017). Studies that have examined suitability-abundance relationships in plants often build ENMs using occurrence data collected according to a political boundary rather than throughout a species' entire distribution, which may yield biased conclusions regarding the ubiquity of observed relationships across the distribution (e.g. Van Couwenberghe et al., 2013;Thuiller et al., 2014).

TA B L E 2
Outcomes of repeated-measures linear regressions and GLMMs of A. utahensis density and demographic rates versus ENMpredicted suitability, study year, and the interaction of suitability and year (where applicable)

F I G U R E 4
Relationship between ensemble ENM-predicted suitability in central and northern study sites only and mean (±SEM) (a) annual per capita seed production, (b) cumulative recruitment from 250 seeds over 3 years of a seed addition study, and (c) annual survival,  (Dallas & Hastings, 2018). Where linear suitability-abundance relationships are weak or nonexistent, this may be attributable to: (a) differences in the identity or scale of climatic conditions that dictate presence compared to abundance (Boulangeat, Gravel, & Thuiller, 2012), (b) a stronger impact of disturbance and land-use history (unrelated to overarching climatic conditions) on abundance (Jiménez-Valverde et al., 2009;Nielsen, Johnson, Heard, & Boyce, 2005), (c) the influence of local biotic interactions on demographic performance and resultant density (Cabral & Kreft, 2012;Davis et al., 1998;Guisan & Thuiller, 2005;Linder et al., 2012), or (d) variation in dispersal ability that allows persistence in low-suitability areas but prevents the colonization of highly suitable areas (Guisan & Thuiller, 2005;Pulliam, 2000;Van Horne, 1983).
Few of these potential explanations for the lack of a suitabilityabundance relationship apply to A. utahensis. For example, while A. utahensis is dispersal limited beyond its northern range boundary (Baer & Maron, 2019), this does not explain the pattern of low density in populations of intermediate suitability that yielded a J-shaped rather than linear suitability-abundance relationship. Moreover, A. utahensis generally grows in sparsely vegetated areas where negative intraspecific density dependence is rare or absent, although interspecific competition may impact density by excluding the species from more benign microsites in some populations (K.C. Baer, personal observation). Predispersal seed predation and pollen limitation decrease the demographic performance (and presumably the density) of A. utahensis, but the extent to which this occurs is similar across populations experiencing varying climatic suitability (Baer & Maron, 2018).
One potential explanation for the observed pattern is that the climatic conditions associated with the suitability of an area for the presence of A. utahensis may have relatively little impact upon its local abundance compared with other aspects of the climate in some areas of the distribution due to local adaptation, which could yield a suitability-abundance relationship such as the one we observed (e.g. Peterson, Doak, & Morris, 2019). For example, the climatic conditions that most strongly constrain A. utahensis' abundance could differ substantially in different parts of its latitudinal range.
It is possible that the conditions contributing most strongly to ENM predictions of suitability for the presence of A. utahensis also play a large role in determining its abundance in central populations (two of the three highest-abundance populations were central populations), but are of relatively minor importance relative to other climatic conditions in other portions of the range. A second possibility is that the observed suitability-density relationship could simply be an artifact of sampling bias: occurrence records could be more numerous where populations are larger and denser, resulting in higher suitability scores for these more visible high-abundance populations.
If bias in the number of presence records only appeared beyond a particular threshold of density, suitability may only exhibit a positive correlation with density once that threshold is reached. A third explanation for the consistently low density of populations below a particular suitability threshold may be that abundance simply does not respond strongly to climatic variation until some threshold of climatic favorability is met, regardless of the conditions that define a favorable climate in a particular region.
The expectation of a linear suitability-abundance relationship assumes that abundance of a species is positively related to its proximity to the "niche center" (as approximated by ENM suitability), much as the abundant center hypothesis predicts abundance to be tied to geographic centrality. Although we did not observe the linear relationship between a proxy for niche centrality (suitability) and abundance, the J-shaped relationship between density and suitability was echoed by a strong correlation between density and the latitudinal position of study populations, with higher density in central than peripheral study sites characteristic of an abundant center distribution (K.C. Baer and J.L. Maron unpublished data). This suggests that estimates of suitability derived from our ensemble ENM (aka niche centrality) are no better a predictor of A. utahensis density across the distribution than the latitudinal position of a study population. This may be because our ensemble model does not adequately describe the entirety of A. utahensis' climatic niche. The ensemble model's specificity value indicates an absence was predicted where a presence had been recorded in roughly twelve percent of cases, meaning that it may often underestimate suitability for this species. Furthermore, evidence of dispersal limitation in A. utahensis (Baer & Maron, 2019) suggests that the current distribution may not fully reflect the range of environmental conditions suitable for its presence and that estimates of suitability based upon climatic conditions within the current distribution may be somewhat biased.

| Suitability versus demographic performance
While some of the demographic parameters we examined were not related to predicted suitability, several important demographic rates did have significant relationships with suitability. The J-shaped relationship between λ s and the suitability of central and northern edge study populations appears attributable to a combination of similar relationships in cumulative recruitment and mean annual survival and a positive linear relationship between recruitment rates from the seed bank and suitability. We previously found that significant declines in λ s between central and northern edge populations were primarily attributable to three demographic rates, of which germination rates from the seed bank and total seed production at the population scale were two (Baer & Maron, 2019). While per capita seed production was not significantly correlated with the suitability of central and northern edge study sites, reproductive plant density was positively associated with suitability in those sites, which likely led to increased total seed production at the population scale in more suitable populations (K.C. Baer and J.L. Maron, unpublished data). It seems likely that the relationship between λ s and suitability is driven primarily by differences in reproductive plant density and germination from the seed bank with smaller contributions from annual survival and overall recruitment rates.
The J-shaped nature of the relationship between suitability and λ s in these populations in central and northern populations and strong positive correlation that we observed between overall A. utahensis density and λ s in these populations, coupled with the J-shaped suitability-density relationship across all southern-northern study populations suggests that a similar J-shaped relationship exists between suitability and λ s across the entire latitudinal distribution. Such a relationship between suitability and demographic performance would indicate that ENM predictions of suitability cannot be directly translated into predictions for the performance of a species across its geographic distribution, as the threshold beyond which demographic performance responds strongly to suitability is not known a priori.
Previous studies that have examined both suitability-abundance and suitability-demography relationships typically find that the former are significant while the latter are not (Bean et al., 2014;Oliver et al., 2012;Pironon et al., 2015;Thuiller et al., 2014). Similarly, studies of solely suitability-demography relationships commonly find that suitability is not correlated with overall demographic performance in any consistent manner (Bacon et al., 2017;Bayly & Angert, 2019;Csergő et al., 2017;Pironon et al., 2018). As with explanations for a lack of a suitability-abundance relationship, the lack of a suitability-demography relationship is generally attributed to the influence of local conditions that yield habitat suitability which differs from climatic suitability. These may include heterogeneity in a site's local physiography, the availability of necessary resources unrelated to climate, or the distribution of facilitative partners, competitors, or predators, all of which may be obscured by single-species SDMs built using solely large-scale climatic predictors (e.g. Bean et al., 2014;Bayly & Angert, 2019). In addition, suitability based upon a model built for the entire distribution may not reflect the climatic optimum for all populations and may poorly predict demographic responses to climate variation across the distribution (Chardon et al., 2020;Peterson et al., 2019).
As with the suitability-density relationship, the nonlinear suitability-λ s relationship we observed in central and northern populations and similar relationship predicted to exist across the entire latitudinal distribution may be attributed to differences in the environmental attributes that most strongly constrain demographic performance among geographic regions. For example, variation in microsite conditions in populations at the southern edge of the distribution may exert stronger control over density and demography in those populations than do larger-scale climatic conditions. This possibility is supported by the fact that three of the four southern edge populations were located in an ephemeral stream bed (South Twin Peak) or beside an irrigation canal (Meadow Creek 1 & 2), which may have led these populations to perform better than predicted by large-scale climatic conditions if the sites in which they occurred were less water-limited than the surrounding landscape. According to the model's suitability threshold for presence of the species (0.203), the suitability of these areas was predicted to be too low for the presence of the species, so supplemental water availability may account for their existence and persistence.
The relatively poor demographic performance of populations with intermediate levels of predicted suitability (of which most were northern edge populations) may indicate that some aspect of the local environment unrelated to climate suppresses demographic performance to a level lower than might be anticipated if the suitability-demography relationship was linear. While A. utahensis tends to grow in sites that support few other species, it is possible that superior competitors in northern edge populations pushed the plants into even less suitable microsites than they inhabit at either low or high suitability.
The shape of the suitability-λ s and suitability-density relationships that we observed was due in large part to high demographic performance and density in two high suitability populations near the center of the distribution (Uinta and Bountiful). As with density, it may be that demographic performance is uniformly poor below a certain threshold of climatic suitability, but responds strongly to increasing suitability beyond this threshold. If Uinta and Bountiful populations existed in areas where climatic conditions exceeded this threshold, this may have allowed them to flourish and help explain the resulting patterns we observed between suitability and both density and overall demographic performance in A. utahensis.

| CON CLUS ION
In this study, we found that projected climatic suitability displayed a J-shaped relationship with the abundance of a species across its geographic range and a similar relationship with overall demographic performance as measured by stochastic population growth rate across the center to north portion of the range. A significant correlation between local density and population growth rate in central and northern edge study populations suggested that the relationship between suitability and population growth rate across the entire latitudinal distribution would be similar to that between suitability and density, meaning that demographic performance may only respond to increasing suitability beyond a threshold that is not known in advance. This implies that suitability is not a reliable metric for estimating either the abundance or demographic performance of a species across its geographic distribution.

ACK N OWLED G M ENTS
The authors thank Peggy Stolworthy and the U.S. Forest Service for permission to perform work on private and public lands in Utah and Idaho. We also thank Zachary Baer, Christina Cain, Becky Fletcher, Allison Klocke, Mike Kimball, Eric Mohrmann, and Jennifer Neville for their assistance with field work and Dr. Thomas Edwards for training in ENM construction.

CO N FLI C T O F I NTE R E S T
None declared.

O PE N R E S E A RCH BA D G E S
This article has earned an Open Data Badge for making publicly available the digitally-shareable data necessary to reproduce the reported results. The data is available at https://osf.io/5ga32/.

DATA AVA I L A B I L I T Y S TAT E M E N T
Data and R scripts associated with this study are publicly available through the Open Science Framework digital repository at https:// osf.io/5ga32/.