Study site

Our study site lies in an extensive fescue (Festuca idahoensis) meadow system located on a broad flat windswept ridge (∼2,332 MASL; 45°46’59′ N, 110°46’40′ W) in the Bangtail Mountains 24 km NE of Bozeman, MT. Our study site (Bangtail Study Area, BTL) was the high elevation grassland site (‘Bridger site’) for the International Biological Program studies (1964–1974) of the grassland biome34. Annual precipitation varies from 700 to1,000 mm, with >80% occurring as snow and the rest as rain in the June through September period. Average July maximum and minimum temperatures are 22 and 8 °C and average January maximum and minimum temperatures are −5 and −10 °C. Although the study plots have not been grazed by livestock since the 1930s, they have been accessible to the full range of large wild herbivores in the area (deer and elk). In 2007, a 1,335 ha portion of the Bangtail Ridge was established as a USDA Forest Service Special Interest Area.

Biomass harvesting

We have harvested total end-of-season (late September to early October) above-ground standing crop in a flat 1,000 m2 portion of this grassland since 1969. It was harvested by clipping all vegetation in two 0.5 m2 quadrats randomly placed within each of five larger (200 m2) subplots (for a total annual harvest of n=10 per treatment). Specific plot locations were not resampled year-to-year. Vegetation was sorted in the field into grass and non-graminoid forb plant functional groups. We did not separate litter from current year’s growth but we took care to avoid litter in the field and previous research showed that carry-over was minimal, particularly in snow treatments. Biomass samples were dried to constant weight and weighed. Early studies14 of vegetative production at BTL showed late growing season biomass to be a good approximation of annual ANPP although it likely represents a slight underestimate. Biomass was harvested in most years from 1969 to 2012 (exceptions were 1975–1978, 1981, 1985, 1986, 1990 and 1995) and continuously since 1996. We used mean ANPP data from the first 6 years of the study reported by Weaver and Collins15. We also compared ANPP during the early part of the record (1969–1974) to ANPP measured over 3 years (1965–1967) in a similar F. idahoensis meadow20 (2,408 MASL) in the GYE ∼90 km southwest of our site. In 1984, we established an independent second set of control plots 60 m northeast of our primary control area and have sampled it since, as described above. For overlapping years, both control areas show qualitatively similar dynamics (Supplementary Fig. 1).

Experimental snow addition

A snow manipulation experiment was established at the site in 1968. The three treatments consist of control plots (described above) and snow-supplemented plots generated with snow fences (n=5) positioned perpendicular to prevailing winds; these double (1.2 m, ‘ × 2’) and quadruple (2.4 m, ‘ × 4’) ambient (<0.6 m) snow levels. Two rectangular areas with continuous snowpack roughly equivalent to the snowfence height over ∼25 m × 30 m area for the × 2 treatment and 50 m × 30 m area for the × 4 treatment resulted.

Chemical and isotopic analysis

Biomass samples were harvested and archived for all treatments for sample years 1987, 1992, 1996, 1998, 2002, 2003, 2007, 2010, and 2012. For these years we randomly selected grass and forb samples (n=3–5) from each treatment and analysed them for chemical and natural abundance isotope composition (n=173). Forbs were not separated into annual and perennial species. The samples were homogenized, ground and analysed for %C (%), %N (%), δ13C (‰) and δ15N (‰) by isotope mass spectrometry at the Woods Hole Biological Station, MA. Biomass %C did not vary significantly (P>0.3) across time, between plant functional types or among treatments. We therefore applied the mean %C (44±2%) to estimate annual C stocks over time. We analysed 13C distributions for patterns of heavy isotope discrimination (Δ=[δ13C air −δ13C plant ]/[1+δ13C plant ]) as a proxy for possible changes in WUE22. To account for the changing 13C signature of atmospheric CO 2 , we used the average 13CO 2 signature for June–September months spanning the 1987–2012 period recorded at La Jolla, CA and Point Barrow, AK, USA (Scripps Institution of Oceanography).

Atmospheric nitrogen deposition

We analysed long-term patterns in regional atmospheric N inputs by first examining the long-term record of annual inorganic N deposition recorded by the National Atmospheric Deposition Program at Yellowstone Park-Tower Falls (WY08) monitoring site located 106 km south of BTL. We also examined March snowpack chemistry from three nearby sites reported by the United States Geologic Survey snow chemistry monitoring network from 1993 to 2012. We analysed data from the three nearest stations: Red Mountain (100 km due west), Big Sky (50 km southwest) and Daisy Pass (145 km south east).

Long-term climate records

We related climate and grassland ANPP using five climate records: (1) hourly precipitation and temperature recorded (1969–2012) at the MSU weather station operated by the Optical Remote Sensor Laboratory in Bozeman (MSU; 24 km southwest of BTL); (2) hourly precipitation and temperature recorded (1969–2012) at Bozeman Yellowstone International Airport (BZN; 29 km west of BTL); (3) daily November–April snow-water equivalent recorded (1974–2012) at Bridger Bowl Ski Area (11 km northwest of BTL); (4) monthly February–May snow-water equivalent and temperature recorded (1993–2012) at Bracket Creek Snow Telemetry station (National Resources Conservation Service, United States Department of Agriculture; 13 km northwest of BTL); and (5) the regional (Bozeman area) Palmer Drought Severity Index (PDSI) for August from 1969 to 2012 (National Atmospheric and Oceanic Administration). The PDSI is based on a supply–demand algorithm of soil moisture that expresses dryness on the basis of recent precipitation and temperature and is effective at determining droughts lasting months and longer. We summarized all weather station data by month and tested for trends in annual, monthly and seasonal (winter (November–April), spring (May–June), growing season (June–September), late-summer (August–September)) values. All the weather stations are located west of BTL and thus record the effect of predominant wind direction and storms that influence BTL. For overlapping years, we evaluated the regional coherence of any secular trends via correlation between data collected at BZN (1,360 MASL) and MSU (1,500 MASL) stations and compared these with the data collected at higher elevation and closer proximity stations (Bridger Bowl: 2,026 MASL; Bracket Creek 1,890 MASL). This analysis indicates strong correspondence among weather stations (Supplementary Figs 2 and 3) suggesting that MSU and BZN records faithfully represent regionally consistent long-term climate patterns. Hence, we used the complete historical records (MSU and PDSI) for analysis ANPP–climate relationships.

Data analysis