Root exudation affects soil biogeochemistry profoundly, yet it is rarely quantified in mature, field-grown trees, and its controls are poorly understood.

We measured rates of carbon (C) exudation in 11 tree species that exhibit divergent root traits, including gymnosperms and angiosperms that associate with either arbuscular mycorrhizal (AM) or ectomycorrhizal (EcM) fungi. Our goal was to explore how tree species, plant functional groups, and root traits collectively influence exudation patterns.

Species-level differences were modest owing to substantial intraspecific variability (among individuals of the same species). However, three of the four highest exuding species were EcM gymnosperms, which exuded ca. two times more C than the other functional groups. Principal component analysis revealed that relationships between root trait organization and exudation were pronounced in EcM-associated trees but weak or absent in AM-associated trees. In EcM trees, exudation rates were negatively correlated with root tissue density (RTD) and positively with specific root area (SRA) and root diameter, driven largely by gymnosperms. In contrast, exudation in AM trees showed only a weak association with specific root length. Consistent with these findings, mixed-effect models also showed that exudation rates were best explained by a combination of tree-mycorrhizal type, phylogenetic group and SRA, though a large portion of unexplained variation suggests an important role for contemporary environmental and local edaphic conditions.

Collectively, our results demonstrate that root exudation is a complex physiological process shaped by interactions among mycorrhizal association, evolutionary history, and root traits, rather than by functional groups or root traits alone. These findings highlight the urgent need for more integrative frameworks and new experimental approaches to incorporate exudation dynamics into plant strategy theory and large-scale ecosystem models.

Site description

This study was conducted in monoculture plots at the Morton Arboretum, Lisle, Illinois (41.81N, 88.05W). The plots were established between 1922 and 1948 to test and study “all the timber trees of the world which might come under consideration for reforestation purposes in this part of the country” (Morton Arboretum Staff, 1929). Soils in the plots are poorly drained Alfisols that form from a thin layer of loess (0.31 m) underlain by glacial till and Mollisols that formed from alluvium (Soil Survey Staff, NRCS, USDA, 2024). The soil series in the plots is primarily Ozaukee silt loams and Sawmill silty clay loam (Midgley & Sims, 2020). The area has a continental climate with temperatures ranging from -6°C in January to 22℃ in July and 800-1,000 mm mean annual precipitation.

Eleven tree species were selected to capture the heterogeneity in root traits among species from distinct functional groups: phylogenetic group (angiosperm (5) vs. gymnosperm (6)), tree-mycorrhizal association (AM (6) vs. EcM (5)), and leaf habit (deciduous (7) vs. evergreen (4)). Within each group, species were chosen based on mean SRL and root tissue N content (root N) - the traits that were found to correlate positively with exudation rates in previous studies (Meier et al., 2020; Sun et al., 2021; Wang et al., 2021). As such, the selected eleven species spanned a wide range of SRL and root N for each group, ensuring that the species captured diverse trait space (Table 1). This allowed for minimizing phylogenetic covariations among traits while maximizing species trait dissimilarities. Out of the eight combinations, only two combinations were absent: evergreen-AM-angiosperms and evergreen-EcM-angiosperms (Table 1).

Root exudation rates

Fine-root exudates were collected during the growing season of 2022 (i.e., from May to July 2022) using an in-situ culture-based cuvette system (Phillips et al., 2008). To mitigate the impact of variable weather, sampling campaigns were conducted under sunny and clear conditions when possible. Each plot was visited twice: once in late May/early June to collect exudates from three individuals of the focal species and once in late June/early July to collect exudates from 2-4 additional individuals of the same species. The terminal roots were excavated carefully from the mineral topsoil below the organic layer. The excavated root segments were examined to ensure that the fine-root system consisted of the first three branching. Organic matter and soil particles adhering to the root system were removed with DDI water with extreme caution while keeping the roots moist with wet paper towels. In cases where the distal fine roots were damaged or broken off, samples were discarded, and a new sample was prepared. The intact root systems were placed in cuvettes (30 mL syringe) filled with sterile, C-free glass beads (>1 mm diameter). The root systems with glass beads were flushed three times with C-free nutrient solution (0.5 mM NH₄NO₃, 0.1 mM KH₂PO₄, 0.2 mM K₂SO₄, 0.15 mM MgSO₄, 0.3 mM CaCl₂) to ensure the root segments and glass beads were well-mixed and to remove any C adhering to the root surface. To ensure the same amount of solution was added to the cuvettes, we added 15mL of nutrient trap solution in the field using a bottle-top dispenser. The cuvette was covered in aluminum foil to allow the root system to equilibrate with the cuvette environment. The same procedure was applied to the control (i.e., no root) cuvette with the same glass beads and nutrient solution. The cuvettes were placed at the exact excavated area and covered with soils and organic matter and incubated for approximately 24 hrs.

After the one-day incubation period, the sampled roots with the cuvette were clipped with care and brought to the laboratory for analysis. Within one hour of clipping, each cuvette was flushed with 15mL of the working nutrient solution three times to remove accumulated exudates in the cuvette. All solutions were filtered immediately through a sterile 0.22 μm syringe filter (Millex-GV 0.22µm PVDF 33mm Gamma Sterilized 50/Pk, Millipore Co., Billerica, MA) and refrigerated at 4°C until analyses (<24 h). All samples were analyzed for non-particulate organic C on a TOC analyzer (Shimadzu Scientific Instruments, Columbia, MD) within a day of sample collection. The total mass-specific exudation rate was calculated with the total C captured from the trap solution minus the total C flushed from the root-free control cuvettes divided by the dry root biomass and day (mg C * g _root^-1 * day^-1).

