Climate change effects on understory plant phenology: implications for large herbivore forage availability

Consistent with a warming climate, the timing of key phenological phases (i.e. phenophases) for many plant species is shifting, but the direction and extent of these shifts remain unclear. For large herbivores such as ungulates, altered plant phenology can have important nutritional and demographic consequences. We used two multi-year datasets collected during 1992–1996 and 2015–2019 of understory plant phenology in semi-arid forested rangelands in northeastern Oregon, United States, to test whether the duration of phenophases for forage species has changed over time for three plant functional groups (forbs, graminoids, and shrubs). Duration of spring green-up was approximately 2 weeks shorter in the later years for forbs (19 ± 3.8 d) and graminoids (13.2 ± 2.8 d), and senescence was 3 weeks longer for graminoids (25.1 ± 5.1) and shrubs (22.0 ± 4.6). Average peak flowering date was 3.1 ± 0.2 d earlier per decade for understory forage species with approximately 1/3 of the species (35%) exhibiting earlier peak flowering dates over time. Variation in late-winter precipitation had the greatest effect on the duration of understory green-up, whereas variation in summer precipitation had a greater effect on duration of the senescent period. Collectively, these results indicate climate-related progression towards shorter periods of peak plant productivity, and earlier and longer periods of plant senescence, the combination of which substantially reduces the temporal window of forage available in growing forms most usable to herbivores. This work adds a needed component to the climate change literature, by describing links between shifting climate variables, multiple phases of understory plant phenology, and possible nutritional consequences for herbivores under a warming climate.


Introduction
Shifting plant phenology is an important biological indicator for understanding the effects of climate warming on biomes worldwide.In addition to being an indicator of climate change, shifts in plant phenology can influence a variety of ecological processes and species interactions (Inouye 2008, Enquist et al 2014, Morellato et al 2016).Changes in key phenological events can result in mismatches between plants and co-occurring species (Kharouba et al 2018) but the direction and extent of these shifts remain unclear.As climate change affects the phenology of plant species around the world (Parmesan andYohe 2003, Wolkovich et al 2012), gaining a better understanding of the onset and duration of these events will be crucial to understand how the timing of resource availability may impose pressures on interacting species.
For many large herbivores, such as ungulates, plant phenology is a key component that influences the availability of forage on the landscape, with green-up being the most studied phenological period (i.e.phenophase).Previous studies have found that spatiotemporal variations in green-up can strongly influence the migratory paths of herbivores (Aikens et al 2017), affecting both fitness and survival (Hurley et al 2014, Middleton et al 2018).Studies that focus on later phenophases (i.e.senescence) or multiple sequential phenophases are uncommon (Augspurger and Zaya 2020) but could have important implications for herbivore nutrition and behavior.Since maternal body stores of capital breeders (e.g.elk, Cervus canadensis) are accumulated late in the previous growing season and are used to compensate the high energetic demands of early reproduction (Gustine et al 2017), phenological changes to forage as it matures can affect the nutritional resources available for reproduction (Cook et al 2004, Monteith et al 2013).Plant senescence can also trigger behavioral cues for migratory ungulates to leave their summer range and begin fall migration (Rivrud et al 2016).Thus, it is critical to study all events in the annual life cycle of plants that represent different nutritional costs and benefits.
Consistent with climate change predictions, temperature changes have been correlated with shifts in plant phenophases (Parmesan and Yohe 2003, Calinger et al 2013, Willis et al 2017, Chmura et al 2019).But exploring additional variables, such as precipitation, may be especially important in water-limited systems where climate change is predicted to result in more frequent droughts (Crimmins et al 2010, McLaughlin et al 2017).As phenophases are sequential, changes in response to climate in one phenophase may affect subsequent phases (Augspurger and Zaya 2020).The duration of each phenophase is delineated by an onset and end date.If shifts in the same direction for the onset and end dates are analogous, little to no change in duration occurs (i.e.consistent shift, figure 1(a)).Alternatively, if the onset or end dates change in opposite directions over time (i.e.inconsistent shift, figure 1(b)), or if the end dates change for early phenophases (figure 1(c)), the duration of subsequent phases can either shorten or lengthen.One or more changes in phenophase duration may result in a specific phase not providing favorable conditions for herbivores seeking to maximize nutrition.For example, warmer springs can result in an earlier onset of plant growth, but also shorten the period of highest quality forage, resulting in reduced net energy gain and overall animal performance (Pettorelli et al 2007, Monteith et al 2015).Understanding the phenology of forage resources across multiple phases provides a comprehensive characterization of the current and future nutritional conditions under which herbivores are exposed.
The goal of our research was to quantify long-term changes to three understory plant phenophases (green-up, vegetative, and senescent) as proxies of forage availability for large herbivores.Our objectives were to: (a) document shifts in the onset and duration of forb, graminoid, and shrub phenophases between two study periods (1992-1996 and 2015-2019); (b) quantify the effects of temperature and precipitation on the duration of understory phenophases; (c) evaluate changes in peak flowering date of preferred forage species for large herbivores (e.g.elk and/or mule deer).We predicted that the onset of all three periods would occur earlier in the season over time (e.g.consistent shift, figure 1(a)) and that the duration of peak green-up would decrease while the duration of peak senescence would increase (e.g.inconsistent shift, figure 1(c)).We hypothesized that flowering dates of forage species would become earlier over time, signaling that forage availability shifts earlier in the season.Finally, we predicted that the onset and duration of peak phenophases would depend largely on precipitation and temperature, with warmer and drier conditions producing earlier and contracted phenophases.

