Seed Availability Does Not Ensure Regeneration in Northern Ecosystems of the Endangered Limber Pine

Research Highlights: When biotic interactions such as disease alter both the seed production capacity of stands, and seedling survivorship, the relative importance of seed availability versus substrate specificity may alter future regeneration opportunities for plant populations. Background and Objectives: We investigated the importance of disease severity, seed availability, and substrate limitation to the regeneration dynamics of the endangered limber pine, Pinus flexilis, in two ecosystems with varying forest composition, and different histories of white pine blister rust infection (WPBR; Cronartium ribicola). Materials and Methods: A total of 17 stands from the montane ecoregion (Alberta, Canada) were sampled for seed production between 2007–2010, seedling density, and age structure. Model selection using an information theoretic approach compared a series of a priori models and their combinations, based on our hypotheses on the role biotic variables play in the regeneration process. Results: Despite higher rates of WPBR infection, 2.3 times more seed was available for avian dispersers in the southern ecosystem. Recent seedling regeneration did not correspond to seed production (83 versus 251 seedlings/ha, in southern versus northern ecosystems, respectively), resulting in a seven-fold difference in seed to seedling ratios between ecosystems. Models suggest that disease and vegetation cover were important factors explaining the absence of regeneration in 79.4% of the plots sampled, while basal area (BA) of live limber pine, rocky substrates, ecosystem, South aspects, and slope enhanced limber pine regeneration. Seedling age structures suggest that recent regeneration is less likely in more diseased landscapes, than it was historically (40% versus 72.8% of seedlings < 20 years old, respectively, in southern versus northern ecosystems). Conclusions: At the northern limits of limber pine’s range, seed availability does not ensure regeneration, suggesting that other environmental or biotic factors hinder regeneration. Regeneration was consistently predicted to be lower in the southern ecosystem than in the northern ecosystem, suggesting that natural regeneration and the potential for population recovery are ecosystem dependent. We recommend that monitoring recent seedlings will aid the identification of biotic and abiotic factors affecting regeneration.


Introduction
One of the fundamental challenges in recovering endangered species is identifying the key factors that allow populations to be self-sustaining despite increased threats. The conservation of plants frequently draws on a rich history of seed ecology studies that test whether species are seed-limited or substrate limited [1]. The tradeoffs between seed size and number of offspring are self-evident in nature with variation ranging 2-3 orders of magnitude in seed size amongst co-occurring species with a similar life form and size at maturity [2]. Seed limitation in coniferous ecosystems may occur Reduced seed availability due to WPBR-induced adult mortality may compromise the seed dispersal mutualism with the Clark's nutcracker [28], leading to seed limitation in landscapes with severe WPBR infestation. While several studies report natural regeneration dynamics without addressing the effects of WPBR [12,26], these studies can be viewed as documenting the potential for "in situ" regeneration, which is of critical importance to population maintenance in disease altered landscapes. Regeneration traits of limber pine make it particularly susceptible to seed limitation, notably its early succession status and occupancy of arid habitats [29], large seed, and dependency on vertebrates for dispersal [30].
We investigated the importance of seed versus substrate limitation to regeneration dynamics of the endangered limber pine in two ecosystems with different histories of WPBR infection. One ecosystem, occurring in Kootenay Plains, AB, represents the northern-most populations of limber pine in North America [16,18], and thus is ecologically significant for the study of demographic characteristics. Our objective was to determine whether seedling densities and age structure differ between ecosystems with different WPBR infection levels, and more broadly whether seed limitation and substrate availability remain useful constructs for maintaining populations of endangered trees in disease altered ecosystems. We hypothesize that severe WPBR landscapes have less regeneration than low WPBR landscapes at the northern periphery due to seed limitation, and greater WPBR induced mortality of seedlings over time. We further hypothesize that seed and substrate limitation processes diminish in importance for maintaining regeneration when other ecological processes (i.e., disease) contribute to species decline. These hypotheses have important implications for recovering an early succession species like limber pine, since fire is presumed to be critical for creating new microsites for regeneration, and substrate controls to regeneration in existing stands have not been investigated. Disease may compromise the importance of seed availability (directly) and substrate availability and competition (indirectly) to the regeneration process by reducing seed production and increasing seedling mortality on favorable microsites.

