Data from: Increases in understory plant cover and richness persist following restoration treatments in Pinus ponderosa forests
Data files
Oct 25, 2023 version files 2.43 GB
-
README.md
25.69 KB
-
Springer_LEARN_archive.7z
2.43 GB
Abstract
A combination of forest thinning followed by prescribed burning is widely applied in the western US to increase ecosystem resistance and resilience to disturbances. Understory plant community responses may be driven both by management treatments and climatic factors. Thus, responses to treatments during a 20-year megadrought have implications for the role of management in fostering adaptive capacity to climate change.
We used a network of five sites (600 plots) spanning an environmental gradient in ponderosa pine (Pinus ponderosa) forests of the American Southwest, an ecosystem that is broadly distributed and actively managed throughout the western US. We used repeated long-term monitoring data to quantify plant community responses to treatment 1-5, 6-10, and >10 years post-implementation. Specifically, we focused on the effects of treatment and abiotic conditions on native and nonnative plant cover and species richness, and on the proportion of native species with northern (cool-mesic) biogeographic affinities.
Overall, thinning and prescribed burning nearly doubled native cover and increased native species richness by about 50% relative to untreated controls. These effects persisted for over a decade after treatment, even under the influence of significant and persistent drought. Cover and richness were also greater on intermediate to wet sites. Finally, native species with northern biogeographic affinities were reduced for up to five years after treatment relative to those with southern (warm-xeric) affinities, and in dry years, indicating that both management and interannual climate variability may foster shifts in plant communities that are more resilient to a warming climate.
Synthesis and applications: In ponderosa pine forests of the American Southwest, tree thinning followed by prescribed burning will generally promote restoration goals of increasing resilience to climate change by enhancing the diversity and abundance of native understory plant species, even during a persistent 20-year megadrought.
Authors
Judith D. Springer, Michael T. Stoddard, Kyle C. Rodman, David W. Huffman, Paula J. Fornwalt, Rory J. Pedersen, Daniel C. Laughlin, Christopher M. McGlone, Mark L. Daniels, Peter Z. Fulé, Margaret M. Moore, Becky K. Kerns, Jens T. Stevens, Julie E. Korb, Sara Souther
Overview
This repository includes field data, R code, and model outputs from “Increases in understory plant cover and richness persist following restoration treatments in Pinus ponderosa forests” by Springer et al. (In Review), J Appl Ecol. Data were archived in .7z format using 7Zip v 23.01 on Windows 11 v 22H2
File Structure
This repository is split into the following subfolders which contain different portions of the project. Descriptions of each file contained in these folders are given within each subsection.
AnalysisOutputs
This subfolder contains generalized linear mixed model (GLMM) objects created in the “Step3_RichnessCoverGLMMs.R” script. Models were fit using the spaMM and glmmTMB packages in R and can be read into R using the “readRDS” function. spaMM and glmmTMB packages should be loaded for easier loading and viewing of these model objects.
- coverGLMM.rds: the fitted model object for the binomial GLMM of understory plant cover
- nAffGLMM.rds: the fitted model object for the binomial GLMM of proportion northern affinity species
- richnessGLMM.rds: the fitted model object for the Poisson GLMM of understory species richness
Code
This subfolder contains R scripts used to process, analyze, or visualize data for the project. Raw data were filtered and/or queried using SQL queries in Microsoft Access, which are not provided here.
- Step1a_HLI_Calcs.R: Calculates a 30-m raster version of heat load index (HLI) using Eq. 3 of McCune and Keon (2002). “Data/SpatialData/SoilsAndTerrain/hliUTM12N.tif”, described below, is the output file. HLI was used in water balance model calculations
- Step1b_TPI_Calcs.R: Calculates a 30-m raster version of Topographic Position Index (TPI) following Weiss (2001). “Data/SpatialData/SoilsAndTerrain/tpi15x15.tif”, described below, is the output file
- Step1c_GIDS-Downscaling.R: Downscales gridded climate data from PRISM to a 30-m resolution using gradient inverse and distance-squared weighted interpolation of Nalder and Wein (1998). Digital elevation models at 800, 270, 90, and 30m were used s ancilary data
- Step1d_CWD-Calcs.R: Uses downscaled climate data from “Step1c…” and heat load index from “Step1a…”, along with soil data, to run monthly water balance models following Lutz et al. (2010). Output datasets are annual climatic water deficit (CWD) and actual evapotranspiration (AET) for each site in “Data/SpatialData/WaterBalance”
- Step2_ExtractSpatial.R: Organizes summarized field data into an analysis-ready format. TPI and CWD are extracted at the location of each field plot and added as additional columns in “Data/plantDataWithSpatial.csv”, which is produced in this script.
- Step3_RichnessCoverGLMMs.R: This script pulls in the “Data/plantDataWithSpatial.csv” file and fits GLMMs for understory cover, species richness, and proportion of northern affinity species. Also makes plots to visualize the results of these models
- Step4_ClimateSeries.R: This script pulls in monthly climate data at each site, and formats it to plot annual precipitation values at each site over time, as part of the study area figure in the paper
Data
This subfolder contains raw and summarized datasets used in this study. Data include field surveys of plant understory communities at each of five experimental sites in Arizona, as well as topographic data, soils data, downscaled climate data, and water balance model outputs at each site. The summarized dataset “plantDataWithSpatial.csv” was used to develop all GLMMs, and is located within the main “Data” folder, rather than subfolders. Each file is described below
- plantDataWithSpatial.csv: This file contains all summarized plot-level data and model covariates included in GLMMs of understory cover, species richness, and proportional dominance of northern affinity species. Columns have a range of data types. It was created using the “Step2_ExtractSpatial.R” script described above. Fields differ in data type and the range of values. Rows represent unique combinations of plots and survey years.
- Site: Two-letter code for each study site. CF = Centennial Forest, FV = Fort Valley, GV = Grandview, MC = Mineral Creek, MT = Mt. Trumbull
- Year: Calendar year for a given survey
- Block: Experimental block ID. Replicated block within a site. Integer from 1 to 4
- Trt: Treatment unit for a given plot. Values of 1 or 10 are “control units” which never received a treatment. Values of 2 are treated units which received thinning and burning treatments at some point during the study. Thus, initial surveys at treated units were still “untreated”
- Plot: Unique plot number. When combined with Site, Block, and Trt, this identifies a unique plot ID in the dataset.
- Long: Longitude (in decimal degrees) of the center of a plot/understory transect. Coordinates are in geographic coordinate system EPSG:4326 (WGS 84; Lat/Long)
- Lat: Latitude (in decimal degrees) of the center of a plot/understory transect. Coordinates are in geographic coordinate system EPSG:4326 (WGS 84; Lat/Long)
- UTMx: Easting (in meters) of the center of a plot/understory transect. Coordinates are in projected coordinate system EPSG:26912 (NAD83 UTM zone 12N).
- UTMy: Northing (in meters) of the center of a plot/understory transect. Coordinates are in projected coordinate system EPSG:26912 (NAD83 UTM zone 12N).
- Year_from_Initial_treatment: For plots that received treatment, the number of years since initial thinning was implemented at a site. Measurements that were prior to initial treatment, or in untreated units have “NA” values. Integer values from 4 to 20
- freq_Tot: Total understory plant cover in a given plot and survey year, as the number of plant “hits” across the 166 line transect points/166 * 100. Because a small number of points had multiple plant “hits”, where separate plants intersected the transect at a given point, the total value of this variable could occasionally exceed 100 for a plot (n = 4 plot measurements). In these cases, these values were truncated to [0,100] for subsequent analyses to allow use of binomial GLMMs, though untruncated values are provided in this archive. Used as a response variable in the GLMM of total plant cover. Floating point values from 0.0 to 131.9277
- freq_I: Understory plant cover of introduced species in a given plot and survey year, as a percentage of the 166 line transect points that had introduced vegetation. See “Data/Raw/speciesList.csv” for descriptions of nativity status for individual species. Floating point values from 0.0 to 65.66265
- freq_N: Understory plant cover of native species in a given plot and survey year, as a percentage of the 166 line transect points that had native vegetation. See “Data/Raw/speciesList.csv” for descriptions of nativity status for individual species. Floating point values from 0.0 to 66.26506
- rich_Tot: Total understory plant richness in a given plot and survey year, as the number of unique plant species on belt transects. Used as a response variable in the GLMM of total plant richness. Integer values from 1 to 74
- rich_I: Richness of introduced species in a given plot and survey year, as the number of unique introduced plant species on belt transects. See “Data/Raw/speciesList.csv” for descriptions of nativity status for individual species. Integer values from 1 to 10
- rich_N: Richness of native species in a given plot and survey year, as the number of unique native plant species on belt transects. See “Data/Raw/speciesList.csv” for descriptions of nativity status for individual species. Integer values from 1 to 68
- propNorthern: Proportion of points with vegetation, where a norther affinity species represented the top plant strata. Used as a response variable in the GLMM of northern affinity proportional dominance. Measurements had no vegetation cover or had no species that could be classified with a biogeographic affinity have “NA” values. Floating point values from 0.0 to 1.0
- vegCount: The number of points with vegetation in line transects. Used as a “n trials” in binomial model of proportional northern affinity species dominance. Integer values from 0 to 147
- tpi15x15: Topographic position index (TPI) of Weiss (2001)*10. Extracted at the center of each plot from “tpi15x15.tif”, which is further described below. Integer values from -11 to 21
- avgCWD_1991_2020: Average (1991-2020) annual climatic water deficit (CWD) at a given plot in mm, calculated using the monthly water balance model of Lutz et al. (2010). Extracted at the center of each plot from “[Site]_CWD_1981-2021.tif” files, which are further described below. Annual values from 1991 to 2020 were extracted separately, and then averaged using the arithmetic mean. Floating point values from 74.53267 to 348.3973
- annCWD_z: Z-score transformed devations of a sampling year from the longer-term a (1991-2020) annual climatic water deficit (CWD) at a given plot in mm, calculated using the monthly water balance model of Lutz et al. (2010). Extracted at the center of each plot from “[Site]_CWD_1981-2021.tif” files, which are further described below. Annual values from 1991 to 2020 were extracted separately, and then averaged using the arithmetic mean.
- treatmentCode: A simplified categorical code representing both treatment status and time since initial thinning treatments for a given plot and year. Character/factor variables with four unique values. Untreated, Treated (1-5 yr), Treated (6-10 yr), Treated (>10 yr). Used as a model covariate in GLMMs
FieldData
This folder includes a single .xlsx file (“LEARN_TransQuadSummaries.xlsx”) with field data obtained from MS Access queries and summaries of raw field data (shown in “Data/Raw”). Descriptions of individual fields and rows are identical to those included in “Data/plantDataWithSpatial.csv” below. This .xlsx file is only included in the archive for reproducibility of the “Step2_ExtractSpatial.R” script described above
Raw
This folder includes four raw data files, which describe understory plant community data collected at each of the five research sites in Arizona, USA. These are likely to be the files that are of the greatest interest to other users. Three of the files are specific to plant communities, while one file gives the years in which treatments were implemented at each site and experimental block
- speciesList.csv: This file summarizes information about each plant species found at our research sites including their taxonomy (as of fall 2023), nativity status in Arizona, biogeographic affinity, common name, and six-letter species code (which helps to link this table to other raw data files). Each row represents a single plant species
- SpeciesCode: Shortened species code (i.e., first three letters of genus and first three level of species). Used to uniquely identify each species in other datasets
- Family: Family-level taxonomy for a given species
- Genus: Genus-level taxonomy for a given species
- Species: Species-level taxonomy for a given species
- NativityStatus: N (native) or I (introduced). Whether a species is native to the study sites in Arizona, USA
- BiogeoAffinity: Biogeographic affinity of each species (as described in supplemental information for the paper). Values are:
- CFP: California Floristic Province
- MAm: Meso-American (warm/xeric affinity)
- MaT: Madro-Tertiary (warm/xeric affinity)
- NA: Non-native species unable to be assigned (species removed from analysis of biogeographic affinities)
- NTm: North-Temperate (cool/mesic affinity)
- SAm: South American (warm/xeric affinity)
- Unk: Unknown affinity, no information available (species removed from analysis of biogeographic affinities)
- WTD: Warm-Temperate/Desert (warm/xeric affinity)
- CommonName: English common name for a given species
- beltTransects.csv: Data used to describe species richness at each plot and survey year. Rows represent each combination of plot, survey year, and plant species. Only species that were “present” in a given plot and survey year are shown in this file. Absences are omitted to reduce file sizes, but species-level absences can be obtained by joining this file with the “speciesList.csv” file. Any species included in “speciesList.csv” but without a row in a given plot/year combination was not observed in that location and survey.
- Site: Two-letter code for each study site. CF = Centennial Forest, FV = Fort Valley, GV = Grandview, MC = Mineral Creek, MT = Mt. Trumbull
- Year: Calendar year for a given survey
- Block: Experimental block ID. Replicated block within a site. Integer from 1 to 4
- Trt: Treatment unit for a given plot. Values of 1 or 10 are “control units” which never received a treatment. Values of 2 are treated units which received thinning and burning treatments at some point during the study. Thus, initial surveys at treated units were still “untreated”
- Plot: Unique plot number. When combined with Site, Block, and Trt, this identifies a unique plot ID in the dataset.
- SpeciesCode: Six-letter species code for which presence is being assessed in a given plot and year. Codes correspond to those in “speciesList.csv”
- lineTransects.csv: Point-line-intercept transects used to describe plant cover and the dominance of northern affinity species. Rows represent individual points (166 per plot and year) along a transect, in a given plot and year. Only points with vegetation are included in this file (to reduce file sizes). Therefore, “missing” points within a plot and year had no cover of any understory plant species. Replicating these missing values (such that each combination of site, block, Trt, Plot, and Year have points ranging from 1-166 and a SpeciesCode of “NONE”) may be needed for some intended uses.
- Site: Two-letter code for each study site. CF = Centennial Forest, FV = Fort Valley, GV = Grandview, MC = Mineral Creek, MT = Mt. Trumbull
- Year: Calendar year for a given survey
- Block: Experimental block ID. Replicated block within a site. Integer from 1 to 4
- Trt: Treatment unit for a given plot. Values of 1 or 10 are “control units” which never received a treatment. Values of 2 are treated units which received thinning and burning treatments at some point during the study. Thus, initial surveys at treated units were still “untreated”
- Plot: Unique plot number. When combined with Site, Block, and Trt, this identifies a unique plot ID in the dataset.
- SpeciesCode: Six-letter species code for the top plant strata at a given point “speciesList.csv”. Multiple species could be identified at each point, however, only the top strata is given here. Belt transects provide better information about the individual species present at a site.
- treatmentYears.csv: File describing the years in which different treatments where implemented at each site and treatment-control unit pair (i.e., block). NOTE: Years in this file correspond to the first growing season in which a block might be considered post-treatment. For treatments in the spring, this was the same year as treatment implementation. for treatments in the summer/fall, this was the year after treatment implementation.
- Site: Two-letter code for each study site. CF = Centennial Forest, FV = Fort Valley, GV = Grandview, MC = Mineral Creek, MT = Mt. Trumbull
- Block: Experimental block ID. Replicated block within a site. Integer from 1 to 4
- ThinningYear: The first growing season in which the treated unit of a block was considered “thinned”. Integer from 1999 to 2007.
- BurningYear1: The first growing season in which the treated unit of a block was considered “burned”. Integer from 2000 to 2012.
- BurningYear2: The first growing season in which the treated unit of a block was considered “burned” more than once. Integer from 2007 to 2008. “NA” indicates that a plot was only burned once
- BurningYear3: The first growing season in which the treated unit of a block was considered “burned” more than twice. Integer from 2015 to 2015. “NA” indicates that a plot was burned one or two times only
- BurningYear4: The first growing season in which the treated unit of a block was considered “burned” more than three times. Integer from 2020 to 2020. “NA” indicates that a plot was burned one, two, or three times
SpatialData
This subfolder includes all spatial datasets used to predict understory plant cover, richness, and dominance of northern affinity species in GLMMs. Files are further divided into subfolders below, which give specific kinds of spatial data used in this study
SpatialData/PRISM_coarse
Monthly climate data (4-km spatial resolution) from the gridded PRISM dataset (PRISM Climate Group, 2022), cropped to the extent of our study sites. Each file is a multi-band raster file (GeoTIFF format) with 489 individual bands. Bands correspond to months in the period of January 1981 (Band 1) to September 2021 (Band 489), and are arranged chronologically (e.g., February 1981 is Band 2). Files are in 32-bit floating point format, to preserve decimals in data, in projected coordinate system EPSG:26912 (NAD 83, UTM zone 12 N).
- PRISM_ppt.tif: Total monthly precipitation (in mm) for the corresponding month and year.
- PRISM_tmax.tif: Average maximum temperature (in C) for the corresponding month and year.
- PRISM_tmin.tif: Average minimum temperature (in C) for the corresponding month and year.
SpatialData/PRISM_downscaled
Files in this subfolder are not included in this archive due to large file sizes. However, these files can be re-created using the “Step1c_GIDS-Downscaling.R” script, which downscales files in the “SpatialData/PRISM_coarse” folder using digital elevation models in “SpatialData/SoilsAndTerrain”. Subfolder is included here for reproducibility of this script. Output of this downscaling script gives statistically downscaled (i.e., 30-m spatial resolution; monthly temporal resoluion) climate data, cropped to the extent of each individual site.
SpatialData/SiteLocations
Files in this subfolder give the locations of all plots included in our analyses. There are four files, all of which are different components of a vector-based GIS point dataset in shapefile (.shp) format. Points give the center location of each understory vegetation transect. This file is in projected coordinate system EPSG:26912 (NAD 83, UTM zone 12 N), with the following fields in the attribute table (i.e., .dbf file).
- Site: Two-letter code for each study site. CF = Centennial Forest, FV = Fort Valley, GV = Grandview, MC = Mineral Creek, MT = Mt. Trumbull
- Plot_ID: A concatenated version of the block, treatment, and plot ID of a corresponding point. The format follows the convention [Block]-[Trt]-[Plot]. This field and the site field can be used to join this dataset with raw plant data in “beltTransects.csv”, “lineTransects.csv”, etc.
- Long: Longitude (in decimal degrees) of the center of a plot/understory transect. Coordinates are in geographic coordinate system EPSG:4326 (WGS 84; Lat/Long)
- Lat: Latitude (in decimal degrees) of the center of a plot/understory transect. Coordinates are in geographic coordinate system EPSG:4326 (WGS 84; Lat/Long)
SpatialData/SoilsAndTerrain
This subfolder includes various biophysical variables included as covariates in GLMMs (e.g., tpi15x15.tif), used to downscale gridded climate data (e.g., demUTM12N.tif), or used to run water balance models (e.g., hliUTM12N.tif, awcUTM12N.tif). All files are raster datasets (GeoTIFF format) in projected coordinate system EPSG:26912 (NAD83 UTM zone 12N). Most files have a 30-m spatial resolution and are aligned to the same reference grid, unless otherwise mentioned below.
- aspectUTM12N.tif: Aspect of a hillslope, with units giving the azimuth at which a slope faces. 0/360 = due North and 180 = due South. Units are in 32-bit floating point format to preserve decimal values. Derived from a 30-m digital elevation model (i.e., demUTM12N.tif)
- awcUTM12N.tif: Water holding capacity (in mm) of the top 200 cm of the soil profile, based on data from POLARIS (Cheney et al. 2016). A combination of soil texture and depth to bedrock characteristics. Units are in 32-bit floating point format to preserve decimal values. Used for water balance modeling
- demUTM12N.tif: Digital elevation model from the USGS 3DEP dataset (USGS, 2021). Originally acquired at 10-m spatial resolution and aggregated to 30-m resolution using “average” resampling in GDAL. Used for statistical downscaling of climate data and to develop aspect and slope grids. 16-bit integer format with pixel values in meters.
- demUTM12N_[resolution].tif: Digital elevation models from the USGS 3DEP dataset (USGS, 2021). Originally acquired at 10-m spatial resolution and aggregated to [resolution] (i.e., 90m, 270m, 900m, or 4000m) using “average” resampling in GDAL. Used for statistical downscaling of climate data. 32-bit floating point format with pixel values in meters.
- hliUTM12N.tif: Heat load index following Eq. 3 of McCune and Keon (2002). Higher values are common in steep, southwest-facing slopes, and low values are common on steep, northeast-facing slopes. Calculated from slope, aspect, and latitude grids in the same subfolder. 32-bit floating point format
- latitudeUTM12N.tif: Latitude (in decimal degrees) of the center of a 30-m pixel. Used to calculate “hliUTM12N.tif” file. 32-bit floating point format
- slopeUTM12N.tif: Slope angle (in degrees) of each 30-m pixel. Used to calculate “hliUTM12N.tif” file. 32-bit floating point format
- tpi15x15.tif: Topographic position index (TPI) of Weiss (2001), which gives the difference between the elevation of a 30-m pixel and the mean elevation of a 450-m square window centered on the 30-m pixel. Units are in meters*10 to convert to integer values and save a decimal of precision. Low, negative values are indicative of valley bottoms, whereas high, positive values are indicative of ridgetops. Used as a covariate in GLMMs. 16-bit signed integer data format
SpatialData/WaterBalance
This subfolder includes outputs of monthly water balance models, calculated in the “Step1d_CWD-Calcs.R” script. Outputs are multi-band raster files (GeoTIFF format), with 41 bands each, a 30-m spatial resolution, and in projected coordinate system EPSG:26912 (NAD83 UTM zone 12N). Bands give the annual total (summed by water year; preceding October 1st to September 30th in focal year) of each metric described below, cropped to the extent of a given site. In each file, band 1 gives the water year total for 1981 (i.e., Oct 1980-Sept 1981), band 2 gives the water year total for 1982 and band 41 gives the water year total for 2021. Most files have a 30-m spatial resolution and are aligned to the same reference grid
- [Site]_AET_1981-2021: Within a [Site] (site codes are identical to those described in datasets above), this dataset gives actual evapotranspiration (AET) in mm/yr, for each water year from 1981 to 2021. Calculated using monthly water balance models of Lutz et al. (2010) and a modified Thornthwaite approach. Not included as a covariate in GLMMs because of high correlation with CWD, and weaker explanatory power. 32-bit floating point format
- [Site]_CWD_1981-2021: Within a [Site] (site codes are identical to those described in datasets above), this dataset gives climatic water deficit (CWD) in mm/yr, for each water year from 1981 to 2021. Calculated using the monthly water balance models of Lutz et al. (2010) and a modified Thornthwaite approach. Used as a covariate in GLMMs. 32-bit floating point format
SpatialData/WorkingDirectory
This folder is intentionally empty, but is the location in which temporary output files are written during processing of some of the R scripts in the “Code” subfolder. Included in this archive for reproducibility of these scripts.
References
- Chaney, N. W., Wood, E. F., McBratney, A. B., Hempel, J. W., Nauman, T. W., Brungard, C. W., & Odgers, N. P. (2016). POLARIS: A 30-meter probabilistic soil series map of the contiguous United States. Geoderma, 274, 54-67.
- Lutz, J. A., Van Wagtendonk, J. W., & Franklin, J. F. (2010). Climatic water deficit, tree species ranges, and climate change in Yosemite National Park. Journal of Biogeography, 37(5), 936-950.
- McCune, B., & Keon, D. (2002). Equations for potential annual direct incident radiation and heat load. Journal of vegetation science, 13(4), 603-606.
- Nalder, I. A., & Wein, R. W. (1998). Spatial interpolation of climatic normals: test of a new method in the Canadian boreal forest. Agricultural and forest meteorology, 92(4), 211-225.
- PRISM Climate Group, Oregon State University. (2022). http://prism.oregornstate.edu.
- United States Geological Survey. (2021). USGS 3D Elevation Program Digital Elevation Model. Accessed from: https://developers.google.com/earth-engine/datasets/catalog/USGS_3DEP_10m
- Weiss, A. (2001). Topographic position and landforms analysis. ESRI user conference, San Diego, CA.
This archive includes long-term plant community data and forest monitoring data from 5 experimental sites (~ 600 plots) in Arizona, USA. These sites were established in the late 1990s to early 2000s to describe the effects of thinning and burning treatments on ecological processes in southwestern US ponderosa pine forests. Understory plant community data describe the presence and abundance of herbaceous and short-statured woody plants, with at least one survey before treatment, and multiple surveys after treatment at each site, as well as data from paired untreated control sites. Plants were identified in the field by trained botanists, or pressed and identified in the laboratory. Most individuals were identified to the species level, but occasioanlly could only be identified to genus. For descriptions of different fiels in the data arrchive, as well as methods of collection, see the "README.md" file