Measurement of phenological events
We conducted research at the Starkey Experimental Forest and Range (hereafter, Starkey;Rowland et al 1997, Wisdom 2005) located in the Blue Mountains of northeast Oregon, USA (45 • × 13 ′ N, 118 • × 31 ′ W).The 83 km 2 Experimental Forest and Range has been the site of long-term studies of wild herbivores and domesticated cattle during the past 50 years (Skovlin 1991, Rowland et al 1997).Wild ungulates including mule deer (Odocoileus hemionus), white-tailed deer (O.virginianus), and elk were present at the time of the study.Cattle (Bos taurus) were grazed 15 June-15 October each year.Elevations ranged from 1120 to 1500 m with rolling uplands and steep canyons.
We recorded plant phenological estimates at 11 sites across Starkey from 1992 to 1996 and 2015 to 2019 (figure S1 in supporting information).Sites were chosen to capture the three major land cover types (grassland n = 4, open pine n = 4, and mixed conifer n = 3) in the study area.Land cover types correspond to generalized vegetation types and associated plant communities in the Blue Mountains Ecoregion (Johnson and Clausnitzer 1992, Johnson and Swanson 2005, Powell et al 2007, Powell 2011) which vary with the wide range of environmental gradients that occur over small spatial scales.The primary drivers for variation in land cover types are related to slope, aspect, soil depth, soil type, and elevation.Forests in the region contained a mosaic of Douglas fir (Pseudotsuga menziesii), grand fir (Abies grandis), ponderosa pine (Pinus ponderosa), lodgepole pine (Pinus contorta), western larch (Larix occidentalis), and Engelmann spruce (Picea engelmannii).Each site consisted of a 50 m permanent transect with four half-meter plots located at 10 m increments.Observations were made at 14 d intervals starting in early April to late October of each year.During some years, site visits were delayed due to snow conditions.On average, phenological data were recorded 14 times throughout the season across all 44 plots (table S1).Phenophase was recorded as green-up, vegetative, or senescent for three understory plant functional groups: forbs, graminoids, and shrubs.The graminoid functional group included all grasses, sedges, and rushes.(a) Green-up (peak nutritional quality of forage for ungulate herbivores; Pettorelli et al 2007) was classified as the onset of leaf greening, young leaves, and increased leaf size; (b) vegetative (peak forage biomass) was classified as fully green leaves, elongated stems, and the presence of flowers or seeds; and (c) senescent (decline in forage quality and quantity) was classified as total leaf senescence, loss of pigment, leaf drop, and cured plant parts.The phenophase of each plant functional group was described using ocular estimates at the plot level by categories of percentage understory canopy cover (1 = 0%-20%, 2 = 21%-40%, 3 = 41%-60%, 4 = 61%-80%, 5 = 81%-100%).
We focused our analysis on the timing of the peak phase (i.e.>60% of plants in a single phenophase) because this metric likely captures the amplitude of the event while minimizing variation across species.Because we assessed phenophase categories at 14 d intervals, the onset and end dates of peak phenophase events were estimated.For 'start' date of peak phenophase, we used the mean interval between the first date on which >60% (categories 4-5) and the date of the most recent prior survey in which a category of <60% was recorded.For example, if green-up for a functional group was first marked as category 4 on 15 April and in the preceding survey on 1 April was marked as category 1, 2, or 3, the estimated onset date of peak green-up was 8 April.Similarly, we used the average of the last date when >60% (categories 4-5) was recorded for a functional group and phenophase, and the date of the next survey in which <60% was recorded for that group and phase to derive the 'end' date for peak phenophase.We summarized the phenophase data in each plot from the start to end date (in number of days).
Linear mixed-effects models were used to test whether the duration of the three phenophases for each functional group changed over time.Fixed effects included land cover type (grassland, open pine, and mixed conifer) and period (early: 1992-1996, late: 2015-2019).We fitted all models with random effects, reflecting the spatial arrangement of the plots nested in transects.Given the sequential (yearly) sampling within each period, we added a residual correlation structure to the model to account for temporal autocorrelation (Zuur et al 2010).We calculated confidence intervals (CIs; 95%) of each model's fixed effects factors and p-values for pairwise comparisons for all combinations (adjusted with the Tukey method) using the multcomp (Hothorn et al 2008) and lsmeans (Lenth 2016) packages in program R (R Core Development Team 2021).