Study Area
Our study was conducted at the northern limits of limber pine range in the Rocky Mountains of North America. Suitable habitat occurs in the montane region, where the climate is primarily Cordilleran, and has a mean daily temperature of 13 • C from June to August, the warmest months, and −10 • C from December to February, the coldest months. A July precipitation maxima occurs and 580 mm precipitation occur annually [31]. We selected two distinct geographical regions of the Montane Ecoregion in Alberta, Canada, that were approximately 400 km apart. Our northern study area (latitude 52.00 • N, longitude 116.50 • W), henceforth referred to as the "northern ecosystem", is situated within the North Saskatchewan River Valley, extending from the Kootenay Plains Ecological Reserve to provincial forests adjacent to Abraham Lake. Our southern study area (latitude 49.60 • N, longitude 114.20 • W), henceforth referred to as the "southern ecosystem", extended east of the Rockies to include the Porcupine Hills formation, and foothills ridges in the Crowsnest Pass ( Figure 1). Climate records from 1999 to 2008 indicate that the northern ecosystem was cooler and drier over this period (mean low, −4.5, mean high 7.2, 469 mm precipitation, 13.09 MJ/m 2 estimated incoming radiation) than the southern ecosystem (mean low, −1.4, mean high 10.4, 525 mm precipitation, 13.65 MJ/m 2 estimated incoming radiation) [32].  (rectangles) in Alberta relative to limber pine distribution in North America (reprinted with permission from [33]).
Limber pine in the North Saskatchewan River Valley (northern ecosystem), is found on a variety of aspects (116°-270°) on rocky ridges, scree slopes, planes, and erosional cliff slopes with alluvial deposits, at elevations of 1348-1455 m. There, limber pine occurred in pure stands, or in association with white spruce, Picea glauca (Moench) Voss, or lodgepole pine, Pinus contorta Douglas Loudon. In the southern ecosystem, limber pine occurred on southwest to west-facing (214°-270°) on the tops of rocky ridges, at similar elevations of 1300-1504 m. There, stands occurred either as limber pine or limber pine/Douglas fir (Pseudotsuga menziesii (Mirb) Franco) mixtures.
The historic fire regime for Alberta's limber pine is uncertain, although as a montane species, it may typically experience a mixed severity fire regime [34]. This regime may have differed between Limber pine in the North Saskatchewan River Valley (northern ecosystem), is found on a variety of aspects (116 • -270 • ) on rocky ridges, scree slopes, planes, and erosional cliff slopes with alluvial deposits, at elevations of 1348-1455 m. There, limber pine occurred in pure stands, or in association with white spruce, Picea glauca (Moench) Voss, or lodgepole pine, Pinus contorta Douglas Loudon. In the southern ecosystem, limber pine occurred on southwest to west-facing (214 • -270 • ) on the tops of rocky ridges, at similar elevations of 1300-1504 m. There, stands occurred either as limber pine or limber pine/Douglas fir (Pseudotsuga menziesii (Mirb) Franco) mixtures.
The historic fire regime for Alberta's limber pine is uncertain, although as a montane species, it may typically experience a mixed severity fire regime [34]. This regime may have differed between the two study regions, with less frequent crown fires predominating in the northern study sites where adjacent subalpine fuel complexes are present, while more frequent, lower intensity fires from adjacent grasslands likely influenced southern study areas. Wildfires have been actively suppressed in both regions for most of the 20th century, and cattle grazing has replaced wildfire in the southern study area. White pine blister rust was detected in 1952 in our southern study area [35], where in 2009 it infected 38% of our study trees (n = 3 study stands; a subsample of stands from Smith [17]). Surveys in our northern study area documented both the occurrence and rapid increase in WPBR infection (1% live tree infestation in 2003, 20% in 2009; n = 3 study stands sampled from Smith [17]).

Sampling Methods
We selected stands based on the presence of prior WPBR infection monitoring data, and their representation of tree species combinations and site attributes associated with limber pine populations (Table 1). Detailed WPBR infection and mortality data was available for three stands from each study area (see, [17]). Study stands were situated from 300 m to 10 km apart from the nearest stand, and >5 ha in size. Marked changes in site types and associated changes in tree species composition frequently occurred in proximal stands. We sampled nine stands in the northern study area, which ranged from pure limber pine, to associations with white spruce, or lodgepole pine. In the southern ecosystem, we sampled eight stands, containing either limber pine, or limber pine and Douglas fir associations.
Stands were sampled for regeneration during the summers of 2008-2009. From a randomly chosen point at the forest edge, we commenced sampling 50 m into the forest. We sampled four sites in each stand along a transect, at random distances of 80-120 m from adjacent sites, due to intervening gaps in limber pine distribution. Each site had a minimum of 10 proximal trees of cone-bearing size. Live and dead BA area by tree species, site-level environmental variables, and annual cone availability for regeneration (2008-2010) from 10 mature trees at each site (40 trees per stand) were recorded at these sites, as previously described [5]. At each of the four sample sites within each stand, three 25 m 2 regeneration plots were situated 15 m from the transect (chosen randomly from the four cardinal directions; n = 12 plots per stand, totaling 300 m 2 ). Thus, 108 and 96 plots were sampled in northern and southern study areas, respectively. Plots were searched systematically for all regeneration ranging from current year germinants to seedlings < 50 cm in size. The height and age of all seedlings were recorded. Seedling age was estimated non-destructively by counting all visible terminal bud scars from the leader tip down to ground level. On all seedlings, terminal bud scars were readily visible for recent years, but could not be reliably counted beyond 20 years for many seedlings; consequently, 20 served as the cut-off value to distinguish between "recent seedlings" and older seedlings.
In order to quantify substrate and seedling regeneration relationships amongst stands and study areas, a 4 m 2 substrate plot was nested within the 25 m 2 plot. The larger plots were oriented to maintain homogeneity of substrate conditions throughout. Substrate abundance was recorded according to the following % cover classes: 1 = 0-5%, 2 = 5.1%-10%; 3 = 10.1%-25%; 4 = 25.1%-50%; 5 = 50.1%-75%, and 6 = 75.1%-100%. During data analysis, we substituted the categorical values with the median value from each cover class range. Substrates were classified as rock, scree (rock < 10 cm diameter), mineral soil, humus, needles, leaf litter, and moss. Additionally, we recorded the actual rooting substrate of each seedling according to the listed substrate categories. Recent litter and humus at the base of each tree were carefully brushed aside. Due to the frequent longevity of seedlings on substrates, many seedlings were assigned two or more rooting substrates (i.e., mineral and needle). Seedlings were classified as healthy (uninfected with WPBR), infected with inactive branch canker, or active branch canker. To evaluate the effects of vegetation cover on seedling abundance, we independently recorded grass cover, herbaceous cover, ground shrub cover (0-0.15 m) and low shrub cover (0.16-1 m), and then summed values to construct an ordinal leaf area index based on these four strata for vegetation up to 1m in height, within a 4 m 2 plot. The percent cover of known WPBR hosts (to genera) was recorded separately as well (i.e., Ribes and Castelleja). Table 1. Stand structure data and site attributes for study sites used to evaluate recent limber pine regeneration relative to seed availability and substrate conditions. Live trees were tallied separately from standing dead trees in basal area counts (BA).

