Local landscape position impacts demographic rates in a widespread North American steppe bunchgrass

landscape position impacts demographic rates in a widespread North American steppe bunchgrass. Abstract. Understanding the environmental drivers of demographic rates and population dynamics over space and time is critical for anticipating how species will respond to climate change. While the in ﬂ uence of temporal environmental variation and large environmental gradients are well recognized, less is known about how local topography and landscape position in ﬂ uence demography over small spatial scales. Here, we investigate how local landscape position (north- vs. south-facing aspects) in ﬂ uence the demographic rates and population growth of a common bunchgrass in western North America, bluebunch wheatgrass ( Pseudoroegneria spicata ), using 6 annual censuses measuring growth, survival, and reproductive output. We found notably lower survival on south-facing slopes, particularly among smaller individuals. In contrast, south-facing slopes maintained comparatively high reproductive output in most years, measured both as spikes per plant and spikelets per spike. When we combined these data in demographic models, we found that lower survival among small individuals and greater reliance on reproduction mean south-facing slopes would generally have to maintain higher recruitment for a stable population. Our results highlight the important in ﬂ uence that landscape position and local topography can have in driving population trends. As conditions warm and dry with climate change (north-faces becoming similar to current south-facing slope conditions), bluebunch wheatgrass may become more reliant on reproduction to maintain viable populations and more sensitive to variability in recruitment.