Temperature and precipitation data
We obtained yearly data for temperature in degrees Celsius (average) and rainfall in millimeters (total amount) from weather stations located on Starkey to evaluate long-term climate trends.To evaluate growing season temperatures and precipitation at each transect location, we used monthly 4 km resolution raster data from the parameter elevation regression on independent slopes model (PRISM Climate Group 2004).PRISM provides a set of time series and long-term mean gridded surfaces of multiple climate variables to reveal short-and long-term climate patterns in the United States (Daly et al 2008).We used monthly PRISM precipitation data for the early period (1992)(1993)(1994)(1995)(1996), late period (2015-2019), and 30 years averages  to estimate the growing season precipitation patterns.We calculated precipitation percent of average (PPA) to identify drier or wetter conditions during the sampling periods (Brown et al 2019).PPA indicated the observed percentage of total precipitation for a month in comparison to the 30 years average total precipitation for that month, where PPA i is PPA in month i of a given year; PPT i is the total precipitation in month i in a given year; and PPT i,mean is the 30 years mean value of the total precipitation in month i for each grid cell.We used this equation to produce monthly PPA layers for each year of the study.
To explore the relationship between PPA and temperature, we intersected the mean monthly temperature ( • C) and PPA at all plot locations.Because the onset of understory green-up and flowering is likely to be affected by winter precipitation in the form of accumulated snow (Lambert et al 2010), we calculated average winter PPA as the period from January to April.For the senescent period, we calculated an average summer PPA as the period from June to August.We averaged temperatures beginning with the month prior to the median day of onset through the month following peak observations for each phenophase (green-up: March-May; vegetative: May-July; senescent: June-August).
We used linear mixed-effects models to predict how temperature and precipitation altered the phenophase duration.Our model included plant functional group, mean temperature, PPA, and all two-way interactions as main fixed effects with plots nested in transects as random effects.We standardized the climate variables by subtracting the mean of the variable and dividing by the standard deviation of the variable.We assessed variance inflation factors among climate variables (Zuur et al 2010).This model did not include land cover or time-period fixed effects because we were incorporating variation across land cover types and years to focus on relationships between phenology and climate variables.We tested climate variables for multicollinearity with variance inflation factors using the performance package (Lüdecke et al 2021) in program R prior to running models.We built models using the lme4 package (Bates et al 2015) and summarized using restricted maximum likelihood t-tests via Satterthwaite approximations for degrees of freedom and estimated significance levels using Wald χ 2 tests in the car package (Fox et al 2019) in program R.