Stand
Live

Model Identification
To test which ecological factors are associated with limber pine seedling regeneration, we used a generalized linear model framework to fit models to our seedling counts from 25m 2 sample plots. Additionally, to account for the hierarchical nature of our sampling design, we used a mixed model approach and nested the random effect of plots within sites and stands via a random intercept in all models [36]. We tested whether a zero-inflated modelling approach was required [37]. Model covariates were only retained when correlations of less than 0.5 were measured (9 of 15 covariates retained; Table 2). The descriptions of model covariates are provided in Table 2. Model selection was conducted using Akaike's Information Criterion corrected for small sample size (AICc) to determine which model was the most parsimonious from our a priori candidate model set [38,39]. While care should be taken using a coefficient of determination estimated for general linear mixed models, we assessed model fit using formulations by Nakagawa and Schielzeth [40] and Johnson [41]. Following the identification of the best candidate model, we inspected model residuals to confirm that errors were homogeneous and independent with respect to model covariates. We analyzed the number of regenerated limber pine seedlings using a two-step procedure. First, we compared four competing model forms using a standardized model with respect to independent variables. Using AICc, we compared Poisson, zero-inflated Poisson, negative binomial, and zero-inflated negative binomial model forms. Secondly, using the model form identified in the previous step, we compared a series of a priori models and their combinations, based on our hypotheses on the role variables play in limber pine regeneration. Our candidate set of models were constructed around the two part nature of the zero inflated model, which contains covariates for both the zero inflated and count component of the model. We included in every model a "null" within the count component related to physical aspects of the sites measured, including ecosystem, slope, and aspect, since there is ample evidence from plant studies that these covariates affect regeneration [1]. Additional to the count component, we added variables related to the hypotheses related to seed and substrate limitation, while the zero component of the model was constructed using variables related to the hypotheses of disease and competition. Our a priori candidate set therefore contained all combinations of the four hypotheses of interest within their respective components of the model. We identified three additional models that include seed limitation, disease, and competition as part of both the count and zero inflated component of the model. All analyses were conducted using the glmmTMB package in R, model selection used bbmle, and ggeffects was used to make model predictions [42][43][44][45][46][47][48]. In all cases we fit the conditional model based on the global model (all covariates listed in Table 2) and included a nested random effects (plots measured at sites in stands) via a random intercept. The zero-inflated portion of the model included only the covariates of basal area of dead limber pine and vegetation cover. Models are shown in decreasing rank, with values for model log-likelihood (LL), number of estimated parameters (K), Akaike's Information Criterion for small samples (AICc), AICc difference (∆AICc), and model weight (W i ).

Age Structure Analysis
We used a one-sample, Kolmogorov-Smirnov goodness of fit test (KS-test) for discrete data [49] to test whether seedling ages in each study region reflected a uniform establishment pattern over the past twenty years (n = 68, and 20, in northern and southern ecosystems, respectively). A uniform age distribution, parameterized from 1 to 20, for the full twenty-year period of bud scar ages was used. Additionally, we used a two-sample KS-test for discrete data [49] to determine whether seedling regeneration in the past 20+ years differed between northern and southern ecosystems. Seedlings greater than 20 years were grouped into a single age category for each site. To minimize the likelihood of older seedlings unduly influencing the age distribution tests between study areas, we performed another KS test where all seedlings were grouped into 5 year age categories, and a 20+ year-old category.

