DOI: 10.5061/dryad.5tb2rbpdx
title: Yellowstone’s free moving large bison herds provide a glimpse of their past ecosystem function

description: While momentum is building to restore bison across North America, most
 efforts focus on small, managed herds, leaving it unclear how large,
 migrating bison shape landscapes and whether their effects enhance or
 degrade ecosystems. We assessed carbon and nitrogen dynamics across the
 migratory landscape of bison in Yellowstone, one of the last large
 migratory populations. Bison stabilized net aboveground production while
 accelerating nitrogen turnover, increasing aboveground nitrogen pools
 while carbon pools remained stable, which improved landscape nutritional
 quality. Effects were strongest in wet, nutrient-rich habitats that
 received higher bison densities and grazing than recommended in rangeland
 management, while soil and plant conditions suggested landscape
 resilience. Restoration should embrace heterogeneity in densities and
 effects across habitats and spatial scales beyond those guiding most
 current recovery efforts.
description: # Yellowstone’s free moving large bison herds provide a glimpse of their
 past ecosystem function This data set includes data files and R markdown
 code to reconstruct the analyses in
 [https://www.science.org/doi/10.1126/science.adu0703](https://www.science.org/doi/10.1126/science.adu0703). The analyses evaluated carbon and nitrogen dynamics in naturally grazed and fenced treatments that prohibited bison grazing in Yellowstone National Park. We used a repeated measures design to assess carbon and nitrogen dynamics across the migratory landscape. The design included two between-subject factors: habitat type (habitat) and grazing treatment (treatment); and one within-subject factor, the month of the growing season (month). We categorized habitats into three types: dry, lawn-forming, and high elevation. Lawn-forming habitats occurred in valley bottoms. Dry habitats occurred on adjacent hillsides and high elevation habitats occurred in topographically wet higher elevation areas. We monitored eight dry habitats, four lawn-forming habitats, and four high elevation habitats. The grazed treatment reflected natural herbivory by bison and other herbivores. We defined “month” as the number of 30-day intervals since snowmelt at each site, which approximated plant phenology, to compare response variables across habitats over similar stages of seasonal plant development. Our response variables were aboveground production, consumption, aboveground plant carbon and nitrogen pools, aboveground plant nutrient percentages, soil carbon and nitrogen pools, soil moisture, soil microbial C mass, ammonia oxidizing bacteria abundance, and soil inorganic N availability. We monitored 16 sites from 2015 to 2021, with each site representing a single habitat level. Sites were approximately 1 ha each. Within each site, grazed and ungrazed treatments were paired, resulting in two samples per site. Each site-treatment combination served as a sample unit, and three to six subsamples were collected and averaged. Sites were selected using stratified sampling informed by GPS collar data to identify bison use hotspots along the migration route. Plots were randomly placed within site boundaries each year. See the main paper and supplement for detailed description of methods. 
 
 ## Description of the data and file structure The files reconstruct the
 analyses presented
 in [https://www.science.org/doi/10.1126/science.adu0703](https://www.science.org/doi/10.1126/science.adu0703). The data files are named based on the table or figure that they reproduce. Files are ordered below in their occurrence in the r Markdown file. The majority of analyses and data files are based on monitoring of 16 sites. Not all sites were not monitored for each variable each year or month.  | **Site** | **Habitat** | **NAP (g m^-2^** **yr^-1^)** | **pH** | **Bulk Density (g cm^-3^)** | **Soil C (Mg C ha^-1^)** | Sand/Silt/ Clay | **USDA Soil Texture** | | :------- | :------------- | :--------------------------- | :----- | :-------------------------- | :----------------------- | :-------------- | :-------------------- | | AGATE | high elevation | 
 
 | 
 
 | 0.95 | 65.3 |  40/15/45 | Loam | | ASPENS2 | dry | 106.4 | 6.8 | 
 
 | 
 
 | 41/39/20 | Clay Loam | | CC1 | high elevation | 254.7 | 5.91 | 0.96 |
 69.3 | 48/27/25 | Sandy Clay Loam | | CC2 | high elevation | 228.4 | 5.85
 | 1.01 | 71.7 | 48/29/23 | Sandy Clay Loam | | CRYS | dry | 65.9 | 6.41 |
 1.15 | 28.9 | 65/21/14 | Sandy Loam | | FIRETR | dry | 76.1 | 6.75 | 1.21
 | 39.6 | 59/20/21 | Sandy Clay Loam | | HRC | dry | 121.2 | 7.58 | 1.09 |
 36.4 | 55/21/24 | Sandy Clay Loam | | HRS | dry | 97.9 | 6.59 | 1.22 |
 64.4 | 51/23/26 | Sandy Clay Loam | | HRT | dry | 106.7 | 6.42 | 1.14 |
 51.8 | 41/26/33 | Loam | | LAMAR | lawn-forming | 310 | 7.66 | 1.02 | 34 |
 49/25/26 | Loam | | LVF | lawn-forming | 205.9 | 7.15 | 0.99 | 38 |
 35/24/41 | Loam | | OPAL | high elevation | 228.2 | 5.38 | 1.02 | 62.3 |
 43/32/25 | Loam | | PEREG | dry | 47.9 | 
 
 | 
 
 | 
 
 | 
 
 | 
 
 | | PLATEAU | dry | 50.4 | 6.58 | 
 
 | 
 
 | 49/25/26 | Sandy Clay Loam | | ROSE | lawn-forming | 201.3 | 7.38 | 
 
 | 
 
 | 38/24/38 | Loam | | SPEC | lawn-forming | 221.3 | 6.05 | 0.92 | 45.7 |
 62/19/19 | Sandy Loam | 
 
 ### Data_Fig2C.csv *Description:* Net aboveground production g m^-2^
 yr^-1^ measured in sites across habitats for grazed and ungrazed
 treatments *Table Headings:* 1. id: unique identifier for repeated
 measures ANOVA  2. site: study site 3. year: growing season year data were
 collected 4. treatment: treatment type - grazed or ungrazed (fenced
 exclosure) 5. nap: net aboveground production g m^-2^ yr^-1^ 6. habitat:
 habitat type - dry, lawn-forming or high elevation ### Data_Fig2D.csv
 *Description:* Plant aboveground carbon: nitrogen ratio measured in sites
 within lawn-forming habitats *Table Headings:* 1. site: study site - note
 there are multiple subsamples for each site and year 2. year: growing
 season year data were collected 3. treatment: treatment type - grazed or
 ungrazed (fenced exclosure) 4. julian: julian day plant sample was
 collected^ ^ 5. ratio: carbon:nitrogen ratio of plant tissue ###
 Data_Fig4A.csv *Description:* Soil carbon Mg Carbon Ha^-1^ measured in
 sites during the 1990s and 2010s *Table Headings:* 1. location: study
 site  2. habitat: habitat type - dry or wet, wet included high elevation
 and lawn-forming habitats in other analyses 3. soil class: USAD soil
 classification 4. MgC_Ha: soil carbon Mg Carbon Ha^-1^  5. period:
 historic for data collected during 1990s and present for data collected
 during 2010s  ### Data_Fig4B.csv *Description:* Net aboveground production
 g m^-2^ yr^-1^ measured in sites during the 1980-90s and 2010s *Table
 Headings:* 1. site: study site  2. year: year of study 3. nap: net
 aboveground production g m^-2^ yr^-1^ 4. habitat: habitat type - dry or
 wet, wet included high elevation and lawn-forming habitats in other
 analyses 5. period: historic for data collected during 1980s-90s and
 present for data collected during 2010s and later ### Data_FigS1.csv
 *Description:* Relative abundance of elk, bison, pronghorn and deer using
 study sites measured by remote cameras *Table Headings:* 1. habitat:
 habitat type - dry, lawn-forming or high elevation 2. species: indicator
 for elk, bison, pronghorn and deer 3. value: relative use identified by
 numbers of animals captured by photography ### Data_TableS3.csv
 *Description:* Summary by site across years of net consumption g m^-2^
 season^-1^ and bison density estimated from gps data *Table Headings:* 1.
 habitat: habitat type - dry, lawn-forming or high elevation 2.
 netcons: net consumption g m^-2^ season^-1^  3. density: bison density
 bison km^-2^ season^-1^ ### Data_TableS4.csv *Description:* Net
 consumption g m^-2^ season^-1^ measured in sites across habitats *Table
 Headings:* 1. site: study site  2. year: growing season year data were
 collected 3. habitat: habitat type - dry, lawn-forming or high elevation
 4. nc: Net consumption g m^-2^ season^-1^ ### Data_TableS5.csv
 *Description:* Daily consumption rate g m^-2^ d^-1^ measured in sites
 across habitats for 30-d intervals of the growing season *Table Headings:*
 1. id: unique identifier for repeated measures ANOVA  2. month: we defined
 “month” as the number of 30-day intervals since snowmelt at each site,
 which approximated plant phenology 3. habitat: habitat type - dry,
 lawn-forming or high elevation 4. cr: daily consumption rate g m^-2^ d^-1^
 ### Data_TableS6.csv *Description:* Net consumption g m^-2^ season^-1^ and
 net aboveground production g m^-2^ yr^-1^ measured in sites across
 habitats for grazed and ungrazed treatments *Table Headings:* 1.
 treatment: The grazed treatment reflected natural herbivory by bison and
 other herbivores and ungrazed reflected no grazing by large herbivores 2.
 site: study site 3. year: growing season year data were collected 4. nap:
 net aboveground production g m^-2^ yr^-1^ 5. nc: Net consumption g m^-2^
 season^-1^ 6. habitat: 1(=dry), 2 (=lawn-forming), 3 (=high elevation)
 indicator variable used in plotting 7. habitatdesc: habitat - dry,
 lawn-forming or high elevation ### Data_TableS7.csv *Description:* Plant
 aboveground carbon pool size g m^-2^ measured in sites across habitats and
 in grazed and ungrazed treatments during 30-d intervals of the growing
 season *Table Headings:* 1. id: unique identifier for repeated measures
 ANOVA  2. habitat: habitat type - dry, lawn-forming or high elevation 3.
 month: we defined “month” as the number of 30-day intervals since snowmelt
 at each site, which approximated plant phenology 4. treatment: The grazed
 treatment reflected natural herbivory by bison and other herbivores and
 ungrazed reflected no grazing by large herbivores 5. meanC_gm2: Plant
 aboveground carbon pool size g m^-2^ ### Data_TableS8.csv *Description:*
 Plant aboveground nitrogen pool size g m^-2^ measured in sites across
 habitats and in grazed and ungrazed treatments during 30-d intervals of
 the growing season *Table Headings:* 1. id: unique identifier for repeated
 measures ANOVA  2. habitat: habitat type - dry, lawn-forming or high
 elevation 3. month: we defined “month” as the number of 30-day intervals
 since snowmelt at each site, which approximated plant phenology 4.
 treatment: The grazed treatment reflected natural herbivory by bison and
 other herbivores and ungrazed reflected no grazing by large herbivores 5.
 meanN_gm2: Plant aboveground nitrogen pool size g m^-2^ ###
 Data_TableS9.csv *Description:* Daily aboveground growth rate g m^-2^
 d^-1^ measured in sites across habitats and in grazed and ungrazed
 treatments for 30-d intervals of the growing season *Table Headings:* 1.
 id: unique identifier for repeated measures ANOVA  2. month: we defined
 “month” as the number of 30-day intervals since snowmelt at each site,
 which approximated plant phenology 3. habitat: habitat type - dry,
 lawn-forming or high elevation 4. gr: daily growth rate g m^-2^ d^-1^ 5.
 treatment: The grazed treatment reflected natural herbivory by bison and
 other herbivores and ungrazed reflected no grazing by large herbivores ###
 Data_TableS10.csv *Description:* Soil percent moisture measured in sites
 across habitats and in grazed and ungrazed treatments during 30-d
 intervals of the growing season *Table Headings:* 1. id: unique identifier
 for repeated measures ANOVA  2. habitat: habitat type - dry, lawn-forming
 or high elevation 3. month: we defined “month” as the number of 30-day
 intervals since snowmelt at each site, which approximated plant phenology
 4. treatment: The grazed treatment reflected natural herbivory by bison
 and other herbivores and ungrazed reflected no grazing by large herbivores
 5. meansoilM: percent soil moisture ### Data_TableS11.csv *Description:*
 Aboveground plant carbon: nitrogen ratio measured in sites across habitats
 and in grazed and ungrazed treatments during 30-d intervals of the growing
 season *Table Headings:* 1. id: unique identifier for repeated measures
 ANOVA 2. month: we defined “month” as the number of 30-day intervals since
 snowmelt at each site, which approximated plant phenology 3. habitat:
 habitat type - dry, lawn-forming or high elevation 4. meanCN: carbon:
 nitrogen ratio 5. treatment: The grazed treatment reflected natural
 herbivory by bison and other herbivores and ungrazed reflected no grazing
 by large herbivores ### Data_TableS12.csv *Description:* 30 d microbial
 mass measured in sites across habitats and in grazed and ungrazed
 treatments near peak of the growing season *Table Headings:* 1. id: study
 site 2. meanmicmass: 30 d microbial mass mg microbial C kg^-1^ soil 3.
 treatment: The grazed treatment reflected natural herbivory by bison and
 other herbivores and ungrazed reflected no grazing by large herbivores 4.
 habitat: habitat type - dry, lawn-forming or high elevation ###
 Data_TableS13.csv *Description:* Abundance of ammonia oxidizing bacteria
 measured in sites across habitats and in grazed and ungrazed treatments
 near peak of the growing season *Table Headings:* 1. id: study site  2.
 month: we defined “month” as the number of 30-day intervals since snowmelt
 at each site, which approximated plant phenology 3. habitat: habitat type
 - dry, lawn-forming or high elevation 4. treatment: The grazed treatment
 reflected natural herbivory by bison and other herbivores and ungrazed
 reflected no grazing by large herbivores 5. AOB: Abundance was measured as
 the number of gene copies per gram of soil on the log scale of the ammonia
 monooxygenase (amoA) gene in Eubacteria quantified by qPCR ###
 Data_TableS14_15.csv *Description:* Mineralized N flux measured by ion
 exchange membranes measured during the 30 days preceding peak of the
 growing season and soil extractable N at the beginning and end of the
 30-day interval *Table Headings:* 1. id: study site  2. month: we defined
 “month” as the number of 30-day intervals since snowmelt at each site,
 which approximated plant phenology 3. habitat: habitat type - dry,
 lawn-forming or high elevation 4. treatment: The grazed treatment
 reflected natural herbivory by bison and other herbivores and ungrazed
 reflected no grazing by large herbivores 5. meanNcm2d: N flux mg N cm^-2^
 d^-1^ 6. meanextrN_strt: soil extractable N at the start of the 30 day
 period mg kg^-1^ soil 7. meanextrN_end: soil extractable N at the end of
 the 30 day period mg kg^-1^ soil ### Data_TableS16.csv *Description:*
 Total soil N pool size measured in sites across habitats and in grazed and
 ungrazed treatments during 30-d intervals of the growing season *Table
 Headings:* 1. id: unique identifier for repeated measures ANOVA  2.
 habitat: habitat type - dry, lawn-forming or high elevation 3. month: we
 defined “month” as the number of 30-day intervals since snowmelt at each
 site, which approximated plant phenology 4. treatment: The grazed
 treatment reflected natural herbivory by bison and other herbivores and
 ungrazed reflected no grazing by large herbivores 5. meanMgn_ha: soil
 total N (organic and inorganic) pool size Mg ha^-1^ ### Data_TableS17.csv
 *Description:* Total soil C pool size measured in sites across habitats
 and in grazed and ungrazed treatments during 30-d intervals of the growing
 season *Table Headings:* 1. id: unique identifier for repeated measures
 ANOVA  2. month: we defined “month” as the number of 30-day intervals
 since snowmelt at each site, which approximated plant phenology 3.
 habitat: habitat type - dry, lawn-forming or high elevation 4. treatment:
 The grazed treatment reflected natural herbivory by bison and other
 herbivores and ungrazed reflected no grazing by large herbivores 5.
 meanMgc_ha: soil total C (organic and inorganic) pool size Mg ha^-1^ ###
 Data_TableS18.csv *Description:* Plant aboveground N pool size when
 reintegrating consumed N measured in sites across habitats and in grazed
 and ungrazed treatments during 30-d intervals of the growing season *Table
 Headings:* 1. id: unique identifier for repeated measures ANOVA  2. year:
 year of study 3. site: study site 4. habitat: habitat type - dry,
 lawn-forming or high elevation 5. month: we defined “month” as the number
 of 30-day intervals since snowmelt at each site, which approximated plant
 phenology 6. treatment: The grazed treatment reflected natural herbivory
 by bison and other herbivores and ungrazed reflected no grazing by large
 herbivores 7. totalN_gm2: Aboveground N pool size g m^-2^when
 reintegrating consumed N ### Data_TableS19.csv *Description:* Plant
 aboveground C pool size when reintegrating consumed N measured in sites
 across habitats and in grazed and ungrazed treatments during 30-d
 intervals of the growing season *Table Headings:* 1. id: unique identifier
 for repeated measures ANOVA  2. habitat: habitat type - dry, lawn-forming
 or high elevation 3. month: we defined “month” as the number of 30-day
 intervals since snowmelt at each site, which approximated plant phenology
 4. treatment: The grazed treatment reflected natural herbivory by bison
 and other herbivores and ungrazed reflected no grazing by large herbivores
 5. totalC_gm2: Aboveground N pool size g m^-2^when reintegrating consumed
 N ### Data_TableS20.csv *Description:* Plant aboveground N pool size when
 reintegrating consumed N measured at the peak of the growing season and
 comparisons to LandSat derived INDVI *Table Headings:* 1. id: unique
 identifier  2. habitat: habitat type - dry, lawn-forming or high elevation
 3. treatment: The grazed treatment reflected natural herbivory by bison
 and other herbivores and ungrazed reflected no grazing by large herbivores
 4. totalN_gm2: Aboveground N pool size g m^-2 ^when reintegrating consumed
 N 5. indvi: LandSat derived integrated INDVI. Refer to supplemental
 material of main paper for curve fitting  ### Data_TableS21.csv
 *Description:* Plant community diversity measures within sites for
 habitats measured during the 1960s and 2010s. *Table Headings:* 1. site:
 study site  2. habitat: habitat type - wet or dry 3. period: present
 represents 2010s and historic represents 1960s 4. simpson: simpson's
 diversity index 5. shannon: shannon's diversity index 6. richness:
 richness 7. pileou_evenness: Pileou's evenness  8.
 distance_to_centroid: distance to centroid ### Data_INDVI.tif
 *Description:* Raster image of integrated INDVI estimated using the
 methods described in 
 [https://www.science.org/doi/10.1126/science.adu0703](https://www.science.org/doi/10.1126/science.adu0703). Raster used in R markdown file to relate INDVI to net aboveground production g m^-2^ yr^-1^  ### Data_UD.tif *Description:* Raster image of utilization distribution for bison estimated using the methods described in  [https://www.science.org/doi/10.1126/science.adu0703](https://www.science.org/doi/10.1126/science.adu0703). Raster used in R markdown file to describe bison density bison km^-2^ ### Data_seasonalUD.tif *Description:* Raster image of utilization distribution for bison estimated using the methods described in  [https://www.science.org/doi/10.1126/science.adu0703](https://www.science.org/doi/10.1126/science.adu0703). Raster used in R markdown file to describe bison density bison km^-2^ during the growing season only. ### Data_supp001.csv *Description:* Site measured net aboveground production g m^-2^ yr^-1^ and INDVI from associated areas *Table Headings:* 1. site: study site  2. year: year of study 3. NAP: Site measured net aboveground production g m^-2^ yr^-1^ 4. INDVI: Integrated NDVI measured using the methods described in  [https://www.science.org/doi/10.1126/science.adu0703](https://www.science.org/doi/10.1126/science.adu0703). ### R Code.Rmd The file R Code.Rmd is an R markdown file that calls the associated data files to reconstruct the tables and figures in [https://www.science.org/doi/10.1126/science.adu0703](https://www.science.org/doi/10.1126/science.adu0703).