Changes in forage flowering dates
For analyses of how flowering dates changed over time, we focused on plant species that have been previously documented in the literature as preferred forage species for large herbivores (e.g.elk and/or mule deer) in the western US (table S2).We recorded species presence and reproductive information at each plot including the date of flowering and total number of flowers present.We counted individual flowers or flowering stems where individual flowers were too small to count.Peak flowering was defined as the date associated with the highest total number of transect-level flowers/species.
We supplemented survey data on flowering with US Forest Service herbarium records collected in the Blue Mountains Ecoregion.The herbarium data (dating back to the 1950s) allowed us to assess flowering trends before the start of the Starkey phenophase surveys.We limited herbarium records to the county (Union), township (4S, 3S) and range (34E) of Starkey and to records that included both a collection date and digitized herbarium image.We confirmed flowering specimens by visually examining herbarium images.We defined forb and shrub specimens with >50% of buds in anthesis as peak flowering at the time of collection (Primack et al 2004).Given the large number of diminutive flowers for graminoid specimens, we could not feasibly count the total number of flowers in digitized records.We confirmed the presence of graminoid flowers in the image and used the collection date of the specimen as a proxy for peak flowering date (Munson and Long 2017).We removed duplicate records (i.e.same date and location) and merged records for synonymous species names.We converted the herbarium collection date and plot survey date to Julian date (1-365), which served as our response variable and a proxy for approximate peak flowering (hereafter 'flowering') date.We used linear regression to analyze change in peak flowering date for each species (response variable) over time.Year was a continuous predictor.The regression was considered significantly different than zero (p ⩽ 0.05), marginally significant (p ⩽ 0.10), or nonsignificant (p > 0.10).

Results
Late period temperature and precipitation differed from the early period and 30 years averages.Average springtime temperatures (1 March-31 May) were slightly higher during the late period (figure 2; 6.6 ± 0.23 • C), but these differences were not significant (early = 6.3 ± 0.24 • C, p = 0.24).By contrast, we observed significantly higher summer temperatures and drier conditions in the late period than in the early period.Summer (1 June-30 August) mean temperatures were 2 • C higher (figure 2, t 459 = −7.80,p = <0.005)and late period precipitation was lower than the 30 years averages (table S3).For example, from 2015 to 2019, average August precipitation (11 mm) was almost half of the 30 years average (22 mm).During the early period, winter (January-April) precipitation totals were similar to the 30 years average, whereas winter precipitation in the late period increased (table S3).
Across all land cover types, mean duration of peak green-up was significantly longer during the early years for graminoids (13.2 ± 2.8 d) and forbs (19 ± 3.8 d) but not for shrubs (figure 3, table S4) Differences in duration were driven by changes in end dates, as the onset dates remained consistent (table S4).The duration of the peak vegetative phenophase was considerably, and significantly, longer during the early versus late period across all land-cover types for all three plant functional groups (figure 3, table 1).Mean duration of the peak vegetative period was ∼1 month (39.3 ± 4.8 d) longer for graminoids during the early years compared to the later years (table S4).For forbs and shrubs, peak vegetative period was 18.2 ± 3.4 and 28.9 ± 3.4 d longer respectively, during the earlier years compared to the later years.In contrast, the peak senescent period was significantly shorter during the early years for graminoids and shrubs, but not for forbs (figure 3, table 1).Peak senescent phenophase was 25.1 (±5.1) days shorter for graminoids and 22.0 (±4.6) days shorter for shrubs during the early years (table S4).
Variation in winter precipitation had a greater effect than variation in mean spring temperature on the duration of green-up.The standardized coefficient for winter PPA (−15) was almost twice as large as the mean spring temperature coefficient (−8; table 2).Green-up duration was longer during years with drier than average winter and spring conditions.Graminoid and forb green-up durations were longer with cooler spring temperatures, whereas shrub green-up duration was associated with warmer spring temperatures (figure 4).Variation in late-spring temperatures (May-July) had a slightly greater effect than variation in winter precipitation on the duration of the vegetative phase; however, we did not find a significant interaction between PPA and mean temperatures (table 2).Plant functional group responses were similar across the three groups with cooler late-spring temperatures and drier than average winters lengthening the vegetative period (figure 4).During the senescent period, variation in summer precipitation had a greater effect than the variation in mean summer temperatures on duration of phenophase (table 2).Summer PPA had a larger standardized coefficient (−14) than the mean summer temperature (−1) and a strong negative correlation across plant functional groups, where senescence was extended during drier than average summers (figure 4).Plant functional group responses to mean summer temperatures differed slightly; graminoid and shrub senescence was extended with warmer summer temperatures, but the duration of forb senescence was unaffected (figure 4).
From the combined herbarium records and phenology plot observations, we examined 1487 specimens from 26 species.Changes in peak flowering date over time showed a high degree of variability among the herbivore forage species (figure 5).On average, peak flowering date occurred 3.1 ± 0.2 d earlier in the year per decade across all species.Nine species (35%) exhibited significant negative shifts in peak flowering date (i.e.flowered at an earlier date over time).Significant negative shifts were observed in 71% of graminoids and 24% of forbs.We did not observe significant shifts for any shrub species.Among significantly responsive species, Danthonia unispicata, a native perennial graminoid, showed the greatest change with flowering occurring 9.4 d earlier per decade.