Modeling Results of Regeneration Density
Recent seedling regeneration (<50 cm tall) was 3.02 times greater in the northern ecosystem than in the southern ecosystem (251 vs. 83 seedling/ha; Figure 2). Recent seedlings (<20 years old) were found in all nine stands in the northern ecosystem (67-533 seedlings/ha), while no seedlings were found in two of the eight stands surveyed in the southern ecosystem (0-233 seedlings/ha).
We recorded 0 seedlings at 79.4% of the plots we sampled and initial model selection suggested that a zero-inflated Poisson was more suitable for the data (w i = 0.737) relative to models without zero-inflation or using the negative binomial distribution (Table 3). Using the zero-inflated Poisson model, we found evidence that a simplified variant of the global model ("final model") was the most parsimonious model (w i = 0. 564, Table 4) and contained at least one covariate from each of the a priori hypotheses (Table 5). Model selection using corrected Akaike's information criterion (AICc) suggests that there is strong support for the final model, which included seed, substrate, disease, and competition, over any less inclusive model containing three or fewer covariates (∆AICc ≥ 3.4; Table 4). In the final model, the zero-inflated model component was comprised of covariates related to the disease and competition hypotheses, while the conditional (count) model component was comprised of covariates related to the hypotheses of seed availability and substrate. Aspect, slope, and ecosystem from the null model also contributed to the count component. There was no evidence that the inclusion of seed limitation, disease, and competition covariates into both the count and zero inflated component of the model was warranted (all w i < 0.05; Table 4) and the covariates were not significant within the models to which they were included (all p < 0.5). The final model adequately described the data, with the marginal model (fixed effects) explaining 68% of the variation in the data, while the conditional model (fixed and random effects) explained 73% of the variation in the data.

Seed Availability
Site BA of live limber pine explained regeneration density better than recent seed production data (2008)(2009)(2010) and was the only seed availability variable retained in the final model (p < 0.001, Table 5). Analyses show that seedlings increased at higher live limber pine BA on average (Table 5). Model predictions suggest that seedlings increase across a range of limber pine BA in both ecosystems, but seedling abundance is greater in the northern ecosystem across a range of seed availability. Interestingly, other seed availability covariates, namely recent cone production, and cone escape, were not retained in the final model.

Substrate
While not abundant (<5% ground surface), the percent of rock was retained in the final model as significant explanatory variable for plots with seedlings (p < 0.001) (Figure 2d). Commonly available substrates, such as mineral soil and humus, which reversed availability between northern and southern ecosystems (57% versus 8%, and 9% versus 42%, mean cover respectively; Figure 3), were not retained in the global model.

Microsite
Seedlings regenerated best on South aspects (135 • -225 • ), particularly in comparison to West aspects (225 • -315 • , p = 0.038; Figure 2b). A variety of aspects occurred in northern ecosystem, while west aspects predominated in the southern ecosystem. Regeneration was associated with increasing slope, though not significantly (p = 0.133), which ranged from 5 • -38 • , and was retained in the final model (Table 5).

Ecosystem
While regeneration density was 3.02 times higher in the northern ecosystem and added significant information to the final model (p < 0.001; Table 5), ecosystems varied in many important attributes. The northern ecosystem had lodgepole pine and white spruce instead of Douglas fir, had a wider variety of aspects (Table 1), markedly more mineral soil (Figure 3), and less understory vegetation. There was a 7.0-fold difference in density between ecosystems when the seed production:seedling ratio is computed between ecosystems (the southern ecosystem produced 1.70 times as many cones per tree × 1.35 times as many trees per hectare, divided by 0.33 times as many seedlings). This difference could be accounted for by the incremental contributions of seed, microsite, substrate covariates which enhanced regeneration, and their greater prevalence in the northern ecosystem, and lower levels of conditions that reduced regeneration (disease and competition). Even after accounting for the variety of processes found to influence regeneration in our final model, regeneration was consistently predicted to be lower in the southern ecosystem than in the northern ecosystem. These findings suggest that other processes are limiting regeneration; amongst the test covariates, we found evidence that both disease and vegetation cover can limit regeneration and may do so disproportionately in the southern ecosystem.

Disease
Site BA of dead limber pine helped explain regeneration absence and added information to our final model. A wider range of dead limber pine BA occurred in the southern ecosystem (Figure 2e). Model selection identified the negative effect of disease, by including it in the zero-inflated component of the model (Table 4). Disease appears to reduce the chances of regeneration on sites where other conditions in the model would facilitate regeneration. Intermediate host species to WPBR, such as Ribes spp. and Castilleja spp., were observed in several stands, but occurred too infrequently in regeneration plots to be included in model selection.
in the northern ecosystem. These findings suggest that other processes are limiting regeneration; amongst the test covariates, we found evidence that both disease and vegetation cover can limit regeneration and may do so disproportionately in the southern ecosystem.