INTRODUCTION
Understanding demographic responses to environmental variation and the consequences of these responses for population dynamics is central to population ecology, evolutionary ecology, and anticipating the response of species and communities to climate change (Donohue et al. 2010, Doak and Morris 2010, Salguero-Gomez et al. 2012, Ehrlén and Morris 2015. The need to understand and anticipate ecological responses to climate change has led to a growing literature linking climate conditions to demographic rates and population dynamics. This work has expanded our knowledge of how plant demographic rates respond to changes in environmental conditions through time (Dalgleish et al. 2011, Chu et al. 2016, Shriver 2016, 2017, Tomasek et al. 2019) and across large spatial gradients (Angert 2006, Eckhart et al. 2011, Gelfand et al. 2013, Merow et al. 2014. For example, Dalgleish et al. (2011) identified winter and summer climate conditions as important for explaining inter-annual variation in demographic rates of steppe plants. Similarly, Eckhart et al. (2011) found that spatial variability in temperature and precipitation across a species' range explained differences in population trends in a California annual plant. While most of this work identifying relationships between demographic rates and environmental conditions has focused on temporal variability or large-scale spatial patterns, the influence of local topography and landscape heterogeneity (e.g., <500 m) on demographic rates is less well understood.
The scarcity of analyses of demographic responses to small-scale topography is surprising given the well-known role topographic position can play in creating microclimate conditions that amplify or ameliorate broader spatial and temporal variability in climate (Bennie et al. 2008, Zellweger et al. 2019. For example, surface temperatures on north-and south-facing slopes can differ 20°C in mountainous terrain at~46°l atitude (Scherrer and Körner 2010). However, the effects of local topographic position on environmental conditions are generally not captured by a single weather station or gridded climate datasets (Zellweger et al. 2019). This lack of small-scale climate data may explain why microclimate and local conditions are often not explicitly considered in demographic analyses.
When demographic analyses do account for microclimate and topographic position, the effects can be striking. For example, Nicolè et al. (2011) found that slope angle was the primary environmental predictor of survival in a rare alpine plant and hypothesized that soil depth and soil moisture drove this effect. Similarly, Dullinger et al. (2004) found that slope was the most important predictor of adult tree mortality at tree line. Finally, Oldfather and Ackerly (2019) found several significant relationships between microclimate conditions and demographic rates, although this did not lead to consistent trends in population growth across microclimates. These examples make clear that if we want to answer pressing ecological questions across spatial scales, understanding how plant demographic rates and population dynamics vary with landscape position and topography will be essential (Gurevitch et al. 2016).
Bluebunch wheatgrass (Pseudoroegneria spicata) is a large perennial bunchgrass common in, and ecologically important to, the Columbia Basin and many other low to mid-elevation areas of the Intermountain West of North America (Rodhouse et al. 2014). In these systems, bluebunch wheatgrass often forms a substantial portion of the herbaceous biomass (Rodhouse et al. 2014). However, these semiarid systems are susceptible to invasion by a number of annual plants, such as cheatgrass (Bromus tectorum), that can come to dominate. Conversion from perennial to annual systems often follows repeated disturbance from livestock overgrazing and fire (Davies et al. 2012). The loss of perennial grasses, including bluebunch wheatgrass, fundamentally alters the structure and function of steppe ecosystems. Communities dominated by invasives provide less forage for livestock, support less wildlife (DiTomaso 2000), and likely sequester less organic carbon (Verburg et al. 2009, Rau et al. 2011. Thus, understanding how patterns of reproduction, growth, and survival of perennial plants vary across landscapes in these systems has important applied implications (Brooks andChambers 2011, Davies et al. 2011). However, in some locations bluebunch wheatgrass can be resilient to disturbances, and when bluebunch wheatgrass or other large deep-rooted perennial bunchgrasses persist at sufficient densities, they can exclude invasive annuals (Rodhouse et al. 2014). The spatial configuration of persistent bluebunch wheatgrass may depend on local topography, with north-facing slopes often supporting robust remnant stands of perennials (Rodhouse et al. 2014). Quantifying differences in demographic rates of bluebunch wheatgrass on north-facing and south-facing slopes are important to understanding how perennial grasslands can resist long-term conversion to invasive annual communities.
Because bluebunch density varies with aspect over small spatial scales, demographic processes are likely to vary over these same scales, but we still lack an accounting of this variability. Variation in demography with aspect is likely to result from divergent temperature and moisture regimes on different aspects. Variation in temperature and moisture correlate with variation in phenology and morphology of bluebunch wheatgrass across larger spatial scales (St. Clair et al. 2013). Further, within a single site, long-term data from this species suggest that annual temperature and moisture drive variation in demographic variables including growth and survival (Dalgleish et al. 2011). Other bluebunch wheatgrass traits may also be impacted by microclimatic variation, for instance, seed quality and offspring fitness (Drenovsky et al. 2016). Thus, we have strong a priori expectations that demographic rates may vary with aspect over small spatial scales in this species.
Here, we use data from six annual censuses of bluebunch wheatgrass populations on northand south-facing aspects from a site in eastern Washington to explore how local landscape position influences demographic rates (survival, growth, and reproductive output) and rates of population growth. Specifically, we addressed the following questions: How does aspect influence rates of individual growth, survival, and reproductive output? What are the consequences of these vital rate differences for annual population growth? Finally, because we lack robust estimates of recruitment (seed-to-seedling transition) from our plots we estimate the recruitment rate required for stable population growth on both north and south aspects and compare this result to available estimates from the literature.

Field data collection
We monitored individual bluebunch wheatgrass plants growing on steep hillsides in Spring Gulch on Whitman College's Wallula Gap Biological Station in the Columbia Basin of Washington (46°00 0 N, 118°54 0 W, 360-420 m elevation). This site receives approximately 350 mm of precipitation per year, primarily during the autumn, winter, and spring (Appendix S1: Fig. S1). Summers are hot and winters are cool. Typically, winter highs are above freezing but average lows are below freezing (Appendix S1: Fig. S2). Spring Gulch runs predominately east to west, so most hillsides face north or south. On the north-facing slopes, there is no exposed bedrock and the silty loam soil is free of rocks and appears >2 m deep. The soil is similarly textured on south-facing slopes and is often >2 m deep, but in places is shallower, and bedrock is exposed rarely. However, there was no exposed bedrock on our monitoring locations. North-facing slopes are primarily vegetated by perennial bunchgrasses (besides P. spicata, also Festuca idahoensis, several species of Poa, and others) though rabbitbrush (Ericameria nauseosa and Chrysothamnus viscidiflorus) is common. Big sagebrush (Artemisia tridentata) is rare on both slopes due to recent (2007,2011,2015) and presumably less-recent fires. South-facing slopes are dominated by annuals, especially cheatgrass (Bromus tectorum), annual fescue (Vulpia sp.), and yellow star-thistle (Centaurea solstitialis), in most places. However, bluebunch wheatgrass stands ranging from 0.03 to 0.17 ha are scattered across the south-facing slopes and cover approximately 5-10% of the surface of these hillsides. Within these south-facing stands of bluebunch wheatgrass, cheatgrass is less common and yellow star-thistle is rare. Only one other large bunchgrass, needle and thread grass (Hesperostipa comata), occurs on south-facing slopes, and it is much less common than bluebunch wheatgrass. When measured over approximately 6 months in 2009, soil temperatures at this site were consistently higher on south-facing slopes than on north-facing slopes (Appendix S1: Table S1), and soil moisture was consistently higher on north-facing slopes to a depth of 30 cm (Appendix S1: Table S2).
In the spring of 2011, we established 14 monitoring sites for bluebunch wheatgrass, seven on north-facing slopes and seven on south-facing slopes. We selected sites where bluebunch wheatgrass was growing. We attempted to locate sites relatively far from other sites on the same slope, and we were able to keep all sites >140 m distant from the next nearest except for one pair of sites on the south face which we located 20 m apart because of the rarity of bluebunch wheatgrass stands on these slopes (mean distance to adjacent sites on the same slope = 180 AE 75 m [SD], min = 20 m, max = 306 m). Each site was v www.esajournals.org 10 m long and 1 m wide and marked at the corners by rebar stakes. We designated 1 × 1 m sampling areas every other meter within each site, such that sites contained 5 non-adjacent sampling areas. During 2011 and 2012, we experimented with methods of placing the 1 × 1 m quadrats and methods of relocating individual plants. In 2013, we marked the corners of each 1 × 1 m quadrat with rebar stakes and marked each bluebunch wheatgrass plant with a small uniquely numbered metal tag held in place by a nail in the soil adjacent to the plant. We used a 1 × 1 m quadrat consisting of a metal frame with wire forming a grid of 10 × 10 cm cells. We anchored this quadrat with the permanent rebar stakes at each visit to the plot and used the grid of 10 × 10 cm cells to construct and update maps of the bluebunch wheatgrass plants within the sampling area.
From 2013 through 2018, a group of two to three researchers visited each plot at least once between late May and late June. During that visit, we searched for every bluebunch wheatgrass plant we could find throughout each plot as well as specifically for bluebunch individuals that had been tagged and mapped on previous visits. We marked any new or previously missed plants we discovered, and we removed tags from any plants that were absent for a second year in a row. Bluebunch wheatgrass at our site is entirely caespitose, and so all the stems of a given individual grow close to each other. Typically, the separation between individuals is clear. However, we may have inadvertently lumped two or more individuals together as a single individual on occasion when the distinction between individuals was unclear. For each bluebunch plant we located, we noted survival and counted the number of living culms and the subset of culms that were flowering (spikes). Both of these counts were highly repeatable (culms R = 0.88, Appendix S1: Fig. S3; spikes R = 0.88, Appendix S1: Fig. S4, see Appendix S1 for further details of repeatability calculations). For all flowering individuals, we haphazardly selected five spikes on which to count the number of spikelets. If there were fewer than five spikes, we counted the spikelets on all spikes. Average number of spikelets per spike was moderately repeatable (R = 0.64, Appendix S1: Fig S5). We also measured the height of the tallest culm, but because this measurement changes over the course of our study season due to growth and damage to plants, we decided not to use this measurement in our analyses. Because new plants were observed at very low rates and the effectiveness of searches for new seedlings varied, we lack robust estimates of the rate at which seeds survived, germinated, and produced visible seedlings present at the next census. We refer to this seed-to-seedling transition as recruitment throughout this paper. We discuss our analytic approach to estimate likely recruitment rates in the absence of robust field data on recruitment in the Population modeling section below. Flowering bluebunch wheatgrass is not easily confused with other species at our site, but observers on this project rarely misidentified small non-flowering individuals of other bunchgrass species as bluebunch wheatgrass. Initially, we planned (Parker and Shriver 2018) to eliminate misidentified seedlings by including only individuals that were flowering or had been observed to flower in a past or future year, in our analyses to minimize the chance of misidentification errors. However, only including young plants later observed to flower would bias survival estimates upwards by selectively eliminating plants that died prior to flowering. Further, including only plants that had flowered in the past would have substantially reduced our sample size for small plants. Thus, we opted to retain all plants identified as bluebunch wheatgrass in the field. Since all identification errors occurred on north-facing plots and such errors were highly unlikely on southfacing slopes where species likely to be confused with bluebunch were rare, identification errors should have only impacted survival estimates from north-facing slopes. Further, those survival estimates would have been biased downwards if misidentified seedling grew to be recognizable as other species, and thus the putative bluebunch seedlings were recorded as missing (dead). This bias downward in survival would make estimates from north-facing plots more similar to south-facing slopes, where survival of small plants was lower than on north-facing slopes (see Results). Finally, observers occasionally failed to detect or record a plant that was present. In plots visited twice to assess sampling reliability, 6.2% of plants detected in at least one of two visits were not detected on the other visit (see v www.esajournals.org Appendix S1). Therefore, we only noted mortality when an individual was absent for 2 or more years, and thus, no survival data are available from the final annual transition (2017-2018) and we only estimate vital rates for four annual transitions over six years.
In the summer of 2014, we collected spikes from 60 bluebunch wheatgrass individuals after seed set to determine the typical number of seeds per spike. We collected five spikes from each target plant, and so did not sample from plants with fewer than five spikes. Otherwise, we sampled with the goal of representing the range of plant sizes encountered at the site. Half of all sampled plants were from north-facing slopes and half were from south-facing slopes. All sampled plants were on the same hillsides as our 14 monitored sites, but were outside of and not immediately adjacent to the monitoring sites to avoid impacting seed rain in those sites. For each plant, we counted the number of spikelets across all five spikes, and the number of seeds produced by those spikelets.
Although we hypothesize that variation in moisture and temperature drive variation in demographic performance of bluebunch wheatgrass in this system, we did not attempt to directly correlate demographic variability with climate variables. First, we lacked a sufficient number of years to compare among-year variation in climate with annual variation in demographic performance. Second, we do not have temperature and moisture data at the level of our survey plots, and so we cannot correlate demographic performance with moisture and temperature at the plot level. However, we are confident that comparing demographic rates between north and south-facing slopes, which have striking differences in moisture and temperature (Appendix S1: Tables S1, S2), provides important insights into the potential roles of environmental conditions in this system.

Population modeling
To understand how individual vital rates combined to influence population growth, we developed integral projection models for north and south-facing slopes. Integral projection models are discrete-time, continuous-state structured population models that allow us to explore how differences in vital rates across aspects impact population growth (Easterling et al. 2000). The model is formulated as where n(x, t, a) is the density of size x individuals at time t on aspect a, K(y, x, t, a) is a kernel describing the transition of size x individuals to size y. In our case, plant size is ln(number of culms). K(y, x, t, a) itself is made up of vital rate functions, Kðy,x, t,aÞ ¼ Gðy,x, t, aÞ Â Sðx, t,aÞ þ Fðx,t, aÞ where G(y, x, t, a) is a kernel describing the size transition of existing individuals, S(x, t, a) is a function of survival rates for individuals of size x, and F(x, t, a) is a function of reproductive output of individuals (i.e., new plants produced per individual). Although F(x, t, a) can include data on both individual reproductive output (seeds produced) and germination and survival of seedlings (recruitment of seeds to seedlings detectable by researchers), as noted above we do not have sufficient field data on recruitment from our plots. Thus, our F(x, t, a) include only the seeds produced per plant. We wish to emphasize that our seed per plant estimates are based on an assumption that the average number of seeds per spikelet counted in 2014 represents the population values across years. Strictly speaking, we expect that this is not true, but our regular seed collections for other purposes suggest that seed number per spikelet do not vary substantially between years.

Vital rate modeling
We estimated vital rates with generalized linear mixed effect models in a Bayesian framework using Stan and the rstan package (Stan Development Team 2020). We parameterized G(y, x, t, a) as a normal kernel, fit to the data with a normal likelihood where the kernel parameters were fit using measured size transitions as Gðy,x, t,aÞ ¼ Normalðyjμ g ðx, t, aÞ,σ 2 Þ μ g ðx, t, aÞ ¼ β 0 þ β 1 x þ β 2 IðaÞ þ γ g ðtÞ þ δ g ðtÞIðaÞ where μ g (x, t, a) is the average size of individual of size x at time t transition to at t + 1. I (a) is an indicator variable for aspect, I = 1 v www.esajournals.org when the aspect is south-facing and 0 when north-facing. γ g (t) is a random effect intercept that varies by year, δ g (t) is a random effect for aspect that also varies by year, β 0 is an intercept, and β 1 and β 2 are fixed effects for size and aspect, respectively.
We estimated survival using a Bernoulli likelihood to observed survival data where the probability of survival, s(x, t, a), is logitðsðx,t, aÞÞ ¼ β 3 þ β 4 x þ β 5 IðaÞ þ γ s ðtÞ þ δ s ðtÞIðaÞ Similar to the size model, γ s (t) is a random effect intercept that varies by year, δ s (t) random effect for aspect that also varies by year, β 3 is an intercept, and β 4 and β 5 are fixed effects for size and aspect, respectively.
Finally, we estimated reproductive output (i.e., seed production) as three parts: (1) the number of seed spikes produced per plant, spike(x, t, a); (2) the average number of spikelets per spike, spikelet (t, a); (3) the average number of seeds per spikelet, seeds.
Fðx, t,aÞ ¼ spikeðx, t,aÞ Â spikeletðt,aÞ Â seeds logðspikeðx,t,aÞÞ ¼ β 6 þ β 7 x þ β 8 IðaÞ þ γ spike ðtÞ þ δ spike ðtÞIðaÞ logðspikeletðt,aÞÞ ¼ β 9 þ β 10 IðaÞ þ γ spikelet ðtÞ þ δ spikelet ðtÞIðaÞ seeds ¼ β 11 Once again, γ spike (t) and γ spikelet (t) are random effects intercept that varies by year, δ spike (t) and δ spikelet (t) are random effects for aspect that also varies by year, β 6 and β 9 are intercepts, and β 7 and β 8 , β 10 are fixed effects for size and aspect, respectively. Because data on the number of seeds per spikelet were not available in each year, we used a single average value for this parameter, β 11 . The number of spikes per plant was fit to the field data using a negative binomial likelihood, while the average spikelets per spike and were fit using a log-normal likelihood. We modeled seeds per spikelet using a normal likelihood, truncated at zero. All random effects were normally distributed with mean 0 and a fit variance parameter.
To implement the integral projection model, we discretized vital rates into 125 size bins using a normalized cumulative normal distribution. Because data on the size distribution of new recruits were unavailable, we assumed all new plants entered the smallest size class.

Model analyses
We quantified a component of the expected annual population growth rate (λ t ) that excludes seed-to-seedling recruitment by calculating the dominant eigenvalue of the discretized integral projection models for each year. It is worthwhile to note that these estimates only represent a component of λ t because our model does not include realistic recruitment rates and instead makes the assumption that every seed becomes 1 established plant. Because the assumption of complete seed survival, germination, and seedling survival is unrealistic, the subcomponent values of λ t are all dramatic over-estimates of true λ t . Because we carried through all of the uncertainty in parameter estimates from our vital rate models to the integral projection model, we present values as posterior distributions.
Finally, we inferred possible seed-to-seedling recruitment rates. We identified the general range of recruitment rates recruitment rates (0-0.02 recruits per seed) that would be needed for both the northern and southern aspects to achieve λ t = 1 (i.e., a stable population) in each year and then calculated population growth (λ t , the dominant eigenvalue) using a numerical approach of calculating λ t across this range of recruitment rates at 0.001 intervals. All newly recruited seedlings entered the smallest size class of the discretized IPM. To determine credible intervals for recruitment rates, we iterated this process over the entire posterior MCMC chains for all vital rates.
To allow transparency in the process of developing our analyses, we pre-registered our analysis plan prior to examination of the data in aggregate (Parker and Shriver 2018). We outline deviations from this plan in the Appendix S1.

Vital rates
Model results indicate that survival and reproductive output consistently differed between north and south aspects, while individual growth did not (Table 1, Figs. 1-3). Both the number of spikes per plant and the number of spikelets per spike were higher on south-facing slope, although the posterior 95% CI for β 8 and β 10 (average difference in spikes and spikelets per plant on south slopes) did overlap with zero. v www.esajournals.org Survival of plants was lower on south-facing slopes. Differences in survival between aspects were most notable in small plants (<20 culms), where the smallest plants on south-facing slopes showed about 50% annual survival and the smallest plants on north-facing slopes showed 70-80% survival. However, large plant survival approached 100% on both aspects (Fig. 1). Although spikes per plant also increased with size, it was in large plants where we observed the greatest difference between aspects. On average, large plants produced about 20-50 more spikes per plant on south-facing aspects than north aspects. Yet, substantially greater reproductive output on south-facing slopes did not occur in all years with individuals on northfacing aspects producing near equivalent numbers of spikes in 1 of the 4 years (Fig. 2). Growth declined on average with size, but showed no consistent differences across aspects. Vital rates also varied substantially from year to year. Most notably, the number of spikes produced per plant tripled on average for large plants from 2015 to 2016 (Fig. 2).

Seed production and recruitment
When we made the assumption that all seeds lead to new recruits, we calculated much higher estimates of a component of population growth rates (i.e., excluding recruitment) on south-facing slopes in 1 of 4 years (2016-2017), and more similar estimates between the two slopes in the remaining years (Fig. 4). Inter-annual variability in vital rates drove up to an almost 2.5-fold variation in the component estimates of population growth rates on south-facing slopes across years; estimates were somewhat less variable across years on north-facing slopes (Fig. 4).
When we set recruitment (seed-to-seedling) rates across a range of plausible values and calculated population growth rates associated with these recruitment rates, population growth rate estimates on south-facing slope were lower across wide ranges of recruitment in three of four years (Fig. 5). As a result, south-facing slopes may need to maintain higher rates of recruitment than north-facing slopes to maintain stable populations (λ = 1; Fig. 5). The recruitment rate required to maintain stable populations varied from year-to-year from <0.001 to > 0.02. However, the recruitment rates required on the southfacing slopes to maintain stable populations were also somewhat more uncertain than those required on north-facing slopes.

DISCUSSION
We found strong evidence that demographic rates in bluebunch wheatgrass, a common bunchgrass in Western North American steppe ecosystems, differed over small spatial scales between north and south-facing slope. While small plants on north-facing slopes survived at much higher rates than those on adjacent southfacing slopes, plants on south-facing slopes had more reproductive output per plant in most years. Despite generally higher seed production per plant on south-facing slopes, the lower survival of small plants means that higher recruitment rates are likely needed for the south-facing populations to remain stable. This may have  (4) 0.01 −0.09 0.14 Notes: Values indicate the difference of south-facing aspects from north-facing ones. Overall fixed effects of aspect for each vital rate (bold) as well as individual year random effects (i.e., yearly deviations from fixed effect) are shown. Parameter values show consistent differences in survival and spikelet number across aspects, modest differences in spike production with considerable annual variability, and no differences in size changes.
v www.esajournals.org major implications for population persistence in the face of climate change given that in semi-arid ecosystems recruitment may be particularly sensitive to increasing temperature and aridity (e.g., James et al. 2019). Although we expect temperature and moisture to be important direct drivers of demographic process, other differences between slopes, such as the degree of competition with invasive annuals, may have also contributed to the divergent demographic patterns we observed. Regardless of the exact ecological mechanisms driving these patterns, our demonstration of demographic variation over small spatial scales contributes to the growing body of work establishing the importance of local topography and landscape position in influencing vital rates and demographic inference.
Our finding that changing landscape position influenced different demographic rates in different ways is consistent with other recent work. For example, Oldfather and Acklery (2019) found notable differences in demographic rates including survival and reproductive output across microclimate conditions. Similarly, others have Fig. 1. Survival probability of bluebunch wheatgrass on north and south aspects across all annual transitions at the Wallula Gap Biological Station, Washington, USA. Points indicate measured field data. Lines represent posterior mean estimates for the average individual survival, with 95% CI shaded regions. For smaller individuals, survival rates were higher on north-facing slopes than on south-facing slopes. Note, points are jittered on y-axis to make it easier to view overlapping points.
found declines in survival with changing slope and aspect, most notably lower survival on steeper, drier slopes (Dullinger et al. 2004, Nicolè et al. 2011. Still, different demographic rates between microsites may not lead to notable differences in population growth if there is compensation among vital rates (Oldfather and Ackery 2019). While we did find evidence of higher reproductive output on south-facing slopes, we hypothesize that this is unlikely to fully compensate for lower survival rates on south-facing slopes compared to north-facing. As we describe below, recruitment rates would probably need to be higher on south-facing slopes than on northfacing slopes to maintain a stable population, but we have reason to expect from ongoing experimental work that bluebunch wheatgrass actually recruits at lower rates on south-facing slopes.
Although we lack estimates of recruitment rates from field data, we do have evidence that changes in recruitment rates could strongly impact population growth. This insight helps us Fig. 2. Spikes produced per bluebunch wheatgrass plant during four growing seasons on north and south aspects at the Wallula Gap Biological Station, Washington, USA. Points indicate measured field data. Lines represent posterior mean estimates for the average number of individual spikes produced, with 95% CI shaded regions. Especially among larger plants (those with more culms), the number of spikes (culms with flowers) increased more rapidly on south-facing than on north-facing slopes in most years.
v www.esajournals.org understand current differences in abundance of bluebunch wheatgrass at our site across aspects and has important implications for the future of these populations. When we made the unrealistic assumption that all seeds led to new seedlings (i.e., the recruitment rate was 1), the incomplete estimates of population growth we generated were similar across aspects, or notably higher for south-facing slopes than north in 1 of 4 yr. But, as we have emphasized, assuming all seeds become new recruits is unrealistic. Preliminary, ongoing experimental work at our site suggests recruitment rates are actually substantially lower on south-facing slopes, where bluebunch wheatgrass is also much less abundant and survival of small plants is lower. Further, based on models in which we systematically varied the recruitment rate, estimates of recruitment needed to maintain stable populations on south-facing slopes were higher than on north-facing slopes in Fig. 3. Size changes for bluebunch wheatgrass on north and south aspects at the Wallula Gap Biological Station, Washington, USA across all annual transitions. Points indicate measured field data. Colored lines represent posterior mean estimates for the average individual size change, with 95% CI shaded regions. Dashed line is 1:1, thus above the dashed line indicates growth, and below the line indicates a reduction in culm number. The pattern is similar on both slopes, with some year-to-year variability, and the best-fit lines are consistent with regression to the mean, with smaller plants more likely to grow and larger plants more likely to shrink. v www.esajournals.org three of four years. This suggests that relatively small differences in recruitment could lead to large differences in population growth in these systems. Since we expect that large differences in recruitment exist between the two slopes, it is plausible that these differences in recruitment rates could be sufficient to lead to lower population growth on south-facing slopes compared to north. For example, when we increased recruitment rates threefold from near 0.005-0.015 Fig. 4. Posterior distributions for an estimated component of population growth rate (λ) on north and southfacing slopes for bluebunch wheatgrass at the Wallula Gap Biological Station, Washington, USA. This is a component of λ because our model does not account for variation in seed survival, seed germination, or seedling survival and instead makes the unrealistic assumption that every seed becomes 1 plant. Because the assumption of complete seed survival, germination, and seedling survival is unrealistic, the depicted values of λ are all dramatic over-estimates of true λ. We depict these values of λ to illustrate the differences between north and south-facing slopes in the components of population growth other than seed survival, germination, and seedling establishment. Each panel corresponds to an annual transition. The variability in estimates within a year and aspect derives from carrying through all uncertainty in parameter estimates from our vital rate models. Estimates of this component of lambda in the absence of germination and establishment data were similar or slightly higher on north-facing slopes in three of four years, and much higher on south-facing slopes in 1 of 4 yr.
v www.esajournals.org recruits per seed in our models, this increased population growth estimates by up to 30% (depending on year and aspect). Given that our ongoing experimental work suggests that the difference in actual recruitment rate could be 10-fold or larger between aspects, the impact of recruitment rate on population growth could be substantial.
The apparent sensitivity of population growth rates to recruitment suggested by our models is consistent with previous work which has established recruitment as a critical process controlling the recovery and restoration of dryland plants and perennial bunch grasses (James et al. 2011), and our estimates of recruitment rates required to maintain stable population are consistent with published rates. For example, James et al. (2019) found consistent probabilities of transitioning from seed to establishment of~0.04. Boyd and James (2013) found similar, but more Fig. 5. Estimated population growth rate in each year from a range of specified recruitment rates for bluebunch wheatgrass on north and south-facing slopes at the Wallula Gap Biological Station, Washington, USA. Lines represent posterior mean estimate for the population growth rate, with 95% CI shaded regions. The confidence intervals for north-facing slopes were moderately narrower than for south-facing slopes. Also, in three of the four years, our population growth estimates for north-facing slopes were higher across the full range of modeled recruitment rates. v www.esajournals.org variable rates. They estimated that germination varied from 20% to 80%, emergence of these germinates was 5-20%, and initial establishment of those that emerged was 20-80%. This would yield recruitment rates ranging from 0.002 to 0.128, which correspond well with our estimated range from 0.001 to >0.02 to maintain stable populations on average. Given the possibility for substantial variability in recruitment rates, the likely sensitivity of recruitment in semi-arid systems to increasingly hotter, drier conditions (James et al. 2011), and invasion of non-native species (Aguirre and Johnson 1991), we expect that variation in recruitment will increasingly impact demography and long-term population persistence of bluebunch wheatgrass in these systems.
Although we do not have detailed measurements of differences in temperature and soil moisture conditions associated with each sample location, other data (Appendix S1: Table S1 and S2) and ongoing experiments from this site as well as other studies indicate that increased temperatures and reduced soil moisture availability, driven by increased solar radiation, may be responsible for vital rate differences. Ongoing experiments performed at this site indicate that when incoming solar radiation is reduced on south-facing slopes (using shade cloth), temperatures decrease and germination and survival of young bluebunch wheatgrass plants increase dramatically (T. H. Parker, unpublished data). Additional studies have also indicated that soil moisture and temperature conditions can substantially alter vital rates of bluebunch wheatgrass. For example, James et al. (2019) found that warmer soil temperatures and reduced precipitation limit the germination, emergence, and early survival of steppe bunchgrasses including bluebunch wheatgrass. Similarly, in a system where snowmelt is an important source of moisture, greater snowpack is associated with greater survival of bluebunch wheatgrass (Dalgleish et al. 2011). Our finding that plot level demographic performance may be quite variable across microsites, even in nearby locations that share similar areal climate conditions, suggests that quantifying microsite climate variation at the scale of plots could provide exciting opportunities to expand our understanding of the relationships between climate and demographic performance across landscapes. Even in the absence of detailed field measurements, models that scale areal climate conditions to local microclimates based on topography are an exciting possibility (Zellweger et al. 2019).
Although north and south aspects differ in temperature and soil moisture, they also support different plant communities, and these community differences are plausible drivers of demographic performance of bluebunch wheatgrass. One of the most obvious differences in plant communities is the higher density of annual plants, including cheatgrass and yellow starthistle, on south-facing slopes. Cheatgrass, in particular, may outcompete bluebunch wheatgrass seedlings and thus hinder their growth (Aguirre and Johnson 1991) and survival, apparently by reducing moisture availability (Harris 1967). In contrast to seedlings, established bunchgrasses apparently outcompete cheatgrass (Rodhouse et al. 2014), thus the negative effects of increased cheatgrass abundance should be most pronounced on smaller, establishing bluebunch wheatgrass individuals that must compete for water in the shallow root zone used by cheatgrass (Young et al. 1987). This scenario is in line with our finding of reduced survival only in smaller individuals on south-facing aspects, although our ongoing shading experiments suggest that direct effects of temperature may be most important. Still, both increased competition and an already warmer, drier climate on southfacing aspects could interact to limit the competitive ability of young bluebunch individuals (Larson et al. 2018).
Although higher temperatures and greater density of invasive annuals may explain the lower survival of small bluebunch wheatgrass plants on south-facing slopes relative to northfacing slopes, the higher reproductive output on south-facing slopes is more difficult to explain. It could be that the lower density of perennial bunchgrasses on south-facing slopes reduces the intensity of competition for soil moisture among well-established individuals and that greater access to resources leads to greater reproductive output. However, various ongoing studies examining the possibility that competition among established bluebunch plants influences the size or distribution of individual plant at this site have revealed little evidence of such competition. Therefore, the driver of higher reproductive output on south-facing slopes remains uncertain.
One of the most striking features of our study site to the casual observer is the lower density of perennial plants, especially bunchgrasses, on south-facing slopes relative to north-facing slopes, and so the differences we observed in bunchgrass demographic rates between these slopes is not a surprise. However, because we lack detailed measurements of recruitment, we remain uncertain whether these differences in population density result from ongoing differences in population growth rate between the two slopes. What we do know, however, is that the hotter, drier conditions on south-facing slopes are associated with reduced survival of small plants on these slopes, and that experimental reduction of temperature on south-facing slopes dramatically increases seedling emergence and survival. Further, our models demonstrate that plausibly low levels of recruitment on southfacing slopes could cause the bluebunch wheatgrass population on those slopes to decline. These observations lead us to hypothesize that lower survival of small individuals coupled with recruitment limitation may be leading to a decline of bluebunch wheatgrass on south-facing slopes and thus may help explain the lower density of bluebunch wheatgrass on these slopes. Further, given forecasts of increasingly warm and dry conditions with climate change in the North American steppe regions (Bradford et al. 2020), we hypothesize that the population decline will worsen, and that as north-facing slopes warm, bluebunch wheatgrass vital rates there may come to resemble what we currently find on south-facing slopes.
We also hypothesize that local landscape position and topography may influence perennial plant populations in semi-arid systems more generally. Given the importance of temperature and moisture to plant growth, survival, and reproduction in these systems (e.g., Badano et al. 2005, Dalgleish et al. 2011) and the large role of topographic position in influencing temperature and moisture (e.g., Bennie et al. 2008), we would be surprised if demographic variability with topographic position were not common in plants in semi-arid systems. Further, we expect variability in recruitment may often play an important role in these systems. Not only did we find clear evidence of the importance of recruitment rate in driving population growth in our system, but there is good evidence that recruitment may be limiting plants in other semi-arid rangelands (e.g., James et al. 2011). Further, when we also consider that topography may drive large differences in recruitment in our system, this suggests that attempting to link variation in recruitment to temperature and moisture, especially as they relate to topographic position, will be important to understanding demographic processes regulating perennial plants in semi-arid systems.

ACKNOWLEDGMENTS
Field work was funded by and the field site provided by Whitman College. Soil samples in 2009 were dried and weighed by THP's undergraduate ecology students. S. Stephan and E. Chin gathered or processed data included in this manuscript. RKS was supported by a NIFA Postdoctoral Fellowship (grant no. 2019-67012-29726/project accession no. 1019364) from the USDA National Institute of Food and Agriculture. Author contributions: THP designed the study and collected data with assistance from EC, HG, AH, SK, MMB, KM, RM, EO, DR, MS, and AW. RKS performed demographic analyses. THP performed repeatability analyses. RKS and THP wrote the manuscript.