Root morphological and chemical traits

Roots originally placed in the cuvette were carefully collected from the cuvette, washed, and stored at 4°C until processing. Fine-root morphology was analyzed for all the fine roots with a transparent flat-bed scanner and the WinRHIZO program (Regent Instruments, Quebec, QC, Canada). Scans were collected at a resolution of 600 dpi. All root samples were dried at 65°C for at least 48 h, and the dried root biomass was used for root trait calculations. Specific root length (SRL, in m g^-1 : the length of the fine roots divided by the corresponding root dry weight), specific root area (SRA, in cm² g^-1 : the area of the fine roots divided by the corresponding root dry weight), root tissue density (RTD, in g cm^-3 : root dry weight divided by root volume), root branching intensity (BI, in the number of tips per total fine-root length), and root diameter (diameter, in cm) were calculated from WinRHIZO. Root N content (% per dry weight) was measured independently in the lab using an elemental combustion system (Costech Analytical Technologies).

Statistical analyses

We used an analysis of variance (ANOVA), linear mixed-effects models, and variance partitioning to characterize the extent to which root exudation rates vary among tree species and across functional groups. To test for differences in exudation rates among tree species, we conducted pairwise comparisons after an ANOVA using a Tukey’s Honest Significant Difference (HSD) test. To test for differences in exudation rates among tree functional groups, we used a mixed-effects model with mycorrhizal type, phylogeny, and their interaction as fixed effects and species-plot as a random effect using restricted maximum likelihood (‘lme4::lmer’ via REML). To evaluate the significance of each nested group in the model after accounting for all other groups, Type III ANOVA with Satterthwaite's Method using the ‘lme4::anova’ was performed to summarize the results of each model. To control the likelihood of false positives in all linear mixed effects models, adjusted p-values from BH Correction (Benjamini-Hochberg) test were performed using the p.adjust function. To quantify the contributions of inter- vs. intraspecific variation to exudation rates in mixed effects models, a variation partitioning analysis was performed using the ‘VEGAN::varpart’. To show co-variations among root traits, a pairwise trait relationships between exudation rates and root traits were also performed using Pearson’s correlations at the individual tree level using ‘corr.test’ function. Root traits and exudation rates were natural-log-transformed prior to analyses to meet model assumptions of residual normality and homogeneity of variance.

To assess how and the extent to which root exudation rates are associated with root trait coordination, we used principal components analysis (PCA; Weigelt et al., 2023) and Redundancy Analysis (RDA). To examine how exudation rates align with major dimensions in the PCA, we created an ordination of RTD, root N, SRL, Diameter, SRA, and BI along with exudation rates using princomp () with standardized PCA. We evaluated the degrees to which each principal component (PC) contributed to the ordination with Broken Stick analysis and tested the statistical significance of each PC with Horn's parallel analysis using ‘PARAN::paran’. To examine the significance of linear relationships between exudation rates and the first four axes, we created a PCA without exudation rates and performed Pearson’s product-moment correlation test between PCs and exudation using cor.test function. To select the best predicting root trait or subset of predictors, we built a PCA with six core variables (RTD, root N, SRL, BI, SRA, and Diameter) and evaluated the relationship between root exudation rates as a trait and the traits that comprise the PCA. We used RDA models for PCs to explain exudation using ‘VEGAN::rda’ and selected the best predicting trait using ‘VEGAN::ordistep’ with both forward and backward stepwise model selection. Given that phylogeny and mycorrhizal association emerged as significant predictors during model selection, we constructed separate PCAs for AM and EcM trees and repeated all analyses to assess whether the organization of root traits predicts exudation differently between the two mycorrhizal types. To test the significant difference among two phylogenetic group in their root trait organization and exudation, we examined the differences in centroid positions in the PCAs through Permutational multivariate analysis of variance (PERMANOVA), using ‘ADONIS::adonis2’ and spread around the centroids using ‘VEGAN::betadisper’.

To identify the functional groups and root traits that collectively predict exudation rates, we used a stepwise model selection approach using linear mixed-effects models by ‘lme4::lmer’ via REML. The fixed effects included six root traits (RTD, root N, SRL, Diameter, SRA, and BI) along with mycorrhizal type or phylogenetic group. Monodominant plot identity (i.e., species-plot) was treated as a random effect. Model selection was based on improvements in Akaike Information Criterion (AIC) and likelihood ratio tests comparing full and reduced models. To highlight the impacts of mycorrhizal type on trait-exudation relationships, we repeated the same analyses to identify the best predicting mixed-effect models for AM and EcM trees. Building on the best-performing model, we further tested interactions between traits and functional groups (e.g., Exudation ~ Trait × Functional Group). We examined the explanatory power of each model by calculating marginal (R²m) and conditional (R²c) R-squared values, where R²m represents variance explained by fixed effects and R²c includes both fixed and random effects (Nakagawa & Schielzeth, 2013). Model assumptions for selected models were verified via checks for residual normality, homoscedasticity, and unbiasedness. All statistical analyses were performed using R v.3.5.3 (R Core Team, 2017).