Disease
Site BA of dead limber pine helped explain regeneration absence and added information to our final model. A wider range of dead limber pine BA occurred in the southern ecosystem (Figure 2e). Model selection identified the negative effect of disease, by including it in the zero-inflated component of the model (Table 4). Disease appears to reduce the chances of regeneration on sites where other conditions in the model would facilitate regeneration. Intermediate host species to WPBR, such as Ribes spp. and Castilleja spp., were observed in several stands, but occurred too infrequently in regeneration plots to be included in model selection. Figure 3. Abundance of substrates by cover classes in mature limber pine stands in northern and southern ecosystems (mean ± 1 s.e., n = 9 and 8 study stands, respectively). Percent cover classes were: 1: ≤ 1%, 2 = 2%-5%, 3 = 6%-10%, 4 = 11%-25%, 5 = 26%-50%, 6 = 51%-75%, and 7 = 76%-100%.

Competition
Seedling absence was associated with increasing understory cover (p = 0.074; Table 5; Figure 2f). Model selection identified the negative effect of vegetation cover, by including it in the zero-inflated component of the model (Table 4). Grass was an important understory component in the southern ecosystem, while ground cover, comprised of bearberry (Arctostaphylos uva-ursi (L.), Spreng.), creeping juniper (Juniperus horizontalis, Moench), and grouseberry (Vaccinium scoparium Leiberg ex Coville), were common in the northern ecosystem. In the northern ecosystem, limber pine seedlings were absent in stands with abundant lodgepole pine.

Age Structure
Recent seedling regeneration patterns differed significantly between the northern and southern ecosystems (dmax = 15.8, p = <0.02; Figure 4). The southern ecosystem has proportionately fewer young seedlings (40% < 20 years old) available to recruit into older age classes than in the north (72.8% < 20 years old). In the most recent 20 year period, regeneration appears to be continuous in both the northern and southern ecosystems (Figure 4), and there was no evidence of age distributions of seedlings differing from a uniform distribution (dmax = 6, p > 0.2, and dmax = 0, p > 0.5, respectively, one-sided KS test).
(72.8% < 20 years old). In the most recent 20 year period, regeneration appears to be continuous in both the northern and southern ecosystems (Figure 4), and there was no evidence of age distributions of seedlings differing from a uniform distribution (dmax = 6, p > 0.2, and dmax = 0, p > 0.5, respectively, one-sided KS test).

Evaluation of Hypothesis Tests
It was timely to investigate the importance of regeneration processes occurring within existing stands for several reasons: (1) a rapid decline of limber pine occurred in the northern Rocky Mountain portion of its range in recent decades [17], and (2) there is a requirement for in situ maintenance of existing populations in Alberta's recovery plan [16]. We based our hypotheses on the role various biotic and abiotic covariates have on limber pine regeneration from a rich history of seed versus substrate research in many systems, and knowledge as a conservation community on both threats and challenges to the natural regeneration of limber pine. Seed availability, substrate conditions, disease, and competition feature prominently in conservation recovery plans, so a full factorial approach to testing all possible combinations of these covariates was warranted. Model selection using an information theoretic approach shows that seed, substrate, disease, and competition (final model) collectively add information to our "null" ecological model, highlighting that limber pine regeneration remains a complex process, best modelled by a variety of processes.

Evaluation of Hypothesis Tests
It was timely to investigate the importance of regeneration processes occurring within existing stands for several reasons: (1) a rapid decline of limber pine occurred in the northern Rocky Mountain portion of its range in recent decades [17], and (2) there is a requirement for in situ maintenance of existing populations in Alberta's recovery plan [16]. We based our hypotheses on the role various biotic and abiotic covariates have on limber pine regeneration from a rich history of seed versus substrate research in many systems, and knowledge as a conservation community on both threats and challenges to the natural regeneration of limber pine. Seed availability, substrate conditions, disease, and competition feature prominently in conservation recovery plans, so a full factorial approach to testing all possible combinations of these covariates was warranted. Model selection using an information theoretic approach shows that seed, substrate, disease, and competition (final model) collectively add information to our "null" ecological model, highlighting that limber pine regeneration remains a complex process, best modelled by a variety of processes. Excluding any one of these four covariates in the factorial design, resulted in a noticeably weaker model. This was an important finding because WPBR features prominently in the literature documenting the decline of high elevation pines, and thus could have overridden other regeneration processes [7,17]. We chose to include ecosystem, and microsite in our "null model", and in all more-inclusive models in the factorial design, since there is ample evidence from plant studies that these environmental covariates affect regeneration [1], and thus provide a reasonable ecological base for modelling the regeneration process. A useful contribution of our final model is the parsing out of the influence of the covariates on the zero and the count components. Interestingly, disease and competition (vegetation cover) helped explain the zero inflated component of the model, while seed availability, site, substrate, and ecosystem are important in the "count" component. A careful evaluation of the covariates in hypotheses is warranted to understand the regeneration processes they model, their relative contribution to the natural regeneration mechanism of limber pine, and threats facing the species.

