Soil nutrient availability affects tundra plant community composition and plant–vole interactions

ABSTRACT

Climate change in the Arctic is predicted to drastically alter carbon and nutrient pools, plant communities, and plant–animal interactions. We examined how four levels of long-term (16 years) nutrient addition in moist acidic tundra affected plant community structure and the abundance of Eriophorum vaginatum, the preferred forage for tundra voles (Microtus oeconomus). We also explored how simulated press and pulse herbivory for four years affected Eriophorum at these different nutrient levels. Shifts in plant community structure, reductions in species diversity and richness, and decreased Eriophorum abundance were pronounced at the highest levels of nutrient addition. Eriophorum abundance was negatively correlated with deciduous shrub abundance, likely caused by light limitation in fertilized plots. Added nutrients generally increased Eriophorum leaf length, but effects varied among years and simulated herbivory treatments. After being clipped once, Eriophorum leaves regrew in fertilized plots to the same length as leaves that had not been clipped, demonstrating recovery. Our results show how interannual variation in plant growth and animal activity might exacerbate or dampen responses over multiple years. Changes to plant community composition and plant–vole interactions associated with increased soil nutrients resulting from warming could have cascading impacts on arctic ecosystems and carbon cycling.

KEYWORDS:

• Arctic
• climate change
• plant-animal interactions
• plant community structure
• plant diversity

Introduction

Arctic ecosystems are critical regions of the Earth for numerous reasons, including storing substantial amounts of carbon (C) in their soils (Gorham 1991; Zimov et al. 2006). Climate change threatens this vital C sink as warming is occurring in the Arctic at three to four times the global average rate (Arctic Monitoring and Assessment Programme 2021; Rantanen et al. 2022), which could have far-reaching consequences for regional and global C cycles. Increased permafrost thaw (Zimov et al. 2006) and subsequent deepening of the active layer (Pithan and Mauritsen 2014) will support greater microbial activity and increase C loss from arctic soils to the atmosphere as a by-product of respiration (Schuur et al. 2009; Natali et al. 2011). Additionally, this increased decomposition of organic material will mineralize inorganic nutrients and thus increase soil nutrient availability (Schuur et al. 2015).

Woody plant species have historically made up a relatively small fraction of plant abundance in low-stature arctic tundra plant communities. However, rising temperatures and the resulting increased thaw depth and soil nutrient availability promote woody shrub growth (Sturm, Charles, and Tape 2001; Stow et al. 2004) leading to increased shrub abundance across many areas of the Arctic (Elmendorf et al. 2012; Kleinhesselink, Magnoli, and Cushman 2014; Myers-Smith et al. 2019). Areas with greater shrub abundance have been shown to trap more snow in the winter, leading to increased winter soil temperatures, and conversely, decreased soil temperatures and potential reductions in permafrost thaw in the summer due to greater shading effects (Myers-Smith et al. 2011; Heijmans et al. 2022). This growth of woody shrubs (C fertilization hypothesis) may increase above- and belowground C storage in the Arctic (Camill et al. 2001). Thus, greater shrub abundance can affect multiple ecosystem processes, including carbon cycling.

Over time in nutrient addition experiments in Alaskan moist acidic tundra (MAT), the community shifts to being dominated by erect deciduous shrubs following an initial increase in abundance of most species (Shaver et al. 2001). Two species [Betula nana, a deciduous shrub (hereafter, Betula), and Rubus chamaemorus, a forb] dominate, reducing species diversity and richness and resulting in greater aboveground net primary productivity (ANPP; Hobbie, Gough, and Shaver 2005; Shaver et al. 2014). The species that decrease in abundance include Eriophorum vaginatum (hereafter, Eriophorum), a tussock-forming sedge that is a vital source of food and shelter for voles and plays an important role in regional C cycling (Curasi et al. 2022).