Discussion
There is strong evidence that global climate change affects the start and end dates of phenological phases across a wide range of plant communities (Oberbauer et al 2013).Our research complements these studies by demonstrating the value of using multiple phenophases across plant functional groups to describe the effects of climate change on understory forage.While investigating the duration of three important phenophases (i.e.green-up, vegetative, and senescence) we found evidence of an inconsistent shift in duration whereby some plant functional groups experienced an expansion of one phase (i.e.senescence) and a protraction of other phases (i.e.green-up and vegetative).On average, the onset of green-up was not significantly different between time periods; however, the end dates of peak green-up for graminoids and forbs have moved forward by approximately 2 weeks over the past 27 years.These results suggest that the timing of green-up has not changed; however, the duration of both green-up and vegetative growth has significantly decreased over time, resulting in earlier senescence compared to a more prolonged period of growth that extended later in the season during the mid-1990s.We found consistent results in the direction of change (e.g.longer senescent period) across land cover types, indicating that the changes in phenophase duration are strong and thus detectable in the face of other variations across the landscape.
Phenophase responses varied according to plant functional groups.Phenological responsiveness to climate variables, such as increased temperature, can be highly variable between related species (Abu-Asab et al 2001).We discovered significant differences between time periods for all three graminoid phenophases, whereas the duration of shrub green-up and forb senescence was not significantly different between time periods.Herbarium data indicated that graminoid species were more likely to flower earlier in the late Table 1.Linear mixed-effects parameter estimates for main effects of phenophase duration for peak green-up, peak vegetative, and peak senescent periods for each plant functional group (FunGroup).Wald χ2 tests and P values P(χ2) were calculated for effects of period (early, late) and land cover type (grassland, mixed conifer, and open pine).Random effects included plot nested in transect.Probabilities <0.05 are designated in bold type.In other semi-arid, montane sites, the onset of early season phenophases, such as flowering, can be triggered by snowmelt date (Diez et al 2012) and warm spring temperatures (Lesica andKittelson 2010, CaraDonna et al 2014).In our system, shorter green-up periods for forbs and graminoids were associated with warmer springs and wetter than average winters.Because most of the precipitation in this system falls in the form of snow during late winter through early spring, a high PPA value denotes higher than average accumulated precipitation.Similar to green-up, the vegetative phase was lengthened with cooler late-spring temperatures and a drier than average winter, whereas variation in summer precipitation had a greater effect than mean summer temperatures on the duration of peak senescence.Future work should focus on improving the understanding of the joint relationships between temperature and precipitation on forage conditions when evaluating the effects of climate change on wildlife (Aikens et al 2020).Besides local temperature and precipitation, long-term, regional climate patterns could impact plant phenology (Ehleringer and Sandquist 2018).Metrics such as the Pacific Decadal Oscillation and El Niño Southern Oscillation could be used to further investigate patterns of understory plant phenology and may be an important tool for predicting forage availability in the future.
The continued increase in mean temperatures expected for our study area (figure S2; Halofsky et al 2018) suggests that we will observe continued declines in the duration of the green-up and a temporal compression of vegetation growth of herbaceous and woody species into even shorter periods, coupled with longer periods of plant senescence.This combination of changes will ultimately result in reduced forage quantity (i.e.shorter growing seasons result in less biomass) and duration of high-quality forage for herbivores, both of which could reduce the body condition and population performance of large herbivores (Cook et al 2004, 2013, Middleton et al 2013).Although migratory species have movement strategies that optimize responses to rapidly changing conditions, for example, by moving to higher-quality habitats, partially migratory or resident species may face nutritional constraints due to a shortened green-up period (Hebblewhite et al 2008).Shifts in the duration of phenophases can also produce mismatches in the timing of nutritious forage and important life history events (Monteith et al 2015, Ross et al 2017), which could have critical implications on the reproductive performance of animal populations (Saino et al 2011, Kerby andPost 2013).During the spring, wet and cold conditions improve forage quality by increasing nitrogen content while dry and hot conditions reduce nitrogen and digestibility (Lenart et al 2002).A contraction of green-up, when forage quality is the highest, can reduce net energy gain and lower the population performance of herbivores (Albon and Langvatn 1992).
Additionally, an earlier onset of plant senescence could truncate the availability of late-summer forage, which is a significant component for building fat reserves going into winter (Cook et al 2004, Rolandsen et al 2017).In this scenario, large herbivores must rely on adaptive behavioral strategies to meet their nutritional demands.Wildlife managers should acknowledge the potential for shorter growing seasons and the associated effects on population performance (e.g.reduced pregnancy rates and recruitment) and animal distributions (e.g.shifts of large herbivores to irrigated agricultural lands) when managing large herbivores.Future work should quantify how altered phenophases influence forage biomass, quality, and temporal availability, and subsequent impacts on herbivore fitness and behavior to assess how such changes might impact herbivore population performance.