Seed Availability
Seed availability is a key factor controlling regeneration density in numerous studies in forest ecosystems [29,30]. Our best performing covariate in our final model, BA limber pine, measured at the site level (i.e., within 15 m of regeneration plots), was a strong predictor of regeneration density. Basal area of seed trees in wind-dispersed species and serotinous species [6,8] is a strong indicator of seed production and seedling regeneration in boreal ecosystems. Dyszoochorous species, like limber pine, which are reliant on the Clarke's nutcracker for long distance seed caching and regeneration [50], would likely have weaker relationships between BA of seed trees and seedling regeneration than other dispersal syndromes, given the reliance on the birds caching habits in determining regeneration patterns. Given the lack of documentation of barochory (gravity dispersal) as a mechanism for successful regeneration of limber pine and the high rates of post-dispersal seed predation by rodents that likely preclude germination, greater seed availability is unlikely to directly contribute to greater seedling density. We suspect that nutcrackers, which are effective dispersers of limber pine [20], responded to increasing BA of seed trees (1.7-11.2 m 2 /ha amongst study stands), and made more caches in the same stands where seed is available, resulting in increases in regeneration density. Within this range of forest cover, mature limber pine stands clearly remain open enough to facilitate regeneration, a fact likely due to the slope and southern aspects that increase solar penetration.
One measure of actual seed available for nutcrackers, the number of cones escaping seed predators in a 2010 mast year (a proxy to historic seed availability patterns), did not provide additional information in the final model. This is not surprising, as recent cone availability may not adequately reflect the variety of species interactions that affect cone availability (i.e., predation, dispersal) that influenced the longer time span over which seedlings regenerated (>20 years). Seed production per tree varied more than 20-fold between sites, and across both ecosystems mast years were critical for limber pine cones to escape red squirrels (Tamiasciurus hudsonicus Erxleben, 1777) interannually [5]. Our research on squirrel predation rates in these landscapes shows that differences in cone availability are conserved even under a wide range of seed predator densities and habitat conditions left in diseased stands [5]. These findings suggest that historically more cones have been available for dispersers in the southern ecosystem, yet this has provided no comparable benefit for regeneration. One possible explanation is that there may be dispersal limitation, not seed limitation in southern sites.
High elevation pines may be disperser limited if cone production and stand health is poor [28], and cone predation by red squirrels limits seed availability in non-mast years [5]. Cone production thresholds of 1000 cone/ha for nutcracker dispersal in whitebark pine [28] suggest that BA of live trees, dead BA, and cone production would all contribute to whether limber pine is dispersal-limited. Such values of cone production were obtained in many of our study stands in the 2010 mast year (Peters, unpublished data), and likely occurred in many of the previous mast years at our study sites. Our variable dead BA of limber pine directly measures the loss of seed producing trees and may be interpreted as the reduced capacity of stands to attract avian dispersers, and a reduction in actual seed cached within seed sources. Research on the seed-dispersal behavior of nutcrackers in northern limber pine ecosystems is necessary to understand their influence on the spatial distribution of limber pine regeneration.

Substrate Availability
Seedling age structures show that germination microsites remain available for recruitment events in mature limber pine stands, particularly in the northern ecosystem. Webster and Johnson [12] used age structures in our study area to suggest that microsites remain available in mature stands as old as 600 years, a fact attributable to the arid environments and niche of limber pine in primary succession. Although seedlings established proportionately more often in plots with mineral soil than expected based on their relative abundance [51], particularly in the southern ecosystem, the availability of mineral soil did not add information to our final models. Coop and Schoettle [25] found that mineral soil and cobble were negatively associated with limber pine seedlings, while leaf litter was positively associated. They attributed this association to the "stability" of microsites, where litter indicates sites where plant growth could be sustained. Our modelling results identified rock, a rarely encountered substrate across both ecosystems (7% mean ground cover), as adding more information than commonly encountered substrates like mineral soil and humus, which differed markedly between ecosystems (Figure 3). Seedlings did not actually regenerate on rock, hence the variable is indicative of a terrain that may function more in nutcracker caching behavior [21], as a nurse object aiding seedling survival [24], or as a site that accumulates soil moisture near crevices, which is a determinant of limber pine regeneration [52]. While substrate helps explain seedling counts in our final model, a finer resolution assessment of substrates with respect to seedling germination and survivorship is needed to help clarify the role of substrate in the establishment phase.