Plant–herbivore interactions can be critical in structuring plant communities in northern regions and might affect how vegetation responds to climate change (e.g., Christie et al. 2015; Gough and Johnson 2018; Roy et al. 2020; Mekonnen et al. 2021; Koltz, Gough, and McLaren 2022). Small herbivores like tundra voles (Microtus oeconomus) or lemmings (Dicrostonyx spp.) can exhibit strong top-down controls on arctic ecosystems by altering plant species composition (Tuomi et al. 2018), and their nonconsumptive activities can affect local nutrient availability (e.g., Roy et al. 2022). Unlike other herbivores in the Arctic, voles can exert year-round impacts on their local environment, as they do not hibernate in the winter. Voles also have high metabolisms and low digestive efficiencies, meaning they must consume several times their body weight in food per day (Ims and Fuglei 2005). These rodents tend to have three to four year boom–bust population cycles (Ims and Fuglei 2005; Olofsson, Tommervik, and Callaghan 2012) that can stochastically and sometimes dramatically impact plant communities and regional C cycling (Min et al. 2021; Rastetter et al. 2021). During a boom, rodents can cause a significant amount of damage to plants through herbivory and other activities. Following a boom, plant biomass can be as much as 12 to 24 percent lower than in the previous year, potentially altering C stocks at a regional scale, the effects of which can be observed from space via satellite imagery (Olofsson, Tommervik, and Callaghan 2012). This short-term, high-level herbivory pressure, followed by several years of sustained low-level herbivory, can be considered a ‘pulse’ event. By contrast, a ‘press’ event would occur if there were several years of sustained high rodent populations and accompanying herbivory pressure.

Climate change is expected to negatively affect small rodent populations (Batzli and Henttonen 1990) and might alter historic boom–bust population cycles as has been seen for some lemming populations (e.g., Schmidt et al. 2012). In the winter, voles and lemmings rely on the safety and insulation of the subnivean zone (the area between the top of the snowpack and the surface of the ground) for survival. The subnivean zone allows these small rodents to freely move around and access their food while under a layer of protection against harsh winter weather conditions like frigid temperatures, wind, and rain. For these reasons, they favor areas that accumulate a thicker snowpack in the winter (Fauteux et al. 2015; Reid et al. 2012); depending on the rodent species, greater shrub abundance may improve winter habitat for the animals by trapping more snow (Steketee et al. 2022). As the Arctic continues to warm, winter conditions will become less favorable by becoming warmer and more variable (Ims and Fuglei 2005). Warming temperatures contribute to increased rainfall that can reduce snowpack depths by melting existing snow, changing snow density, and decreasing snow accumulation. Reduced snowpack depth increases the risk of small mammals being exposed to extreme winter weather conditions, which could be detrimental to their survival (Aars and Ims 2002; Pitelka and Batzli 2007).

Tundra plant species vary greatly in their susceptibility to herbivory by small mammals. Eriophorum is a well-studied tussock-forming sedge that accounts for a significant portion of the biomass in MAT (Shaver, Chapin, and Gartner 1986) and contributes to regional C stocks (Curasi et al. 2022). Tundra voles rely on this species for food and shelter, often selectively foraging for sedges over woody shrubs (Batzli and Lesieutre 1991; Christie et al. 2015). Eriophorum tolerates herbivory by increasing leaf production after a defoliation event (e.g., Archer and Tieszen 1983). However, tolerance to herbivory and recovery after defoliation might be reduced when subjected to fertilization (Johnson and Gough 2013), likely because of a shift in competitive interactions with other plant species (although also see Grogan and Zamin 2018). While Eriophorum can tolerate both nutrient addition (Heskel et al. 2012) and warmer temperatures (Tissue and Oechel 1987), this species can be outcompeted by other plants, especially woody shrubs, under those conditions (Heskel et al. 2013; Shaver et al. 2014), potentially because of reduced light availability (Fetcher 1985). Fertilization enriches Eriophorum leaf nitrogen (N) content, making it more palatable, and increases foraging benefit, which can lead to increased herbivory rates (Gough, Ramsey, and Johnson 2007; Gough et al. 2008; 2012; Tuomi et al. 2018). This suggests that climate change might lead to compounding negative effects on Eriophorum’s photosynthetic capacity and overall success by increasing interspecific plant competition for light availability and herbivory rates.

The goal of this study was to examine how MAT plant community composition and Eriophorum’s response to simulated rodent herbivory are affected by soil nutrient availability. Unlike most previous studies in Alaskan MAT, here we examined multiple levels of soil nutrients. Based on earlier studies, we hypothesized (H1) that with increasing levels of nutrient addition, plant communities will have more pronounced shifts toward the dominance of deciduous shrubs (specifically Betula) and reduced Eriophorum abundance. We quantified this by examining relative cover, species richness, and diversity over 16 years in plots receiving one of four different levels of nutrient addition. We also hypothesized (H2) that fertilization would not affect the ability of Eriophorum to recover from simulated vole herbivory. We tested H2 by implementing a 4-year simulated herbivory experiment across this nutrient addition gradient with both pulse (clipped once) and press (clipped annually) herbivory treatments.

