Mycorrhizal legacy mediates seedling success following timber harvesting in Northeastern forests
Data files
Jun 17, 2026 version files 42.07 MB
-
2020Initial_seedlingdata.csv
54.69 KB
-
August_2021_seedlingdata.csv
49.30 KB
-
August_2022_seedlingdata.csv
49.42 KB
-
belowground.data.csv
27.74 KB
-
Corinth__ITS_Processing.Rmd
7.80 KB
-
Corinth_data_analysis.Rmd
85.17 KB
-
Corinth_fungal_community.Rmd
28.84 KB
-
Corinth_resindata.csv
12.73 KB
-
Corinth_seedling_data.csv
483.04 KB
-
Corinth_Soil_Test_Results_9-24-24.csv
1.01 KB
-
Corinth.alldata.csv
8.87 KB
-
Corinth.DBH.species.csv
41.83 KB
-
Corinth.project.alldata.xlsx_-_Site_info.csv
808 B
-
corinth.taxa.2021.csv
3.80 MB
-
corinth.taxa.2022Summer.csv
1.91 MB
-
DNA_ID_fall2020_summerfall2021_AFitch.csv
4.59 KB
-
feature_table_2020.csv
4.74 MB
-
feature_table_2022.csv
4.40 MB
-
feature_table_2022Summer.csv
1.99 MB
-
feature_table_all_2021.csv
4.74 MB
-
FUNGuild_2021_trimmed
4.27 MB
-
FUNGuild_2022_trimmed.csv
3.61 MB
-
FUNGuild_2022_trimmedSummer.csv
1.71 MB
-
FUNGuild_F20_trimmed
4.37 MB
-
ph_postharvest.csv
2.03 KB
-
PLFA_NLFA_2022.csv
4.70 KB
-
PLFA_NLFA_condensed.csv
32.07 KB
-
plfa.model.csv
36.56 KB
-
README.md
19.77 KB
-
seedling_position_index.csv
40.05 KB
-
Soil.moisture.correction.csv
6.67 KB
-
taxon_table_2021.csv
3.67 MB
-
taxon_table_2022Summer.csv
1.85 MB
-
Tree.info.csv
14.10 KB
Abstract
Timber harvesting, a human-caused disturbance, can interact with other disturbance mechanisms such as invasive pests and climate change to speed forest community change and affect forest productivity and species composition. Overstory trees rely on two types of mycorrhizal fungi, arbuscular mycorrhizae (AM) or ectomycorrhizae (EcM), for establishment and growth, but we know very little about how timber harvesting affects these mycorrhizal fungi and subsequent downstream effects on seedling success. In this field experiment, we investigated how the legacy of AM and EcM forest stands affected seedling survival and growth, soil fungal communities, and soil chemistry including nutrient availability and pH following timber harvesting. We established sixteen plots, half of which were AM and EcM dominated. Eight plots were harvested and then planted with four AM- and four EcM-associated seedling species in a split plot design. The eight unharvested control plots allowed us to compare changes in fungal biomass, and soil pH and nutrients. Our results show that AM fungal biomass decreased in both AM and EcM plots only in the first fall after timber harvesting, while EcM fungal biomass continued to decrease in both legacy types over three growing seasons. The mycorrhizal legacy also affected seedling success: AM-associated seedlings had higher survival when planted in AM legacy plots, but there was no effect for EcM-associated seedlings. Patterns in foliar nitrogen (N) indicated that although AM seedlings closer to plot edges may have acquired more N through mycorrhizae relative to EcM seedlings at plot edges, they also had lower percent N. This pattern indicates a tradeoff between increased access to mycorrhizae for colonization and competition with proximity to live roots. Higher soil calcium (Ca) in AM legacy and control sites provides evidence for an alternative mechanism behind our observed AM legacy effect, where some AM tree species are sensitive to soil Ca levels. Together, these results illustrate a more complicated and potentially contrasting set of mechanisms that drive seedling success than is postulated by the “mother tree” hypothesis through common mycorrhizal networks.
Dataset DOI: 10.5061/dryad.18931zd75
Manuscript published in Ecological Applications: EAP24-0875.R2
Description of the data and file structure
Our research takes place in the Clement Woodlot, a northern mixed-deciduous forest in Corinth, Vermont (44°02'48.4" N, 72°19'03.0" W). We established 16 0.1 ha plots where half of the plots were dominated (>80% basal area) by EcM trees and half by AM trees. Four 0.1 ha plots per mycorrhizal type (referred to as “EcM legacy” and “AM legacy” plots) were harvested in the winter of 2020-21 and four were left as unharvested controls, for a total of eight 0.1 ha plots in each control and treatment group. To assess how harvesting affected soil conditions, we measured soil moisture, pH, fungal biomass, and nitrogen (N) availability before timber harvesting and after timber harvesting. To measure the relative differences in nutrient availability, we measured available soil nitrogen (NH4 and NO3) and phosphate (PO4) using cation and anion exchange resins. To assess harvest effects on soil microbial community composition, we combined an estimate of fungal biomass with internal transcribed spacer (ITS) sequencing to calculate biomass of fungal guilds within both the harvested and control plots. From May 14-17, 2021, we planted 20 seedlings each of eight different species per plot, four AM species (Acer saccharum ACSA, Acer rubrum ACRU, Nyssa sylvatica NYSY, and Prunus serotina PRSE), and four EcM species (Quercus rubra QURU, Carya cordiformis CACO, Tilia americana TIAM, and Betula lenta BELE). To further investigate how access to available nitrogen could affect seedling survival and growth, we measured ratios of δ15N in seedling leaves after the first season of growth in fall 2021.
We used R software for all data formatting and analyses. We used linear mixed effects models 1) to investigate potential changes in AM, EcM, pathotrophic, and saprotrophic fungal biomass, soil nutrients, pH, and soil moisture from before and after timber harvesting and between control and harvested plots using the BACI design, 2) to test for the interactive effects of mycorrhizal type (MT), or legacy, and seedling mycorrhizal association on seedling survival and growth, and 3) to test for differences in soil Ca between plot mycorrhizal types following findings that AM seedlings had a survival benefit in the AM legacy plots.
Files and variables
Seedling Data
seedling_position_index.csv
Index file linking seedling tags with species identity, plot location, and spatial position within plots.
Variables:
- Unit – Forest unit identifier
- Plot – Plot identifier, within unit
- Species – Seedling species code
- MYC – Seedling mycorrhizal association (AM or EcM)
- MT – Plot mycorrhizal legacy type (AM or EcM)
- Num. – Seedling number within plot
- tag – Seedling tag identifier (unique to each seedling)
- survival – Seedling survival status (1 alive, 0 dead)
- Prox – Seedling location within the plot (Center or Edge)
Corinth_seedling_data.csv
Complete seedling dataset that combines multiple files. See accompanying code.
Variables:
- Year – Measurement year
- Measurement.date – Date of measurement
- Unit – Forest unit identifier
- Plot – Plot identifier, within unit
- Species – Seedling species code (4 letter)
- Number – Seedling number (1 - 20, within plot)
- ID1-ID4: Combinations of information (unit, plot, species, number) for the purposes of downstream models and plot visualization
- Health – Visual health assessment (L live, D dead, V low vigor)
- Herbivory – Herbivory damage indicator (0 none - 3 intense herbivory)
- Condition – Seedling condition classification (0 none, 1 top dead or removed)
- HT.cm – Seedling height (cm); NA entries correspond to dead seedlings with no measurement
- RCD.mm – Root collar diameter (mm); NA entries correspond to dead seedlings with no measurement
- HT.cm.initial – Initial seedling height (cm); NA entries correspond to dead seedlings with no measurement
- RCD.mm.initial – Initial root collar diameter (mm); NA entries correspond to dead seedlings with no measurement
- Soil – Soil condition classification at seedling level (binary; 1 = Scarified, 0 = unscarified); NA entries correspond to dead seedlings with no measurement
- Slash – Slash cover presence classification at seedling level(0 none - 4 > 90% cover); NA entries correspond to dead seedlings with no measurement
- VSMxtree – volumetric soil moisture (%); NA entries for missing data
- Planting.date – Seedling planting date
- MT – Plot mycorrhizal legacy type (AM or EcM)
- Treatment – Harvest treatment (control or cut)
- SGgw – Species wood density coefficient (unitless, ratio of wood to water weight)
- Time.months – Months since planting
- Biomass.g – Estimated seedling biomass (g); NA entries correspond to dead seedlings with no measurement
- Biomass.initial.g – Initial biomass estimate (g); NA entries correspond to dead seedlings with no measurement
- RG.g.month – Relative growth rate (g month⁻¹); NA entries correspond to dead seedlings with no measurement
- survival, bi.cond – Seedling survival status (1, alive, 0 dead)
- HT.growth – Height growth since planting (cm); NA entries correspond to dead seedlings with no measurement
- MYC – Seedling mycorrhizal association (AM or EcM)
- placeholder - used in calculation code
2020Initial_seedlingdata.csv
Initial measurements of seedlings taken immediately after planting.
Variables are defined under "Corinth_seedling_data.csv"
NA values are for missing data in the soil.moisture column and for seedlings that died following the initial planting prior to first measurement.
August_2021_seedlingdata.csv
Seedling measurements collected in August 2021.
Variables are defined under "Corinth_seedling_data.csv"
Blank cells for "HT.cm" and "RCD.mm" are for seedlings that died and did not have these measurements.
August_2022_seedlingdata.csv
Seedling measurements collected in August 2022.
Variables are defined under "Corinth_seedling_data.csv"
Blank cells for "HT.cm" and "RCD.mm" are for seedlings that died and did not have these measurements.
Soil Chemistry and Environmental Data
Corinth_resindata.csv
Soil nutrient availability measured using ion exchange resin probes.
Variables:
- NH4 – Ammonium concentration (mg / L); NA entries for missing data
- NO3 – Nitrate concentration (mg / L); NA entries for missing data
- PO4 – Phosphate concentration (mg / L); NA entries for missing data
- Unit – Forest unit identifier
- Plot – Plot identifier, within unit
- Direction – Sampling direction (cardinal directions)
- dir_code – Direction code (a - d that represent normalized cardinal directions)
- moisture – Soil moisture (volumetric, %)
- Season – Sampling season (Fall or Summer)
- Year – Sampling year
- MT – Plot mycorrhizal legacy type (EcM and AM)
- Treatment – Harvest treatment (Control or Cut)
Corinth_Soil_Test_Results_9-24-24.csv
Laboratory soil chemistry analysis from soil samples.
- Site - ("Unit"_"plot")
- ph - 1:1 soil to water measured in lab
- P, K, Ca, Mg, Fe, Mn, B, Cu, S, Zn, Na, Al: soil extractions for elements (mg/Kg)
- SOM%: The percent soil organic matter
- CEC: Cation exchange capacity
Soil.moisture.correction.csv
Gravimetric soil moisture calculations used to calibrate field soil moisture measurements.
Variables:
- Sample – Sample identifier unique
- tin – Mass of empty tin (g)
- wet + tin – Mass of wet soil plus tin (g)
- dry + tin – Mass of dried soil plus tin (g)
- grams of water – Water mass lost during drying (g)
- moisture – Calculated gravimetric soil moisture
belowground.data.csv
Belowground environmental variables in a processed dataset from multiple files - see accompanying code for PCA.
PC1 - principal component 1
PC2 - principal component 2
- NH4 – Ammonium concentration (mg / L)
- NO3 – Nitrate concentration (mg / L)
- PO4 – Phosphate concentration (mg / L)
- dir_code – Direction code (cardinal directions)
- moisture – Soil moisture (volumetric, %)
- Season – Sampling season (Fall or Summer)
- Year – Sampling year
- MT – Plot mycorrhizal legacy type (EcM and AM)
- Treatment – Harvest treatment (Control or Cut)
- Site - ("Unit"_"plot")
- Date - "season"_"year"
- ID2: "unit"_"plot" _ "dircode"
- time_f: factor for months since harvest
- time: fime_f + 1
- CN: carbon to nitrogen ratio of the bulk soil (unitless ratio)
- ph: 1:1 ratio of soil to water measured in lab
ph_postharvest.csv
Post-harvest soil pH measurements from experimental plots.
- Date - "season"_"year"
- Number: unique number
- ph - 1:1 soil to water measured in lab
Corinth.alldata.csv
Comprehensive dataset that includes variables from multiple files - see accompanying code.
- Unit – Forest unit identifier
- Plot – Plot identifier, within unit
- Direction – Sampling direction (cardinal directions)
- ID: unique ID used in processing code
- number: used in code to combine to other data
- Horizon: OA (organic A) or B
- MT (Am or EcM)
- Treatment – Harvest treatment (Control or Cut)
- Root.wt.g: Root weight (g)
- Bulk density (g / cm^3)
- ph: 1:1 ratio of soil to water measured in lab
- per.C: Soil percent carbon
- per.N: Soil percent nitrogen
- GWC: gravimetric water content (g water/ g soil)
- Biomass: microbial biomass carbon (ug C / g soil); NA entries for missing data
Microbial Biomass and Community Composition
plfa.model.csv
Phospholipid fatty acid (PLFA) and neutral lipid fatty acid (NLFA) biomarker dataset used to estimate microbial biomass and community composition.
- ID_Only — unique sample identifier (unitless)
- AMF_NLFA — arbuscular mycorrhizal fungi neutral lipid fatty acid biomarker concentration (nmol g⁻¹ dry soil)
- AM_PLFA — arbuscular mycorrhizal fungi phospholipid fatty acid biomarker concentration (nmol g⁻¹ dry soil)
- NLFA.PLFA — ratio of AMF_NLFA to AM_PLFA (unitless)
- Total.NLFA — total neutral lipid fatty acids (nmol g⁻¹ dry soil)
- Total.PLFA — total phospholipid fatty acids (nmol g⁻¹ dry soil)
- Gram.Negative — Gram-negative bacterial biomarkers (nmol g⁻¹ dry soil)
- Eukaryote — eukaryotic microbial biomarkers (nmol g⁻¹ dry soil)
- Fungi — fungal biomarkers (nmol g⁻¹ dry soil)
- Fung_2 — alternative fungal biomarker group (nmol g⁻¹ dry soil)
- Gram.positive — Gram-positive bacterial biomarkers (nmol g⁻¹ dry soil)
- Actinomycetes — actinomycete biomarkers (nmol g⁻¹ dry soil)
- AMF.RA — relative abundance of arbuscular mycorrhizal fungi (% or proportion of total PLFA; unitless)
- gram.neg.RA — relative abundance of Gram-negative bacteria (% or proportion; unitless)
- Euk.RA — relative abundance of eukaryotes (% or proportion; unitless)
- Fungi1_RA — relative abundance of fungal biomarker group 1 (% or proportion; unitless)
- Fungi2_RA — relative abundance of fungal biomarker group 2 (% or proportion; unitless)
- Fun1_2_RA — relative abundance of fungal biomarker groups 1 and 2 (unitless)
- gram.pos — relative abundance of Gram-positive bacteria (% or proportion; unitless)
- actinomycetes.RA — relative abundance of actinomycetes (% or proportion; unitless)
- Fungi.Bacteria — fungal-to-bacterial biomass ratio (unitless)
- Predator.Prey — predator-to-prey microbial ratio (unitless)
- Gram..Gram. — Gram-positive to Gram-negative bacteria ratio (unitless)
- Sat.Unsat — saturated-to-unsaturated fatty acid ratio (unitless)
- Mono.Poly — monounsaturated-to-polyunsaturated fatty acid ratio (unitless)
- GNeg.Stress — Gram-negative bacterial stress index (unitless)
- Unit — Forest unit identifier (unitless)
- Plot — plot identifier (within unit)
- Direction — sampling transect or plot direction (categorical; cardinal)
- Number — sample number or replicate identifier (unitless)
- Year — year of sampling
- Season — season of sampling (Fall or Summer)
- months_post_harvest — time since harvest (months)
- ID2 — secondary sample identifier (unitless)
- MT — AM or EcM plot type
- Treatment — Control or cut
- Before.After — sampling period relative to treatment or harvest (categorical; unitless)
PLFA_NLFA_2022.csv
Variables listed and described metadata for "plfa.model.csv"
PLFA_NLFA_condensed.csv
Variables listed and described metadata for "plfa.model.csv"
Fungal Community Sequencing Data
feature_table_2020.csv
ASV abundance table from ITS sequencing of fungal communities from 2020 samples. The column "number" corresponds to the sample ID, and all following columns are ASVs.
feature_table_all_2021.csv
ASV abundance table for fungal ITS sequences from 2021 samples with the same structure as "feature_talbe_2020.csv"
feature_table_2022.csv
ASV abundance table for fungal ITS sequences from 2022 samples with the same structure as "feature_talbe_2020.csv"
feature_table_2022Summer.csv
ASV abundance table for fungal ITS sequences from summer 2022 samples with the same structure as "feature_talbe_2020.csv"
taxon_table_2021.csv
Taxonomic classification of fungal ASVs detected in 2021 sequencing. Each ASV (first column) has corresponding taxonomic information in the subsequent columns from Kingdom to Species. NA values indicate no information from the FUNGUILD databse.
taxon_table_2022Summer.csv
Taxonomic classification of fungal ASVs detected in 2022 summer sequencing with the same structure as "taxon_table_2021.csv"
corinth.taxa.2021.csv
Taxonomic classification of fungal ASVs for 2021 samples with the same structure as "taxon_table_2021.csv"
corinth.taxa.2022Summer.csv
Taxonomic classification of fungal ASVs for 2022 summer samples with the same structure as "taxon_table_2021.csv"
FUNGuild_2021_trimmed
Functional guild assignments for fungal ASVs from 2021 sequencing using the FUNGuild database, with the same structure as "taxon_table_2021.csv" for taxonomic information, and additional columns:
- taxon — taxonomic name assigned to the fungal group (e.g., family, genus, species; unitless)
- guid — globally unique identifier for the taxon (unitless)
- mbNumber — identifier from the MycoBank database (unitless)
- taxonomicLevel — taxonomic rank represented by the record (e.g., family, genus, species; coded value, unitless)
- trophicMode — primary ecological trophic mode (e.g., Saprotroph, Symbiotroph, Pathotroph, or combinations; categorical; unitless)
- guild — ecological guild assignment describing resource use or function (categorical; unitless)
- confidenceRanking — confidence level of the guild assignment (e.g., Possible, Probable, Highly Probable; categorical; unitless)
- growthForm — fungal growth morphology (e.g., Clavarioid, Agaricoid, Corticioid; categorical; unitless)
- trait — additional ecological or functional trait information (categorical; unitless; may contain missing values)
- notes — explanatory notes regarding taxonomy, ecology, or guild assignment (text field; no units)
- citationSource — literature source supporting the ecological classification (text field; no units)
FUNGuild_2022_trimmed.csv
Functional guild assignments for fungal ASVs from 2022 sequencing with the same structure as FUNGuild_2021_trimmed"
FUNGuild_2022_trimmedSummer.csv
Functional guild assignments for fungal ASVs from summer 2022 sequencing with the same structure as FUNGuild_2021_trimmed"
FUNGuild_F20_trimmed
Functional guild assignments for fall 2020 fungal sequencing data with the same structure as FUNGuild_2021_trimmed"
Forest Structure Data
Corinth.DBH.species.csv
Tree diameter measurements used to calculate basal area and characterize canopy composition within plots.
- Unit — experimental unit identifier (unitless)
- Plot — plot identifier (unitless)
- Species — tree species code or species name (categorical; unitless)
- myc.type — mycorrhizal association type (e.g., AM = arbuscular mycorrhizal, ECM = ectomycorrhizal; categorical; unitless)
- Multiplier — expansion factor used to scale individual trees to a per-area basis (unitless)
- DBH — diameter at breast height (cm)
- radius — tree radius at breast height, calculated as DBH / 2 (cm)
- area — stem cross-sectional area at breast height, calculated as π × radius² (cm²)
- BA — basal area, typically expressed on an area basis after applying the multiplier (m² ha⁻¹)
Tree.info.csv
Plot-level forest community metrics from a larger experiment where we use some of their plots. This file was used to evaluate a plot's % basal area of FRAM and ACSA species.
- Plot.ID — unique plot identifier (unitless)
- Plot — plot identifier or plot number (unitless)
- Treatment — experimental treatment category (categorical; control, resilience, demo, exploitative)
- plotrichness — number of species present in the plot (species)
- plotabundance — total abundance of individuals in the plot (count)
- totalplotBA — total plot basal area (m² ha⁻¹)
- FRAMIV —% basal area of FRAM (Fraxinus Americana - American ash) tree species
- FRAMIVbin — binned FRAM presence (categorical; 0 - 5)
- ACSAIV — % basal area of ACSA tree species
- ACSAIVbin — binned or categorical version of ACSAIV (categorical; 0 - 5)
- plotdiversity — plot-level species diversity index (typically Shannon diversity, H′; unitless)
- plotevenness — plot-level species evenness index (unitless)
- Moisture — soil moisture content (%); NA values are for missing data
- Temp.F — soil temperature (°F); NA values are for missing data
- Worm.stage — earthworm invasion stage or earthworm impact class (categorical; 1 - 5)
Corinth.project.alldata.xlsx_-_Site_info.csv
Metadata describing experimental plots, subplot locations, and geographic coordinates.
- Unit — Forest unit identifier (unitless)
- Subplot — subplot identifier within unit (unitless)
- Percent — percent cover, occurrence, or measurement value (%) or AM or EcM trees
- Type — Am or EcM plot type
- long — longitude (decimal degrees)
- lat — latitude (decimal degrees)
- Treatment — experimental treatment category (categorical; control or cut)
- Notes — supplementary notes or comments associated with the observation (text field; no units); Some sites do not have notes.
DNA Metadata
DNA_ID_fall2020_summerfall2021_AFitch.csv
Metadata linking soil samples to sequencing plate wells used during DNA sequencing.
- Site — site identifier or site name (categorical; unit, plot)
- Direction — sampling direction, cardinal
- Date — sampling or observation date (season ,year)
- ID — unique sample identifier (unitless)
- Well — well identifier for 96 well plate; NA indicates missing data.
Code for Data Processing and Analysis: R and R Studio
Corinth__ITS_Processing.Rmd
R Markdown workflow used to process raw ITS sequencing data and generate ASV tables.
Corinth_fungal_community.Rmd
R Markdown analysis of fungal community composition and guild structure.
Corinth_data_analysis.Rmd
R Markdown file containing statistical analyses and figures used in the manuscript.
Access information
Other publicly accessible locations of the data:
- NA
Data was derived from the following sources:
- NA
Site description and experimental setup
Our research takes place in the Clement Woodlot, a northern mixed-deciduous forest in Corinth, Vermont (44°02'48.4" N, 72°19'03.0" W), dominated by overstory tree species sugar maple** (Acer saccharum) and** white ash (Fraxinus americana) across AM-dominated areas, and by eastern hemlock (Tsuga canadensis), yellow birch (Betula alleghaniensis), and American beech (Fagus grandifolia) in EcM-dominated areas. Corinth experiences a mean annual temperature and total annual precipitation of 5.3° C and 1139 mm, respectively (PRISM Climate Group, 2014).
We established 16, 0.1 ha plots, where half of the plots were dominated (>80% basal area (BA)) by EcM-associated trees and half by AM-associated trees (Fig 1). To determine plot mycorrhizal dominance, we measured the diameter at breast height (DBH) for all trees that were > 10cm in diameter, and recorded the species and mycorrhizal type (see Table S1 for plot characteristics including dominant species). Four 0.1 ha plots per mycorrhizal type (referred to as “EcM legacy” and “AM legacy” plots) were patch cut in the winter of 2020-21 and four were left as unharvested controls, for a total of eight 0.1 ha plots in each control and treatment group (Fig. S1). To reduce potential differences among plot locations, we chose locations for control and treatment (harvest) plots with similar slopes, aspects, and dominant tree species within a mycorrhizal type. We measured plot characteristics before harvesting in 2020 across the 16 sites including soil moisture, pH, nutrient availability, and fungal biomass (described below), to employ the Before-After Control-Impact (BACI) framework with an equal number of paired control and treatment sites (Underwood, 1992).
Soil sampling and plot characterization under the BACI design
To assess how harvesting affected soil conditions, we measured soil moisture, pH, fungal biomass, and nutrient availability before timber harvesting and after timber harvesting. We measured volumetric soil moisture with a probe (METER, Teros 10) at three locations per plot when collecting soil samples (October, 2020, June 2021, June 2022). To assess soil pH and fungal communities, we collected 10 cm of the OA horizon prior to timber harvesting (October 2020) and following timber harvesting in both the early growing season (June 2021 and 2022) and at the end of the growing season (October 2021 and 2022), using a 5 cm diameter soil probe in triplicate at each plot after removing the leaf litter (Oi layer) if present. Soils were transported back to the lab on ice, where they were frozen at -20°C until processing. We analyzed soil percent carbon (C) and N on an elemental analyzer (EA-IRMS; Thermo Scientific) to calculate the soil C:N ratio on pre-harvest soils. To measure soil pH, we used a 1:1 ratio of 5 g of MiliQ water and 5 g of thawed soil on a Mettler Toledo S220 pH/ion meter. Finally, we measured extractable calcium (Ca) (Agricultural and Environmental Testing Lab, UVM Extension), but only using samples from the same time point. B horizon data, combined three subplots, sampled in June 2024.
To measure the relative differences in nutrient availability, we measured available soil nitrogen (NH4 and NO3) and phosphate (PO4) using cation and anion exchange resins (Thermo Scientific), held in bags constructed using pantyhose (L’eggs Womens Everyday). We collected these data prior to timber harvesting (October 2020) and following timber harvesting in both the early growing season (June 2021 and 2022) and at the end of the growing season (October 2021). In each plot, we deployed triplicates of each resin bag type at the OA horizon interface and at least 3 cm below the soil surface, for one month. Following field collection, we extracted resin bags using a 2M KCl solution, following the protocol from Lajtha, 1988, and then measured extracts for respective nutrients on an auto analyzer (Lachat instruments). For both soil sampling and resin deployment, we established permanent locations to resample near the same location when possible.
Soil microbial community characterization
To assess potential changes in the microbial community after harvest and in comparison to no harvest controls, we combined an estimate of fungal biomass with internal transcribed spacer (ITS) sequencing to calculate biomass of fungal guilds. First, we measured the relative biomass of fungi by phospholipid fatty acid (PLFA) and neutral lipid fatty acid (NLFA) analyses, performed at the Soil Health Assessment Center at the University of Missouri (Buyer & Sasser, 2012). This method is applicable for quantifying the amount of fungal (PLFA) and arbuscular mycorrhizal fungal (NLFA) biomass (Lewe et al., 2021).
To characterize fungal community composition, we extracted DNA from OA horizon soils using Qiagen PowerSoil Kit (Qiagen Science Inc., Germantown, MD). DNA was sequenced at the Hubbard Center for Genome Studies, University of New Hampshire, for the fungal ITS2 region using ITS4-FUN (AGCCTCCGCTTATTGATATGCTTAART) and 5.8S-FUN (AACTTTYRRCAAYGGATCWCT) primers. The first round of PCR amplified the ITS2 marker along with priming regions for the TruSeq read primers. We performed PCR in 12 µl reactions using 6ul µL of a hot-start, high fidelity polymerase (Kapa 2X HIFI HotStart ReadyMix), 0.7 µL of each primer (at 5 uM), and 2 µL of extracted DNA. Thermal cycling conditions included a 3 minute hot start at 94˚C, 35 cycles of denaturing (94˚C, 0:30), annealing (58˚C, 0:30), and extension (72˚C, 0:30) and a final extension of 7 minutes at 68˚C. Successful amplification was verified using agarose gel electrophoresis (1.5-2.0%). The first PCR reaction was cleaned with a 1:10 dilution. The second round of PCR added the P5 and P7 flow cell adapters to prepare the library for sequencing on an Illumina MiSeq, along with an external set of sample barcodes located between the flow cell adaptors and read primers. PCR was performed in 15 μl reactions using 6ul µL of a hot-start, high fidelity polymerase (Kapa 2X HIFI HotStart ReadyMix), 5 µL of each primer (at 5 nM), and 2 µL of diluted product from the first PCR as template. Thermal cycling conditions were a 3 minute hot start at 94˚C, 15 cycles of denaturing (94˚C, 0:20), extension (72˚C, 0:15) and a final extension of 7 minutes at 72˚C. Amplicons were cleaned with the Qia-quick PCR purification kit. Purified products were quantified using a Qubit 2.0 fluorometer with the Qubit dsDNA HS assay (Thermo Scientific, Grand Island, NY). Amplicons were pooled at equal concentration and sequenced on an Illumina MiSeq using V3 chemistry using paired-end sequencing (300 cycles). Sequences were separated by barcode at the Hubbard Center for Genome Studies, then filtered for quality and assigned to amplicon sequence variants (ASV’s, equivalent to 100% identity operational taxonomic units) using the DADA2 program (Callahan et al., 2016) in R software. Non-singleton ASVs were identified to the lowest confident taxonomic level using the naïve Bayesian classifier RDP using the UNITE database for ITS reads (Nilsson et al., 2019). Since exact sequence variants are not appropriate biological units for fungi (Tedersoo et al., 2022a), we agglomerated ASVs to species level based on taxonomic identification.
While other regions such as LSU and SSU have been preferentially used over the ITS region for identifying AM fungi, some studies have shown that the ITS2 region has successfully characterized potential soil AM communities as well as other regions (Heeger et al., 2019; Kohout et al., 2014; Tedersoo et al., 2022b). We matched ASVs with taxonomic and functional guild information using the FUNguild database (Nguyen et al., 2016) and assigned all fungal species to one of six broad putative categories: arbuscular mycorrhizal (AM), ectomycorrhizal (EcM), plant pathogens/endophytes, saprotrophs, other, and unassigned. The “other” category consisted of species assigned to rarer categories (lichenized fungi, animal pathogens, etc.). We made no attempt to separate plant pathogens from other, non-pathogenic plant endophytes because even individual fungal strains can switch between these lifestyles based on plant host and environmental conditions. For species assigned to multiple categories, for analytical purposes we assigned them to a single guild with the following priority (EMF > plant pathogen/endophyte > saprotroph) for analytical purposes. Because sequencing data of communities are relative (Gloor et al., 2017), we used the fungal fatty acid signature 18:1 w9 from the PLFA data (Kaiser et al. 2010) to estimate EcM, saprotrophic, and pathotrophic fungal biomass by multiplying these data by the summed relative abundances for guild of interest for each sample (Lewe et al., 2021). Measured NLFA values (16:1 w5c) represented the AM fungal biomass (Kaiser et al., 2010; Olsson & Lekberg, 2022).
Seedling survival and growth
From May 14-17, 2021, we planted 20 seedlings each of eight different species per plot, four AM species (Acer saccharum, Acer rubrum, Nyssa sylvatica, and Prunus serotina), and four EcM species (Quercus rubra, Carya cordiformis, Tilia americana, and Betula lenta). We chose to include an equal number of AM and EcM-associated species that aligned with species of interest from regional forest managers, as well as potential climate projections for the region. For planting locations, we split the plot into equal halves for AM and EcM seedlings (Fig 2), and then randomly assigned locations to the four species within a grid. Within the grid, the seedlings were planted at a minimum distance of 1 m from other seedlings and were caged to prevent deer browsing effects on seedling growth and survival. We then measured seedling height and root collar diameter, noted whether seedlings were planted on scarified soil (disturbed from logging equipment), and categorized the level of slash present within half a meter around the seedling on a scale from 0-4 on increments of 25% (0 representing no slash, and 4 representing 75-100% slash coverage). At the end of each growing season, August 30-September 1, 2021 and August 29-August 31, 2022, we measured seedling growth (root collar diameter and standing height to the apical bud), documented seedling survival (dead = 0, or alive = 1), documented the “condition” of the terminal bud as intact or missing (1 = missing, 0 = intact), and noted herbivory on a scale from 0-3 (0 = no herbivory, 1 = slight herbivory, 1 = moderate herbivory, and 3 = major herbivory and/or missing apical bud). If a new leader had sprouted due to a missing terminal bud, we measured the height to the top of the new leader.
To estimate green biomass accumulation, we followed equations from Clark et al., (2022), where seedling volume was assumed to be conical and calculated based on the volume of a cone and calculated green biomass using constants from Miles & Smith (2009). Finally, we calculated the relative growth rates as the natural log of change in biomass from initial measurements (May 2021) to the end of the second season of growth (2022), divided by time from seedling planting to measurements at the end of the second season (Hunt & Cornelissen, 1997).
Leaf 15N stable isotopes
To further investigate how access to available nitrogen could affect seedling survival and growth, we measured ratios of δ15N/14N (‰) in seedling leaves after the first season of growth in fall 2021. The ratio of δ15N/14N gives an indication of the nitrogen source, where a lower ratio of δ15N/14N can be an indication of fractionation of N by mycorrhizal fungi as they transfer N to seedling roots (He et al., 2009; Hobbie & Högberg, 2012). We took duplicate hole punches from two seedlings per species within 8 m of the plot center and within 8 m of the plot edge in each plot. We foil-balled duplicate ¼ cm diameter leaf punch samples, and measured them on an isotope ratio mass spectrometer (IRMS, Thermo Scientific Flash EA with Delta V Advantage).