Processes Inhibiting Regeneration
Disease-related factors appear to reduce the chances of successful limber pine regeneration at the northern limits of the species range. Several of our findings support factors identified in other studies that influenced limber pine regeneration. Firstly, regeneration was more abundant in areas with low blister rust infection and fewer WPBR-killed snags (252 and 83 individual seedlings/ha, <50 cm tall in our northern and southern ecosystems, respectively). Smith et al. [17] found nearly identical seedling cluster densities in 2009 in our study areas, with 250 and 100 seedling/ha (<1.3 m tall) in the Kootenay Plains (northern ecosystem) and Porcupine Hills/Whaleback (southern ecosystem). In range wide monitoring studies in the central and southern Rocky Mountains, seedlings <1.37 m in height averaged 141 stems/ha, but declined with WPBR infection levels (tree cankers, branch dieback, and dead BA) [53]. While earlier limber pine stand dynamics studies have suggested that age structures reflect natural stand dynamics in disease-free stands [12,[25][26][27], our study suggests WPBR is changing regeneration patterns. We suspect that the severity and duration of WPBR infection in the southern ecosystem is contributing to limited regeneration. Currently, WPBR is increasing more rapidly at the northern limits of both limber [17] and whitebark pine's range [54] than it is in central regions, which will likely alter the age structures of regeneration in the future.
Our cone production studies indicate that in our study system, disease lowers seedling regeneration by influencing the seedling mortality process, rather than through an overriding negative effect on seed availability [5]. We have to be cautious in identifying what disease-related effects dead BA models. This variable not only reflects reduced seed availability (see Section 4.2.1), but may also contribute to zero values by reducing survivorship of seedlings. Proportionately fewer recent seedlings in the southern ecosystem may suggest that younger seedlings are succumbing to WPBR infections. Given their larger size, older seedlings may persist for several years following infection. This interpretation is hard to reconcile with other studies that have found that taller seedling have higher rates of infection and mortality than shorter seedlings in limber [17,53] and whitebark pine [54], respectively. In fact, only 0 and 1% of seedlings < 0.5 m tall were found to be infected in 2003 and 2009 range-wide monitoring surveys in Alberta [17]. We suspect that some of these differences in findings may be attributable to our attention to smaller seedling classes and new germinants that may be under-represented in health monitoring protocols, and which would not persist following infection or mortality. Provincial monitoring of WPBR in Alberta show that infection levels and mortality in adult trees increased from 1% to 20% between 2009 to 2014 (within a subset of our northern study stands), and more broadly across the northern part of limber pine range, from 2% to 11% [17]. This increase could have contributed to increased seedling mortality even in our healthier northern ecosystem, underscoring the "resolution" challenges in inferring disease and regeneration relationships from retrospective age analyses of living seedlings. Surveys of whitebark pine regeneration at a regional level in Alberta show similar patterns to our study ecosystems, with consistent regeneration in northerly sites where WPBR infection is low, and very limited regeneration in the southern portion of the province where WPBR infection is high [54]. Collectively, these studies suggest that disease effects are relatively new in the northern ecosystem, hence many seedlings may die in coming decades, reducing their potential to sustain in situ regeneration.

Conservation Biology of Limber Pine
The natural regeneration processes that we identified in our analyses show that prior limber pine regeneration dynamics do not reflect current or future dynamics. Surviving seedlings in the southern ecosystem reflect historic regeneration events, culled by disease over several recent decades, while the northern ecosystem reflects continuous regeneration patterns, less affected by WPBR, a recent player on the landscape [17]. Mature limber pine stands that are disease-free are no longer present, even in the northern ecosystem; such ecosystems may still where they provide benchmarks for "expected regeneration" densities for the southern part of the species range. In mature stands, with no evidence of recent fire, range-wide surveys of seedlings vary from 172 clusters/ha in the northern Rocky Mountains (derived average from [17]) to 141 stems/ha in the central and southern Rocky Mountains [53]. Post-fire studies have reported similar densities, and while seedling regeneration is highly variable between fires, averages range from 38-314 [25], and 290-508 stems/ha [27]. As an early succession species in arid environments, "natural" rates and densities of seedling establishment are quite low, particularly compared to other conifer species that invade limber pine stands in more mesic sites [55]. Growing season moisture and soil temperature can limit first year survivorship of limber pine seedlings [52]. Climate differences between our two ecosystems may have contributed to differences in regeneration patterns; however, greater annual moisture in the southern ecosystem (525 versus 469 mm) may compensate for the added heat stress of its warmer climate (−1.4 vs. −4.5 • C, mean annual temperature). Climate differences also likely contributed to the differences in codominant tree species between ecosystems, and the ensuing difference in vegetation cover. It is important to note that both of our ecosystems support fewer seedlings than are considered necessary to maintain self-sustaining populations of limber pine in Alberta, given the threats that limber pine face provincially [16].
Age distributions of seedlings clearly show that the southern ecosystem has proportionately fewer recent seedlings available to recruit into older and larger size classes than the northern ecosystem. Depictions of age versus height suggest that seedlings take about 20 years to reach 25 cm in height (Peters, unpublished data). Difficulties with identifying terminal bud scars increase as seedlings age [56], limiting our ability to comment on recruitment dynamics after this point. Age structure differences suggest either higher rates of attrition occur in recent southern seedlings, or that barriers to regeneration in the past 20 years are disproportionately affecting the southern ecosystem. The northern ecosystem appears to reflect "healthy" limber pine regeneration dynamics, bearing a similar pattern to several studies in stands that are unaffected by WPBR, where recent seedlings are most abundant followed by saplings [12,25,26]. In contrast, the southern ecosystem was more recruitment limited over the past 20 years. While probable, we do not have the data to conclude that WPBR is directly contributing to the observed age structure, as it was beyond the scope of this study to collect seedling infection levels and relate infection to seedling survivorship. Our analyses identified that including disease in model selection improved the model fit and helped explain the absence of seedlings in heavily diseased stands. In range wide monitoring of WPBR in Alberta in 2009, seedling mortality was commonly found in infected stands [17], and may prevent successful regeneration.
One of the most intriguing and alarming findings of our study is that seed production processes do not ensure regeneration in the southern ecosystem. We found that 3.02 times more seedlings occurred in the northern ecosystem, despite only 43.5% of the net cone production observed in the southern ecosystem [4], leading to a 7-fold seed:seedling ratio difference between ecosystems. While greater rates of limber pine cone predation by red squirrels occur in the southern ecosystem, more cones are still available in this ecosystem because the habitat supports fewer squirrels, and squirrel populations display temporal fluctuations that allow cone escape in mast years [5]. The remaining difference in seedling abundance could be accounted for by the incremental contributions of substrate and microsite covariates which enhanced regeneration, and their greater prevalence in northern ecosystems, and lower levels of conditions that reduced regeneration (disease and competition). Even after accounting for the variety of processes found to influence regeneration in our final model, regeneration was consistently predicted to be lower in the southern ecosystem than in the northern ecosystem ( Figure 2). These findings suggest that other processes are limiting regeneration, which could range from differences in grazing regimes (the southern ecosystem has cattle grazing), post dispersal seed predator interactions, or other abiotic factors. Clearly, self-maintaining regeneration is possible in mature stands in the northern ecosystem, even with low seed production. While our final model indicates that regeneration will increase on average with increasing seed availability (BA limber pine), model predictions in Figure 2a suggest that greater seed production in the southern ecosystem does not appear to be adequate to compensate for a combination of disease and vegetation, or site conditions that are unfavorable in many plots in the southern ecosystem.
Given the higher seed production in our southern ecosystem, we may infer from other studies [28] that nutcrackers would be more attracted to these high cone producing regions, so the low levels of regeneration observed were unlikely to be attributable to either seed or disperser-limitation. We observed many seemingly favorable microsites that remained unoccupied, despite considerable caching activity. Further declines in seed production through WPBR induced mortality and reduced seed output in infected trees will only exacerbate in situ conservation of existing stands. Such limber pine stands may be beneficial for nutcracker assisted dispersal potential to new habitats created through fire; however, similar explorations in diseased whitebark pine systems have not led to adequate regeneration in some studies [57,58]. The relationship between fire and limber pine regeneration needs further study, particularly in northern portions of the species range, where it has only been tested in one study [59]. In our severely infected southern ecosystem, we suspect that this is a short-term recruitment opportunity with no lasting population-level benefit, given the prevalence of disease and the mortality that seedlings may experience in ensuing years. Our results suggest that when in situ conservation is prioritized, diseased stands will require greater restoration planting of WPBR-resistant seedlings because natural regeneration may not adequately restock declining stands.

Conclusions
We found support for an inclusive regeneration model that identified seed availability, substrate, microsite, and ecosystems as contributing covariates to regeneration counts, and disease and competition as covariates that helped explain the absence of regeneration. Across both ecosystems, limber pine regeneration patterns remain responsive to a variety of favorable biotic and abiotic factors. As hypothesized, our more diseased southern ecosystem had less regeneration, but there was little evidence of seed limitation caused by white pine blister rust as a contributing factor given the higher cone productivity per tree and higher seed tree basal area. Thus, despite seed availability, other biotic or environmental factors are hindering regeneration processes in our southern ecosystem. While age structure differences between ecosystems suggests that WPBR induced mortality of seedlings contributes to this finding, further study of the timing and causes of seedling mortality is needed to confirm this interpretation. Dispersal limitation may still occur if Clark's nutcrackers use southern landscapes less than the northern ecosystem; however, this is unlikely, given the much higher seed production in the southern ecosystem. Differences in regeneration may be partially attributable to suitable microsites and germination substrates being more available in the northern ecosystem, and less competing vegetation cover. Overall, seedling densities in our northern ecosystem remain comparable to studies in other parts of limber pine range, where WPBR has not been documented. Age structure of seedlings suggests that our northern populations may be self-sustaining at present, while southern populations are not. Recent increases in WPBR at the northern limits of WPBR, combined with lower natural seed production, suggests that less natural regeneration and more WPBR-induced seedling mortality will occur in the future. This leads us to process-based recommendations for both monitoring the status of limber pine populations, and the implementation of recovery plans in jurisdictions where it is endangered. In ecosystems where natural regeneration remains responsive to a variety of biotic and abiotic factors, we recommend that attention to monitoring recent regeneration (<50 cm) is warranted, as we believe it may be a sensitive indicator of changes in the capacity of an ecosystem to sustain natural regeneration. In ecosystems where low seedling densities and the ratios of young to older seedlings suggest recruitment limitation, we recommend that recovery efforts prioritize restoration with WPBR-resistant seedlings because natural regeneration cannot be relied on to restock declining stands.
Author Contributions: V.S.P. conceived, designed, acquired funding, administered, and collected the data for the study. V.S.P and D.R.V. analyzed and visualized the data, and wrote the paper.