Methods

A nutrient addition gradient experiment in MAT was established by the Arctic Long-Term Ecological Research project in 2006 at Toolik Field Station, located in the northern foothills of the Brooks Range, Alaska, USA (68.63, 149.60). Soils in this region are moist and mainly glacial till and outwash from glacial advances approximately 60,000 years before the present (Hobbie and Kling 2014). The region is snow-covered for two-thirds of the year, and the maximum thaw depth ranges from ~0.25-1.0 m during summers. The growing season in this region is short, typically ranging between 120-126 days, starting around mid-May until mid-September (Cherry et al. 2014). Over the past several decades there has been no consistent warming trend documented in this area (https://www.uaf.edu/toolik/edc/). Tussock-forming sedges (Eriophorum vaginatum), deciduous shrubs (Betula nana), ericaceous evergreens (Rhododendron tomentosum and Vaccinium vitis-idaea), and mosses (Sphagnum spp.) dominate the vegetation at this site (Shaver et al. 2014, Table S1). Most of the previous nutrient addition experiments conducted on MAT have included only one high level of fertilization (N 10 g/m²/yr and P 5 g/m²/yr) (Shaver et al. 2014). The MAT nutrient addition gradient used in this study was established to monitor the response of the plant community to a range of nutrient addition levels described below.

Nutrient addition

The MAT nutrient gradient experiment consists of blocks of 5 × 20 m plots that are assigned different fertilization treatments (Figure 1). Control plots (CT) receive no fertilizer, F2 receives 2 gN/m²/yr + 1 gP/m²/yr, F5 receives 5 gN/m²/yr + 2.5 gP/m²_/yr, and F10 receives 10 gN/m²/yr + 5 gP/m²/yr (Prager et al. 2017; Suchocki 2020). Plots are fertilized annually with granular sodium nitrate, ammonium chloride, and triple superphosphate at the beginning of the growing season when the plots are free of snow (prior to 2012, N was applied solely in the form of ammonium nitrate). All these nutrient amendments are likely greater than the amount of additional N being mobilized by increased rates of decomposition resulting from warming in this plant community [estimated by Prager et al. (2017) as ~0.3 gN/m²/yr].

Figure 1. Nutrient addition gradient study with simulated herbivory experimental design in three blocks, each containing four 5 × 20 m fertilization plots (CT, F2, F5, and F10). Each fertilization plot has two no-herbivory (NH), two pulse (PU), and two press (PR) individual tussocks.

Plant community sampling

In mid to late July from 2008-2021 (except for 2009), relative vascular plant abundance was estimated to species or genus in the four nutrient level treatments in blocks 1-3 of this experiment. The relative abundance of mosses, lichens, and nonliving groups, including plant litter, bare ground (including frost boils and standing water), and woody standing dead plant tissues (mainly deciduous shrub die-off), were also estimated. A permanent 8 × 1 m long monitoring area was established in 2008 on the northern side of each plot. Along this area, eight 1 m² quadrats were aligned roughly 0.2 m in from the edge of the plots, and aerial percent cover was estimated for each species/group (only four quadrats were used per plot in 2008, and five quadrats per plot in 2015 and 2020). To determine relative abundance, the percent cover of individual species/groups was divided by the sum of all percent cover within a quadrat. Plants were identified to species, then consolidated into growth forms (deciduous shrubs, evergreen shrubs, graminoids, and forbs; see Table S1) except for Eriophorum, which was excluded from the graminoid group and evaluated independently.

Diversity was calculated with the data at the species level using the Shannon Diversity Index (Equation 1(1) $\mathit{H}=-\sum _{\mathit{i}=1}^S{p}_i\mathit{ln}{p}_i$(1) ). Species richness and diversity were calculated at a 1 m² level from the cover quadrats.

Equation 1(1) $\mathit{H}=-\sum _{\mathit{i}=1}^S{p}_i\mathit{ln}{p}_i$(1) : Where H is diversity, i is one species, p_i is the proportion of the entire community made up of species i, S is the total number of species in the plot.

(1) $\mathit{H}=-\sum _{\mathit{i}=1}^S{p}_i\mathit{ln}{p}_i$(1)

Simulated herbivory

In 2018, we established a simulated herbivory experiment within the fertilization gradient to examine the effects of leaf removal on Eriophorum (Suchocki 2020; Williamson 2022). Six Eriophorum tussocks with no discernable signs of vole herbivory were randomly selected in every plot, except in the F2 treatment in Block 3 where there were no tussocks in the sampling area (N = 66). Two Eriophorum tussocks were randomly assigned to one of three herbivory treatments (NH: no herbivory, PU: pulse or PR: press). This provided a total of six replicates for each of the twelve possible combinations of fertilization and simulated herbivory treatments, except for the F2 treatments that only had four replicates (Figure 1). In late June 2018, the leaves (and any other species of plant growing on the tussock) of each tussock in a pulse or press treatment were cut to the base of the leaf and removed. This process was repeated in 2019, 2020, and 2021 for the press treatment, while the pulse treatment was allowed to recover.

Six tillers (ten in 2021) from each Eriophorum tussock were randomly selected for measurement. Selected tillers were fitted with metal identification tags, given unique identifiers, and monitored for the duration of the growing season. If a tag fell off, a new tiller was randomly selected to maintain a total of six tillers (ten in 2021) per Eriophorum tussock. Tags were removed from press treatments before leaves were cut, and new tillers were selected following each simulated herbivory event. The length of the longest green leaf within a tiller was measured for all treatments every 2 weeks for the duration of the growing season (about 8 weeks), except in 2020, when leaf lengths were collected only once at peak season (in late July) because of staffing and travel limitations during the COVID-19 pandemic.

Statistics

Statistical analysis was performed with R version 4.0.3 (2020). Statistical significance was evaluated at an alpha level of 0.05. A multivariate analysis of variance was used to test significant differences among mean relative abundance of live plant growth forms and Eriophorum (excluding plant litter, bare ground, and woody-standing dead) among fertilization treatments for 2008 and 2021 (run independently) to provide a comparison between the start of the experiment and 16 years of treatment. Similarly, species richness and diversity were analyzed separately for 2008 and 2021 as a one-way analysis of variance (ANOVA) to compare treatment effects in these 2 years. Linear regressions were used to test the relationships between the relative abundance of Eriophorum and deciduous shrubs among the fertilization treatments. Two-way ANOVA with interaction was used to test differences in peak season leaf lengths among simulated herbivory and fertilization treatments (each year tested independently of the others). If the ANOVA showed significant differences among treatments or the interaction, the means were compared using Tukey’s honestly significant difference.

Results

Plant community response to nutrient addition

In 2008, after 2 years of nutrient addition, there was little difference in mean relative abundance for most functional groups across nutrient addition treatments (Figure 2). However, moss and Eriophorum had significantly different mean relative abundances across the fertilization treatments (Table S2). Eriophorum abundance was greater under the highest levels of nutrient addition (F5 and F10) compared to the CT and F2 treatments, while moss had the opposite trend (Figure 2).

Figure 2. Mean (n = 3) relative cover of growth forms/species among years (2008–2021) and fertilization treatments (see Figure 1 for nutrient levels). Relative cover was not measured in 2009, and F2 treatments were not sampled in 2010 and 2011.

Fourteen years later (2021), much of the plant community had shifted in response to fertilization, most dramatically in the F10 treatment (Figure 2, Table S2). Although we did not measure productivity or biomass in this study, NDVI (Normalized Difference Vegetation Index) data suggest that aboveground biomass of vascular plants was greatest in the F10 treatment (Williamson 2022; Figure S1), like results from older experiments at this site in which ANPP significantly increased in response to F10 treatments (Shaver et al. 2014). As in 2008, the relative abundance of moss in 2021 continued to be significantly lower in the fertilized treatments compared to the control (Table S2). Eriophorum, however, had the opposite trend, seen in 2008, tending to have a lower relative abundance in the fertilized treatments compared to the control. Evergreen shrubs also had significantly lower abundance as fertilization increased, while lichen and graminoids (not including Eriophorum) were not affected (Figure 2). Deciduous shrubs (of which Betula comprised more than 70 percent with almost 100 percent in F10) and forbs (dominated primarily by Rubus) were the only living growth forms to increase in abundance in response to nutrient addition, although this was only significant for deciduous shrubs. By 2021, deciduous shrubs made up over half of the relative abundance in the F10 plots, compared to about a quarter in F10 plots in 2008 (Figure 2). No differences in relative abundance of litter, woody standing dead, or bare ground were detected across fertilization treatments or years.

In 2008, mean species richness was similar across the control and all fertilization treatments (Table 1; fertilization treatment not significant). By 2021, richness had declined in the F5 and F10 treatments, while the control and F2 treatments remained relatively unchanged (main effect of fertilization: F_3,8 = 25.2, p < .001; Table 1). This effect was most pronounced within the F10 treatment, which saw a loss of over half the number of species documented in 2008. Several plant species went locally extinct in the F10 treatment by 2021, all of which occurred in low relative abundance in control plots (Table S1), including two deciduous shrub species (Arctostaphylos alpina and Vaccinium uliginosum), one evergreen shrub (Cassiope tetragona), and one forb (Pedicularis sp.). Shannon diversity had similar trends where there was no difference in diversity across the nutrient addition treatments in 2008, but diversity was significantly lower in the F5 and F10 treatments compared to CT and F2 in 2021 (Table 1, main effect of fertilization in 2021: F_3,8 = 36.1, p < .0001).

Table 1. Mean (n = 3; ± SD) species richness and Shannon diversity index per 1 m² by fertilization level (see Figure 1 for treatment abbreviations) for the years 2008 and 2021. Superscript letters indicate means were significantly different within year.

	CT	F2	F5	F10
Richness
2008	12 ± 1	13 ± 2	12 ± 2	12 ± 2
2021	12 ± 1^a	11 ± 2^a	8 ± 1^b	5 ± 1^b
Diversity
2008	2.24 ± 0.10	2.23 ± 0.18	2.13 ± 0.20	2.05 ± 0.22
2021	1.98 ± 0.13^a	1.76 ± 0.08^a	1.18 ± 0.24^b	0.63 ± 0.21^b

Negative correlations between deciduous shrub (dominated by Betula) and Eriophorum relative abundance occurred across all treatments; however, they were only significant under the highest levels of nutrient addition, F5 and F10 (Figure 3). Notably, the range of relative deciduous shrub abundance increased with the level of fertilization (0.15–0.75 in the F10 treatment compared to 0.14–0.34 in the CT). In the F10 treatments, relative deciduous shrub abundance in 2008 (solid triangles; Figure 3) was between 0.2 and 0.3, while in 2021 (solid circles) shrub abundance was always over 0.5. A similar, but less dramatic trend was seen in the F5 treatment. This increase in relative shrub abundance over time in the F5 and F10 treatments is responsible for driving the negative relationship as Betula abundance increased and Eriophorum abundance decreased. In the CT treatment, values from 2008 and 2021 are clustered together, indicating little to no change in deciduous shrub or Eriophorum cover over time.

Figure 3. Correlation between deciduous shrub and Eriophorum relative cover in each fertilization treatment. Each point represents the average for one block (n = 8 quadrats within block) for 3 points per treatment per year (2008–2021). Confidence interval around regression line is SE.

Eriophorum recovery from herbivory

In general, leaf growth was measured across the growing season as either an increase or no change compared to the first measurement. In some cases, mean leaf length fluctuated within a growing season because of new leaves being selected if tags fell off or the tip of a leaf started to senesce, reducing the live length (Figure 4).

Figure 4. Mean (± SD, n = 3 except in F2 n = 2) Eriophorum tiller heights (cm) among simulated herbivory and fertilization treatments (see Figure 1 for treatment abbreviations) over four growing seasons (2018–2021). Note that tiller heights were only measured once in 2020 due to COVID-19 restrictions.

The seasonal pattern of growth in the 3 years for which multiple dates were sampled generally did not differ across fertilization treatments. The tussocks that were not subjected to clipping (NH) showed a positive effect of fertilization on leaf length in 2018, but in the other years leaf lengths were lower overall, and there was little effect of fertilization. The PU treatment had longer tillers in the F5 and F10 treatments across all years, whereas the PR tussocks were generally not able to regrow following clipping in any of the fertilization treatments (Figure 4).

The effects of the treatments on leaf length at the final measurement from each growing season differed in 2018 and 2021. In 2018, a marginally significant interaction effect was found between clipping and fertilization (F_6,53 = 2.21, p = .06, Figure 5). This occurred because the NH leaves grew longer with greater nutrients, the PU tillers in F5 and F10 also grew longer than the other two fertilization treatments, but the PR leaves were not affected by fertilization. Leaf lengths were significantly different among the simulated herbivory treatments as well (main effect of clipping: F_2,53 = 60.0, p < .0001). PR and PU were significantly shorter than NH (p < .001 for both treatments) but were not different from each other (p = .26). Leaf lengths were also significantly affected by fertilization (main effect of fertilization: F_3,53 = 5.82, p < .001). F10 and F5 fertilization treatments had significantly longer final measurement leaf lengths compared to the control (p = .01) but were not significantly different from each other or the F2 treatment (Figure 5).

Figure 5. Final tiller heights (n = 3 except n = 2 in F2) among simulated herbivory and fertilization treatments for four growing seasons (2018–2021; treatment abbreviations in Figure 1). Xs indicate means; letters indicate significant differences among the means caused by a significant interaction between herbivory and fertilization in 2021.

In 2021, the fourth year of the simulated herbivory experiment, a significant interaction effect was found between clipping and fertilization (F_6,51 = 3.30, p = .008, Figure 5). Fertilization resulted in different leaf lengths within the PU treatment but not the NH or PR treatments. In the PU herbivory treatment, Eriophorum in the F10 fertilization plot produced significantly longer leaves compared to CT and F2, but none were significantly different from F5. In addition, the PU leaves in F5 and F10 grew to similar lengths as the NH leaves in all fertilization treatments, demonstrating recovery (Figure 5). PR Eriophorum produced significantly shorter leaves compared to all other treatments because of annual clipping. Leaf lengths were also found to be significantly different among simulated herbivory treatments in 2021 (main effect of clipping: F_2,51 = 90.6, p < .0001). PR and PU were significantly shorter than NH (p < .001 for both treatments) but were not significantly different from each other (p = .59). Leaf lengths were also found to be significantly affected by fertilization (main effect of fertilization: F_3,51 = 2.98, p = .04). Eriophorum in F10 fertilization plots produced significantly longer leaves compared to the CT (p = .04), but not when compared to any other fertilization treatment.

Discussion

Plant community composition shifted toward Betula dominance and decreased Eriophorum in response to fertilization

The MAT plant community was significantly altered by nutrient addition over the course of 16 years, with a notable decrease in relative Eriophorum abundance. The magnitude of the shift depended on the level of fertilization, and thus supports our first hypothesis. Our results agree with previous studies that used the F10 level of nutrient addition in MAT and showed an increase in ANPP resulting from the increase in Betula (e.g., Shaver et al. 2001; Mack et al. 2004; Hobbie, Gough, and Shaver 2005). By 2021, the relative abundance of deciduous shrubs (driven by the increase of Betula) was significantly greater in the F5 and F10 treatments compared with the CT and F2 treatments. Added nutrients resulted in significantly more soil ammonium and nitrate availability in the F5 and F10 plots in 2021 compared with CT and F2 (M. Wilcots, unpublished data), and this difference in available N may be responsible for the more dramatic vegetation shifts in the higher nutrient treatments. Betula can quickly accumulate a significant height advantage over most other MAT species when fertilized because of rapid branching and sturdy wooden structures that support increased canopy area (Bret-Harte et al. 2001). Betula also appears to have a leaf-level physiological advantage over Eriophorum under both ambient and increased soil nutrient conditions (Heskel et al. 2013; Prager et al. 2020).

In response to fertilization, the only other growth forms/groups to increase in relative abundance were forbs (primarily Rubus) and standing dead (also dominated by Betula, particularly in the F10 where Betula comprised 97 percent of standing dead cover in 2021), although this change was not statistically significant for forbs. Rubus has been found to respond positively to nutrient addition in other MAT experiments (Shaver et al. 2014), creating an understory layer beneath the taller Betula shrubs. Greater standing dead Betula in the F5 and F10 plots likely resulted from winter frost exposure as the shrubs extend above the snow or disease, although the exact cause remains uncertain. This die-off suggests that Betula may not be able to maintain long-term dominance even with greater soil nutrient availability. Most of the other growth forms (evergreen shrubs, Eriophorum, lichen, moss, and bare ground) decreased significantly in response to fertilization likely caused by increased shading from deciduous shrubs, resulting in the species richness and diversity decline in MAT over time.

Although only significant in the F5 and F10 treatment, Eriophorum abundance was negatively associated with Betula abundance, and this relationship intensified as the nutrient addition continued over time. Eriophorum abundance was higher in the F10 treatment in 2008 compared with the control, reflecting a short-term positive growth response to added nutrients before the shrubs were able to grow tall enough to outcompete this species (also see Shaver et al. 2001). Eriophorum growth was reduced when artificially shaded (Chapin and Shaver 1985), and this species increased productivity and biomass more under fertilization when Betula was removed (Bret‐Harte et al. 2004, 2008). At a landscape scale, a recent remote sensing study found that across Alaska and part of northwestern Canada, deciduous shrubs are increasing at the expense of graminoids (Macander et al. 2022), and another study determined that Eriophorum is decreasing in abundance across much of Alaska (Curasi et al. 2022), which could also be correlated with increased shrub abundance.

At the level of the individual tiller in the NH treatment, the effect of fertilization on Eriophorum leaf length was mixed. In the first year of the experiment (2018), fertilization caused an incremental increase of leaf length with added nutrients as expected (Figure 5), consistent with similar studies that reported greater leaf lengths under fertilization (Shaver, Chapin, and Gartner 1986; Shaver and Laundre 1997; Johnson and Gough 2013) and in warming treatments (e.g., Hobbie, Shevtsova, and Chapin 1999). However, in the other three years, growth of NH tillers was not consistently greater with fertilization. Eriophorum might have been increasing its production in other ways, such as increasing tiller number (Shaver and Laundre 1997). In addition, because of the heterogeneity of tussocks occurring on the landscape, biomass of Eriophorum often does not significantly differ between control and fertilized plots (e.g., Shaver et al. 2014), so, although leaf length is increasing in fertilized plots, without the establishment of new tussocks, biomass of this species at the community level may not significantly change while relative abundance declines because of the dramatic absolute and relative increase of Betula.

Another explanation for the interannual variation in leaf length we observed is that variables such as weather, herbivory, or interspecific plant competition are interacting with fertilization to prevent a positive growth response to additional nutrients. For example, Eriophorum can experience reduced growth when exposed to higher air temperatures but can have increased growth when exposed to both elevated temperatures and fertilization (Chapin and Shaver 1985), suggesting that growth responses might differ depending on the environmental conditions. As discussed in the previous section, woody shrubs often grow to dominate over other species in fertilized plots and can lead to understory shading that limits Eriophorum growth. In addition to increased competition for light, natural vole herbivory on Eriophorum has been found to increase with fertilization (Gough et al. 2012). Estimates of vole activity and vole damage to Eriophorum tussocks in the treatment plots sampled here did not detect significant differences across soil nutrient levels (data not shown); effects were heterogeneous and patchy even in relatively high vole abundance years. Snap trap and mark-recapture estimates suggest that vole populations in this area were low in 2017 and 2018 and higher in 2019–2021 (J. Steketee, personal communication, 10 October 2023). This suggests that natural vole activity could potentially have affected our results. The combined negative impacts on Eriophorum growth from plant competition and increased vole activity could lead to suppressed leaf growth responses in fertilized plots.

Nutrient addition does not always support Eriophorum recovery from vole herbivory

In the fourth year of the simulated herbivory experiment, there were inconsistent responses in Eriophorum leaf length to simulated herbivory and nutrients. In addition, the pulse and press treatments resulted in different responses. For the press treatment, our results supported our second hypothesis that added nutrients did not affect Eriophorum recovery from repeated clipping. In 2018, after a single year of clipping, there was a suggestion that leaves on tussocks in the higher fertilization plots were able to recover marginally faster and achieve slightly longer lengths by peak season, mimicking the pattern in the no herbivory treatment, but there were no significant differences among treatments. After 4 years of simulated herbivory in the press treatment, there were no differences in tiller length among fertilization treatments, and on average tillers remained less than 10 cm in length compared to around 15 cm in the no herbivory plots. The results of the press treatment demonstrate that the negative impact of repeated clipping could not be overcome by increased nutrients. In other studies where leaves were clipped multiple times within a growing season, a threshold has also been demonstrated beyond which Eriophorum does not have the resources available to regrow (Chapin and Slack 1979) even when receiving nutrients at the F10 level (Johnson and Gough 2013). Thus, if herbivore pressure were to become more consistent across years, Eriophorum abundance would likely decline regionally.

In the PU treatment, however, our second hypothesis that greater nutrients would not affect Eriophorum growth following simulated herbivory was not supported after 3 years of recovery. In 2018, following the initial clipping treatment, leaves in the F5 and F10 treatments recovered faster and achieved greater lengths than those in the F2 and CT plots. In 2021, the pulse F5 and F10 leaves had recovered to approximately 15 cm in length, similar to leaves not subjected to clipping. However, the CT and F2 leaves were still shorter, demonstrating that leaves that were not able to take advantage of the higher levels of nutrients could not recover from biomass removal 4 years earlier. A similar result was seen 1 year after clipping when Eriophorum was fertilized with just N, and shoot production was still lower than leaves that had not been clipped (Grogan and Zamin 2018). These results suggest that Eriophorum requires significant nutrients to recover from a single clipping event after 4 years, which is an approximate cycle for voles in this region. In a similar study conducted in MAT at the F10 level after 10 years of fertilization, Eriophorum was also not able to recover from clipping after 3 years (Johnson and Gough 2013). In that study, adding nutrients without clipping did consistently increase leaf length. These results suggest that weather, as well as the surrounding plant community, are likely affecting the ability of Eriophorum to regrow.

Conclusion

The composition of MAT plant communities and specifically the abundance of Eriophorum are affected by soil nutrient availability (Shaver et al. 2014). Under the greatest levels of nutrient addition, Betula dominates the community, resulting in an overall loss of diversity and reduction in Eriophorum relative abundance. The levels of nutrient addition applied here are higher than would be likely to naturally occur with warming (Prager et al. 2017). Our results suggest a similar community trajectory toward dominance by Betula at lower levels of soil nutrients than F10, mirroring results from a passive warming experiment in MAT (Sistla et al. 2012). As the Arctic continues to warm, and soil nutrient concentrations subsequently rise, the MAT plant community will likely become woodier, as is already occurring in many regions (e.g., Myers-Smith et al. 2019). Depending on the scale and magnitude of these plant community composition shifts, they could potentially impact the regional carbon balance as well as small mammal herbivores like voles.

It is unclear how a warming Arctic will impact tundra vole populations. Voles rely on a thick snowpack to remain insulated from harsh winter weather conditions. Greater woody shrub abundance could therefore benefit voles in the winter, as the taller canopy provided by the shrubs can trap snow and create a thicker snowpack (also see Steketee et al. 2022). Voles might initially experience a dietary benefit as increased soil nutrient availability contributes to greater forage quality. However, in this study we documented that at the same time as deciduous shrubs increase in abundance, tundra voles’ preferred forage species (Eriophorum) decreases. Our results demonstrate that Eriophorum growth does not rapidly recover from herbivory events even when soil nutrients are drastically increased. Therefore, if greater forage quality promotes increased herbivory rates, it could further accelerate the decline of Eriophorum, and thus negatively affect vole populations over the long term. The role and degree of impact that small mammals will have on arctic plant communities and the C cycle are complex but have the potential to be significant.

Acknowledgments

The authors thank many field assistants for their help in collecting these data and maintaining these experiments, particularly through the COVID-19 pandemic restrictions, including Ruby An, Lilly Conteh, Audrey Fatone, Savannah Kjaer, Peter Lamousin, Megan Muller-Girard, Nicholas Patel, Miriam Ritchie, Sunniva Sheffield, and Julian Sheppy. We thank Jim Laundre and Gus Shaver for their long-term maintenance of the Long-Term Ecological Research experiments, and Rebecca Rowe, Jennie McLaren, Kevin Griffin, Natalie Boelman, and Ed Rastetter for their collaboration and comments on this article. We also thank Megan Wilcots for sharing unpublished data, and the anonymous reviewers for comments that improved the article.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Data availability statement

The data that support the findings of this study are openly available in the Arctic Long-Term Ecological Research Data Catalog at https://arc-lter.ecosystems.mbl.edu/arctic-data-catalog.

Supplementary material

Supplemental data for this article can be accessed online at https://doi.org/10.1080/15230430.2024.2356276

Additional information

Funding

This work was supported by the U.S. National Science Foundation under Grants 9810222, 1026843, 1637459, and 1603760.