Conclusions
Future climate projections in the Pacific Northwest indicate sustained warming in the Blue Mountains region (Halofsky et al 2018), and may result in increased intensity and duration of summer drought conditions (McLaughlin et al 2017, IPCC 2021).Sustained alterations in hydrological conditions can alter the growing regimes of understory plant communities and decrease productivity in dry upland forests in the Pacific Northwest (Latta et al 2010, Kerns et al 2018).Given the strong relationship between PPA/summer temperature and duration of forage senescence, we can expect that these plant functional groups' phenologies will shift even more in the future.By expanding our understanding of how different plant groups respond over time to changing conditions we are better able to predict which fauna may be affected and through what mechanisms.

Figure 1 .
Figure 1.A schematic framework of potential changes in duration of phenophases for understory plants over time.If both onset and end dates shift in the same direction, we expect to see (a) little to no change in duration.If onset or end dates change in the opposite direction, or if only onset or end date changes, we expect to see inconsistent shifts in duration of phases (b) and (c).(Dotted lines represent future, and solid lines indicate present conditions.).

Figure 3 .
Figure 3. Estimates and 95% confidence intervals of duration of each phenophase period by major land cover type.Letters represent statistically significant differences for time period (early, late) across each phenophase (green-up, vegetative, and senescent) with the same letter indicating that time period contrasts are not significant (Tukey contrasts, p-value threshold of 0.05).Shrubs were only present at the mixed conifer and open pine sites.

Figure 4 .
Figure 4. Interaction plots of plant functional groups for the duration of (a) green-up, (b) vegetative, and (c) senescent phenophases.Predicted values of plant functional group by mean temperatures and plant functional group by PPA.Confidence bands correspond to a 95% width.

Figure 5 .
Figure 5. Changes in peak flowering date of preferred forage species (see tableS2in supporting information) over a 64 years period.Symbols represents the significance of the slope of a line (days/year) from linear regressions of peak flowering over time.Significant shifts are symbolized with squares (p ⩽ 0.05), nonsignificant shifts with triangles (p > 0.10), and marginally significant shifts with circles (0.05 < p ⩽ 0.10).Whiskers represent 95% confidence intervals.The dashed vertical line represents no change over time.Species are listed alphabetically within plant functional groups (forbs = green, graminoids = orange, and shrubs = purple).Sample size (n) of observations (i.e.herbarium specimens, field observation) is given following species name.
, with over half of the graminoids exhibiting significant advancement.Our results suggest that grazers (i.e.ungulate herbivores that primarily consume grasses) are especially vulnerable to the altered growing seasons of different graminoids species, some of which are important summer forage species for large period

Table 2 .
Linear mixed-effects parameter estimates for main effects and all two-way interactions to predict duration of peak phenophase at Starkey Experimental Forest and Range.Wald χ2 tests and P values P(χ2) were calculated for plant functional group (FunGroup), mean temperatures (MeanT), and percent precipitation average (PPA).Random effects included plot nested in transect.Probabilities <0.05 are designated in bold type.