Quantifying the effects of switchgrass ( Panicum virgatum ) on deep organic C stocks using natural abundance 14 C in three marginal soils

Abstract

Perennial bioenergy crops have been shown to increase soil organic carbon (SOC) stocks, potentially offsetting anthropogenic C emissions. The effects of perennial bioenergy crops on SOC are typically assessed at shallow depths (<30 cm), but the deep root systems of these crops may also have substantial effects on SOC stocks at greater depths. We hypothesized that deep (>30 cm) SOC stocks would be greater under bioenergy crops relative to stocks under shallow‐rooted conventional crop cover. To test this, we sampled soils to between 1‐ and 3‐m depth at three sites in Oklahoma with 10‐ to 20‐year‐old switchgrass (Panicum virgatum) stands, and collected paired samples from nearby fields cultivated with shallow rooted annual crops. We measured root biomass, total organic C, ¹⁴C, ¹³C, and other soil properties in three replicate soil cores in each field and used a mixing model to estimate the proportion of recently fixed C under switchgrass based on ¹⁴C. The subsoil C stock under switchgrass (defined over 500–1500 kg/m² equivalent soil mass, approximately 30–100 cm depth) exceeded the subsoil stock in neighboring fields by 1.5 kg C/m² at a sandy loam site, 0.6 kg C/m² at a site with loam soils, and showed no significant difference at a third site with clay soils. Using the mixing model, we estimated that additional SOC introduced after switchgrass cultivation comprised 31% of the subsoil C stock at the sandy loam site, 22% at the loam site, and 0% at the clay site. These results suggest that switchgrass can contribute significantly to subsoil organic C—but also indicated that this effect varies across sites. Our analysis shows that agricultural strategies that emphasize deep‐rooted grass cultivars can increase soil C relative to conventional crops while expanding energy biomass production on marginal lands.

INTRODUCTION

Soil horizons deeper than 30 cm contain the majority of Earth's soil organic carbon (SOC)—possibly holding well over 1,000 Pg of C globally (Harrison, Footen, & Strahm, 2010; Jobbágy & Jackson, 2000). While the bulk of deep soil C tends to exchange slowly with the atmosphere (Mathieu, Hatté, Balesdent, & Parent, 2015; Trumbore, 2009), SOC losses from deep soil horizons following land use change have been substantial—accounting for the majority of the 133 Pg of SOC lost following the global expansion of agriculture (Sanderman, Hengl, & Fiske, 2017). By implication, successful attempts to reverse SOC loss in agricultural lands must restore SOC in deep horizons. Furthermore, C concentrations at depth are relatively low—implying that subsoils have a large capacity to store C and thus might sequester a significant amount of additional atmospheric CO₂ (Lorenz & Lal, 2005; Minasny et al., 2017; Paustian et al., 2016; Rumpel & Kögel‐Knabner, 2011).

A range of processes introduce C to subsoils, including dissolved C transport in percolating water, burial of aboveground litter via physical mixing, and C fluxes from root exudates and root turnover at depth (Rumpel & Kögel‐Knabner, 2011). Deep roots in particular have been identified as a potentially useful conduit for increasing subsoil C stocks in managed landscapes (Kell, 2012; Lynch & Wojciechowski, 2015). A large fraction of SOC is root derived, and the depth distribution of SOC correlates with rooting distributions across biomes in natural ecosystems (Grayston, Vaughan, & Jones, 1997; Jobbágy & Jackson, 2000; Rasse, Rumpel, & Dignac, 2005). Dead roots and root exudates fuel production of microbial biomass, which subsequently becomes a primary source of mineral‐associated C that can persist over long timescales (Sokol, Kuebbing, Karlsen‐Ayala, & Bradford, 2019). Deeply rooted bioenergy crops can also enhance production of microbial extracellular polysaccharides, cementing soil aggregates that may protect SOC (Sher et al., 2020). In theory, increasing SOC via deep roots might be achieved without displacing conventional food crops if bioenergy crops are grown on marginal lands—which are otherwise not ideal for food production due to low fertility or environmental sensitivity (Gelfand et al., 2013; Lemus & Lal, 2005; Robertson et al., 2017).

While cultivation of perennial bioenergy crops and restoration of perennial grasslands have been widely shown to increase SOC stocks relative to conventional crops, the majority of studies have focused on the top 30 cm of soil (Anderson‐Teixeira, Davis, Masters, & Delucia, 2009; Beniston, DuPont, Glover, Lal, & Dungait, 2014; Conant, Cerri, Osborne, & Paustian, 2017; Harris, Spake, & Taylor, 2015; Monti, Barbanti, Zatta, & Zegada‐Lizarazu, 2012; Qin, Dunn, Kwon, Mueller, & Wander, 2016). Furthermore, the magnitude of the difference in SOC stocks following conversion to perennial grassland is highly variable (Conant et al., 2017). Predicting the effect of deep roots on subsoil C across different soil types will ultimately require more field studies spanning edaphic gradients that sample deeply (i.e., ≥1 m).

Evaluating the effects of deep roots on subsoil C in the field is challenging, however, because differences in SOC stocks between different land use types are often small relative to total SOC stocks (Syswerda, Corbin, Mokma, Kravchenko, & Robertson, 2011). Ideally, changes in SOC under different plant types would be quantified in long‐term experiments in which initial conditions are controlled and quantified (Liebig, Schmer, Vogel, & Mitchell, 2008; Sanford et al., 2012). An alternative is to sample opportunistically using a paired design (Fisher et al., 1994; Liebig, Johnson, Hanson, & Frank, 2005); in this case, the plant cover of interest is compared to a neighboring “reference field” representing the conventional management practice and initial conditions are assumed to be the same across the two plots. This approach cannot detect net change SOC over time given that SOC stocks in the reference plot may not be at steady state—but it can detect divergence in SOC stocks under different management scenarios (Sanderman & Baldock, 2010). Furthermore, the paired design can be applied rapidly in locations where initial data are unavailable, enabling wider sampling of edaphic gradients.

Naturally occurring C isotopes (¹³C, ¹⁴C) can be used as sensitive tracers of C fluxes (Jones & Donnelly, 2004), and are useful for constraining the effect of deep roots on subsoil C when the paired sampling approach is applied (Balesdent et al., 2018; Marin‐Spiotta, Silver, Swanston, & Ostertag, 2009; O'Brien, Jastrow, McFarlane, Guilderson, & Gonzalez‐Meler, 2013; Richter, Markewitz, Trumbore, & Wells, 1999). For instance, ¹³C is commonly used to quantify the fraction of SOC derived from recent plant inputs in cases where the photosynthetic pathway of the plant cover is replaced, changing the ¹³C signature of the inputs (Balesdent et al., 2018; Balesdent, Mariotti, & Guillet, 1987; Garten & Wullschleger, 2000). However, ¹³C‐based mixing models require a clear transition between C₃ and C₄ vegetation (Balesdent & Mariotti, 1996), and are thus challenging to apply in agricultural systems with complex cropping histories.

In systems where no clear transition between C₃ and C₄ vegetation has occurred, the radioisotope ¹⁴C provides an alternative to ¹³C. Atmospheric radiocarbon concentrations are sustained by production of ¹⁴C in the stratosphere, and were elevated by introduction of ¹⁴C from atomic weapons testing during the 1950s and 1960s (Hua, Barbetti, & Rakowski, 2013). Deep soil C exchanges slowly with the atmosphere and thus becomes naturally depleted in ¹⁴C as it undergoes radioactive decay (Trumbore, 2009). Consequently, recently fixed C introduced to subsoils via increased belowground inputs should have an elevated ¹⁴C signature relative to the preexisting subsoil C pool (Richter et al., 1999). ¹⁴C can thus provide upper limits on the magnitude of differences in SOC that emerge after replacing conventional crops with deeply rooted crops.

In this paper, we explore C storage in marginal lands cultivated with switchgrass (Panicum virgatum L.), a deeply rooted perennial grass grown as forage and as a cellulosic bioenergy feedstock. We used a paired sampling design at three sites in Oklahoma with different soil textures that experienced soil degradation during the American Dust Bowl and were planted with switchgrass in either 1998 or 2008 and sampled in 2018. Given that 10 years is typically sufficient to measure C stock differences at shallow depths (<30 cm) when comparing switchgrass to conventional cropland (Anderson‐Teixeira et al., 2009), we hypothesized that C stocks at greater depths (>30 cm) would also diverge between switchgrass and paired reference plots over this timespan. Identifying rates of SOC divergence in subsoils under perennial bioenergy crops is important because the majority of existing studies on land conversion to perennial crops still deal with relatively shallow sampling depths: increasing the number of studies that sample deeply is an imperative for improving regional‐ to global‐scale prediction of perennial crop effects on SOC (Ledo et al., 2020). We tested our hypothesis by quantifying both total C and ¹⁴C, which we used to develop sensitive estimates of the component of the total C stock that could be attributed to switchgrass.

MATERIALS AND METHODS

Field sites

Sampling took place in 2018 at three sites in Oklahoma, USA. At each site, we sampled deep soil cores in >10‐year‐old switchgrass plots and compared these with paired cores collected from nearby fields cultivated with annual crops. The two sites in Southern Oklahoma; Red River farm, Burneyville (hereafter the “Sandy Loam” site; Lat: 33°53ʹ20.52ʺN, Lon: 97°17ʹ7.13ʺW) and Pasture Demonstration Farm, Ardmore (hereafter the “Clay” site; Lat: 34°13ʹ11.00ʺN, Lon: 97°12ʹ36.96ʺW) had been planted with “Alamo” switchgrass in 2008. The location in Northern Oklahoma, near Stillwater (hereafter the “Loam” site; Lat: 36°8ʹ0.16ʺN, Lon: 97°6ʹ15.42ʺW), was planted with “Kanlow” switchgrass in 1998. At the Sandy Loam site, switchgrass was uncut, whereas at the Clay and Loam sites switchgrass was mowed and harvested annually (Loam) or 1–2 times annual (Clay). The switchgrass stands at each site were unfertilized, although the stands at the Clay site were originally established as part of a short‐term P response study and thus received fertilizer initially after planting. All three sites were near the outer geographic boundary of the American Dust Bowl during the 1930s and likely experienced wind erosion at that time. Before European settlement, the region likely hosted tall‐grass prairie dominated by C₄ grasses (Cotton, Cerling, Hoppe, Mosier, & Still, 2016). After European settlement in the 19th century, soils in the region were cultivated with C₃ cereal crops (Paulsen & Shroyer, 2008). The three sites have a broadly similar mean annual climate (Table 1).

1 TABLELocation and climate of study sites. Mean annual temperature (MAT) and mean annual precipitation (MAP) were obtained using gridded PRISM climate data (Prism Climate Group, Oregon State University, 2011)

Site	Latitude	Longitude	MAT (°C)	MAP (mm)
Sandy Loam	33°53ʹ20.52ʺN	97°17ʹ7.13ʺW	17	954
Loam	36°8ʹ0.16ʺN	97°6ʹ15.42ʺW	16	933
Clay	34°13ʹ11.00ʺN	97°12ʹ36.96ʺW	17	959

At the Sandy Loam site, the reference field had been cultivated with the C₃ grass rye (Secale cereal L.) in winter and the C₄ plant crabgrass (Digitaria sanguinalis (L.) Scop.) in the summer for at least the last 15 years under no‐till management. Nitrogen fertilizer was applied at approximately 150 kg N/ha in the reference field annually at this site. At the Clay site, the most recent species grown in the paired reference plots was wheat (Triticum aestivum L.) with a winter cover crop mix; this site was managed with conventional tillage, N was applied at an average rate of 67 kg N/ha annually, and fields were grazed by cattle in winter. At the Loam site, the reference field was typically planted with wheat and managed with conventional tillage—although during several years throughout 1998–2018 the reference field was planted with the C₄ grass sorghum (Sorghum bicolor (L.) Moench); N was applied at an average rate of 72 kg/ha annually. To our knowledge, none of the sites were limed.

The sites spanned a soil texture gradient driven by parent material composition. The Sandy Loam site featured coarse alluvial soils (Coarse‐loamy, mixed, superactive, thermic Udic Haplustolls; National Cooperative Soil Survey, 2020). The Loam site featured soils derived from alluvial and eolian deposits (Fine‐loamy, mixed, superactive, thermic Fluventic Haplustolls; National Cooperative Soil Survey, 2020). Notably, the soils at this site included a buried soil (paleosol) at >1 m depth. The Clay site included a range of relatively fine‐textured soils weathered from Permian shales and sandstones (Fine‐loamy, mixed, active, thermic Udic Argiustolls; National Cooperative Soil Survey, 2020). Soils at the third site varied between clays, clay loams, and sandy clay loams based on the USDA texture classification system; we chose the label “Clay” for this site because it was the most common texture class.

Field sampling

At each site, three soil cores were collected under switchgrass and three cores were collected in an adjacent reference field. We treated the three cores taken in each field as replicate samples, but we acknowledge that these cores are “pseudo‐replicated” in that they were collected from the same field, and that a larger sample size would have been ideal (Kravchenko & Robertson, 2011). The low sample size was necessitated by the larger amount of labor required to process >3 m soil cores and the costs of radiocarbon analyses. Cores were spaced apart in each field so that they would capture within‐field variation to the extent possible: cores at the Sandy Loam and the Clay sites were collected in June 2018, from a 20 m² area within each field, and at the Loam site cores were collected October 1, 2018, also within a 20 m² area. The reference fields at the Clay and Loam sites were approximately 50 m distant from the switchgrass fields, and the reference field at the Sandy Loam site was approximately 500 m distant but situated in the same soil series. Cores at the Sandy Loam and Cores were taken using a Giddings probe with 10.16 cm (4 inch) inner‐diameter tooling and sampling 120 cm intervals. Sampling proceeded to a depth of 3 m unless the probe reached refusal at a shallower depth (this occurred at the Clay site at a depth of 120–150 cm, likely due to calcium cementation at depth). Each core was photographed and divided into 30 cm intervals in the field. The reference plots were chosen to match the soil properties of the switchgrass plots based on field observations.

At all sites, bulk density was estimated by weighing a 4 cm subsample from the center of each core interval in the field and correcting for the gravimetric water content of the subsample to obtain the subsample dry mass. This mass was then divided by the volume of the subsample to calculate bulk density for that interval. Compression during sampling was on‐average 7 ± 5% at the Clay site and <1% at the Loam and Sandy Loam sites. To account for compression during sampling, volumes were linearly corrected over each sampling interval by scaling the observed core length to the expected length (Parfitt et al., 2010). Particles >2 mm comprised a negligible fraction of the total mass of each interval, and so no correction for rock fraction was performed. C stock calculations were later performed on an equivalent mineral soil mass basis to minimize sensitivity to bulk density estimates (see below).

Roots were removed from the bulk density subsample by hand; approximately 20 person‐minutes were spent removing roots per interval. Roots were washed and dried to obtain the root mass in each interval and scaled by the volume of the interval to obtain root biomass estimates. Soil used for total C and C‐isotope analysis was sampled from the remainder of the core interval after removing 1 cm from its exterior to exclude soil from upper horizons that might have contaminated the interval during sampling. Soil sampled from the interior of the core was sieved to 2 mm and dried at 105°C before being subdivided for physical and chemical analyses.

Laboratory analyses

Soil physical and chemical analyses were conducted at the Oregon State University Crop and Soil Science Central Analytical laboratory (https://cropandsoil.oregonstate.edu/cal). Total C and N were quantified by combustion at 1,150°C using an Elementar Macrocube analyzer. Soil texture analysis, soil pH, and exchangeable cations were also quantified on samples from the 0–30, 30–60, and 60–90 cm depth intervals and select intervals at greater depths. Texture was quantified by the sieve and pipette method after removal of organic matter and carbonates (Soil Survey Staff, 2014). Soil pH was measured by electrode in a 1:1 soil:water slurry. Exchangeable cations were quantified by 0.1 M barium chloride extraction and analysis by ICP‐OES (Soil Survey Staff, 2014).

Inorganic C was quantified at Lawrence Livermore National Laboratory by treating finely ground subsamples of each sample with 1 M phosphoric acid in a sealed jar and measuring CO₂ evolved using a LI‐850 infrared gas analyzer (Robertson, 1999). Where carbonates were present, total organic C was obtained by subtracting inorganic C from total C.

C isotopes were quantified on a subset of the soil that was ground to a fine power by hand. Soils that contained carbonates were treated with 1 M HCl to remove inorganic C before isotope analysis. Direct addition of dilute (~1 M) HCl has measurable but relatively small (<1‰) effects on ¹³C and ¹⁴C in soils and sediments (Brodie et al., 2011; Komada, Anderson, & Dorfmeier, 2008) and appears to be no more biased than alternative treatment approaches (Brodie et al., 2011). HCl was added to each sample until effervescence ceased and then was allowed to evaporate to prevent leaching of acid‐soluble C. Acid‐treated soil was analyzed for ¹³C at the University of California Berkeley Center for Stable Isotope Biogeochemistry (https://nature.berkeley.edu/stableisotopelab/). The ¹³C content of each sample (δ¹³C) was reported in per mil (‰) relative to the V‐PDB isotope standard. Radiocarbon values were measured on the NEC 1.0 MV Tandem accelerator mass spectrometer (AMS) or the FN Tandem Van de Graaff AMS at the Center for AMS at Lawrence Livermore National Laboratory. Samples were prepared for ¹⁴C measurement by sealed‐tube combustion to CO₂ in the presence of CuO and Ag and then reduced onto iron powder in the presence of H₂ (Vogel, Southon, Nelson, & Brown, 1984). The ¹⁴C content of each sample (∆¹⁴C) was reported in ‰ relative to the absolute atmospheric ¹⁴C activity in 1950. We report ∆¹⁴C here rather than mean residence times because reporting ∆¹⁴C does not require an assumption that SOC pools are at equilibrium; negative ∆¹⁴C values generally indicate less interaction between SOC and the atmosphere and longer SOC residence times. To calculate ∆¹⁴C, measured δ¹³C values were used to correct for mass‐dependent fractionation to yield ¹⁴C activity at a reference δ¹³C of −25‰ (Stuiver & Polach, 1977). Radiocarbon analyses were conducted in late 2018–2019 (exact dates are listed for each sample in Table S1). Because collection and analysis occurred within a short period, no correction was performed for decay of ¹⁴C between sampling and analysis. The average instrument uncertainty for ∆¹⁴C was ±4‰, and the average precision estimated from a set of six duplicate samples was ±5‰.

C stock calculations

We used measured C stocks to directly estimate the net difference in C between the switchgrass and reference fields. We also used ¹⁴C measurements to develop an indirect estimate that was independent of the measured C stock in the reference field. The C stock calculations were carried out on an equivalent soil mass (ESM) basis using the cumulative coordinate approach (Gifford & Roderick, 2003; Rovira, Sauras, Salgado, & Merino, 2015). We used this approach because it is robust to differences in bulk density, and thus better suited to comparing C stocks under different land uses (Wendt & Hauser, 2013). Calculations were performed separately on the surface soil layers—which we defined as the top 500 kg/m² of soil—and the subsoil—which we defined as the 1,000 kg/m² of soil directly below the uppermost 500 kg/m² of soil.

We obtained C stocks using linear interpolation to predict cumulative C mass from cumulative soil mass (Gifford & Roderick, 2003). The mineral mass of each depth interval was used as the basis for developing mass coordinates (Rovira et al., 2015). Mineral mass was obtained by multiplying the mass of the interval by the 1 minus the soil organic matter fraction [soil organic matter fraction = % organic carbon × (1/100) × 2; Pribyl, 2010]. We then used linear interpolation to develop a piece‐wise function defining cumulative OC mass as a function of cumulative mineral soil mass (Gifford & Roderick, 2003):[Image Omitted. See PDF]where C(t) is the cumulative C mass at the target cumulative soil mass M(t), C(z_a) and C(z_b) are the cumulative C masses at the upper and lower boundaries of the sampling interval containing M(t), and M(z_a) and M(z_b) are the cumulative mineral masses at those boundaries (Gifford & Roderick, 2003). Using this approach, we estimated topsoil C contained in the first 500 kg/m² of soil, and then obtained subsoil C by calculating the total C stock to 1,500 kg/m² and subtracting the topsoil C stock. Isotopic values for the topsoil and subsoil were calculated by weighting the values associated with each sampling layer by the contribution of that layer to the C stock. When the lower boundary of the topsoil or subsoil occurred within a layer, isotopic values from that layer were weighted by the C mass that contributed to the topsoil or subsoil.

Isotope calculations

We initially explored the use of ¹³C as a quantitative tracer of switchgrass inputs in our system. The mixed history of C₃ and C₄ vegetation at all three sites—and in particular the recent history of periodic C₄ cropping at the Sandy Loam and Loam sites—suggested that our sites did not experience a clear transition between vegetation types. Depth weighted average δ¹³C values for the subsoil (defined over 500–1,500 kg/m² ESM) in the reference plots at our sites ranged between 16.1‰ and −14.9‰, which is at the higher end of the C₄ plant range (O'Leary, 1988). We measured the δ¹³C of switchgrass roots at the three sites and obtained a range of −13.73‰ to −13.34‰—indicating that the difference between isotopic end‐members in a potential ¹³C‐based mixing model in the subsoil was only 2–3‰. This range is comparable to ~2‰ fractionation effects that apply to plant‐tissue end‐members in isotopic mixing models and are a possible source of uncertainty (Menichetti et al., 2015; Werth & Kuzyakov, 2010). Given these clear limitations, we concluded that δ¹³C—while useful for qualitative interpretation of the SOC depth profiles at our sites—could not be used for identifying switchgrass contributions to SOC quantitatively.

Instead of ¹³C, we used ¹⁴C to develop estimates of the amount of C introduced to subsoils by switchgrass that were independent of the observed C stocks in the reference plots. The ¹⁴C signature of plant inputs depends on the composition of the atmosphere, and is thus identical in switchgrass and reference plots. Consequently—while root‐derived inputs are presumably lower under the reference vegetation—some atmospheric ¹⁴C is introduced into the subsoil in both cases, and ¹⁴C can be used to identify net differences in C when comparing the two plots. This contrasts with ¹³C, which is typically used to estimate gross contributions of recently fixed C in the context of paired sampling (Balesdent & Mariotti, 1996).

We did not carry out ¹⁴C‐based calculations for the uppermost 500 kg/m² of soil (approximately 30 cm depth) because the ∆¹⁴C values of the uppermost 500 kg of soil in the reference plots were similar to the range of ∆¹⁴C value of the recent atmosphere at two of the sites. Specifically, we obtained empirical 95% confidence intervals for the ∆¹⁴C value of the uppermost 500 kg of soil using Monte Carlo sampling (see Section 2.6) spanning [−74, 14]‰ at the loam site, [−165, −71] at the Clay site, and [−160, −11]‰ at the Sandy Loam site. At two sites, these intervals approached or overlapped the ∆¹⁴C of the recent atmosphere (assumed to be −7‰ in 2018; Hua et al., 2013), indicating little separation between the isotopic end‐members at the surface. This suggests that ¹⁴C may only be a useful tracer of increased root inputs at depth, where SOC tends to be ¹⁴C depleted and contrasts strongly with recent inputs.

We divided the subsoil SOC stock under switchgrass (C_s, kg C/m²) into two parts: (a) a component equal to the C stock under the reference plot (C_r, kg C/m²), representing the initial C stock plus the mass of C equal to what was accrued or lost under the reference vegetation since 1998 or 2008; and (b) a component equal to the additional or “new” C accrued under switchgrass since 1998 or 2008 (C_n, kg C/m²). By definition C_s = C_r + C_n. Each of these components was assigned an accompanying ¹⁴C signature: ∆_r and ∆_s, which represented the measured ∆¹⁴C of the reference and switchgrass plot soils, respectively, and ∆_n, which represented the assumed ∆¹⁴C of C_n. These values were related via an isotopic mixing equation:[Image Omitted. See PDF]

This mixing relationship was used to obtain the fraction (f_n) of the C stock under switchgrass comprised by C_n and to solve for C_n:[Image Omitted. See PDF][Image Omitted. See PDF]

The ¹⁴C‐based isotopic mixing model thus provided an estimate of the C stock difference based on the observed C stock in the switchgrass plot and the shift in ¹⁴C values between the two plots.

Parameterizing Equation (3) required three ∆¹⁴C values: ∆_s, ∆_r, and ∆_n. We estimated ∆_s and ∆_r as the stock‐weighted average ∆¹⁴C values of the subsoils in the switchgrass and reference fields, respectively. In contrast, ∆_n could not be assigned a fixed value because the ∆¹⁴C of the atmosphere changes over time and there can be lags between root production and integration of root‐C into SOC. However, ∆_n could be constrained within relatively narrow range based on the known atmospheric ∆¹⁴C and plausible decomposition rates for root‐derived SOC since planting. To constrain this range, we modeled the ∆¹⁴C of SOC produced since 1998 or 2008 using a one‐pool soil C model.

The one‐pool C model was implemented in SoilR (Sierra, Müller, & Trumbore, 2012) using the function “OnepModel14” and a published atmospheric CO₂ record for northern hemisphere, extended to 2018 by assuming a 5‰ annual decrease in atmospheric ∆¹⁴C (Hua et al., 2013). The model was initiated in 1998 or 2008 with zero initial C. Inputs were fixed at an arbitrary, constant, nonzero value as the modeled ∆¹⁴C value was independent of the input rate. While a varying input rate would influence the modeled ∆¹⁴C value of the SOC, we had no basis for parametrizing a varying rate and the effect of varying inputs was small (e.g., halving litter inputs for the first 4 years reduced the final ∆¹⁴C by 4‰). The decomposition rate constant was set to two extreme scenarios: either zero (no decomposition) or ln(2) (a 1‐year half‐life). The modeled ∆¹⁴C value of the SOC pool in 2018 under each scenario was used to define a range for ∆_n. This range spanned from 0‰ to +15‰ for the Sandy Loam and Clay sites (planted in 2008), and from 0‰ to +44‰ for the Loam site (planted in 1998).

Statistical analyses

We evaluated C stock differences between the reference and switchgrass plots by propagating statistical uncertainties using Monte Carlo simulations. Simulations were used to obtain distributions for each estimate of the difference in C stocks between plots given the uncertainties in the input parameters. We obtained 95% confidence intervals from the Monte Carlo distribution of each estimate by computing quantiles of the final distributions (Buckland, 1984), and we obtained empirical p values from the Monte Carlo intervals to test the hypothesis that the difference in stocks was greater than 0. p Values were obtained using the formula p = (r + 1)/(n + 1), where r was the number of Monte Carlo replicates less than 0 and n was the total number of simulations (Davison & Hinkley, 1997). The error in each of the field‐measured properties (C stocks and isotope signatures) was modeled by generating normal distributions with the standard deviation and mean obtained from the replicate cores (Huang, 2019). To generate the normal distributions, estimated standard deviations were corrected to account for sample size by dividing them by a correction factor (c₄) which equals 0.886 when n = 3 (Huang, 2019). The distributions were assumed to vary independently. In the case of ∆_n, we assumed a uniform distribution that ranged between the limiting cases defined in Section 2.5 above. Parameter sets were drawn from the distributions 100,000 times. For each parameter set, we calculated one of two quantities: an estimate of C_n from the observed stock difference (C_s − C_r) or the ¹⁴C‐based stock difference (f_n × C_s).

RESULTS

Soil physicochemical characteristics

The three sites varied in texture, pH, and exchange properties (Table 2). Clay content and exchangeable cation concentrations were lowest at the Sandy Loam site and highest at the Clay site (Table 2). Ca was the dominant exchangeable cation at the Sandy Loam and Loam sites, whereas Mg and Ca were approximately equal contributors at the Clay site (Table 2). Soil pH values were mildly acid to mildly alkaline across three sites, and exchangeable Al concentrations were below detection, or less than 1% of the total cation pool at all sites, and thus not reported.

2 TABLESoil texture and exchange properties. Data are from three replicate cores sampled under switchgrass and paired “reference” annual crops at three sites in Oklahoma characterized by different soil textures. Values represent means of all six cores sampled at each site calculated on averages of the top three depth intervals sampled (0–30, 30–60, and 60–90 cm). Standard deviations are listed in parentheses

Site	Particle size (%)	Exchangeable cations (meq/100 g)	pH
Sand	Silt	Clay	Ca	Mg	Na	K
Sandy Loam	63 (3)	28 (3)	9 (1)	3.9 (0.5)	1.5 (0.2)	N.D.	0.1 (0.04)	6.2 (0.2)
Loam	41 (6)	37 (5)	22 (2)	8.8 (0.8)	3.1 (0.4)	N.D.	0.2 (0.03)	7.1 (0.5)
Clay	46 (10)	15 (7)	39 (13)	7.3 (1.4)	7.4 (4.0)	0.8 (1.0)	0.2 (0.02)	6.5 (0.4)

Root biomass

Root biomass values and rooting depth under switchgrass differed substantially between sites. Rooting profiles were deepest at the Sandy Loam site and comparatively shallower at the Loam and Clay sites (Figure 1). Root biomass was much greater under switchgrass at all sites (Figure 1). However, the reference plots were sampled after harvest, and the small number of cores collected (n = 3) may mean that we bypassed roots. Consequently, these differences are likely not representative of growing season conditions.

1 FIGURE. Root biomass versus depth. Data are from three replicate cores sampled under switchgrass and paired “reference” annual crops at three sites in Oklahoma characterized by different soil textures (Sandy Loam, panel a; Loam, panel b; Clay, panel c). Data from each replicate core are shown individually. Cores taken under switchgrass are shown in blue, and cores taken under the reference plot are shown in red. The soil at the Loam site (panel b) featured a buried profile, which is shown with a dashed gray line

Organic C

Total organic C concentrations were lowest throughout the soil at the Sandy Loam site, intermediate at the Loam site, and highest at the Clay site (Figure 2). At the Sandy Loam site, organic C concentrations were highest in the three cores sampled under switchgrass throughout the uppermost 200 cm of soil (Figure 2a). At the Loam site, organic C concentrations were higher in the cores sampled under switchgrass in the top 100 cm of the soil, with the largest difference in the top 30 cm (Figure 2b). We also observed a substantial “bulge” in organic C below 200 cm at the Loam site, which matched the top of the buried paleosol that we identified both in the soil series description and in our field observations. The organic C content of the buried soil was higher in the cores sampled under the reference vegetation (Figure 2b). In contrast to the Sandy Loam and Loam sites, at the Clay site organic C concentrations were generally similar under both vegetation types (Figure 2c).

2 FIGURE. Organic C concentrations versus depth. Data are from three replicate cores sampled under switchgrass and paired “reference” annual crops at three sites in Oklahoma characterized by different soil textures (Sandy Loam, panel a; Loam, panel b; Clay, panel c). Data from each replicate core are shown individually. Cores taken under switchgrass are shown in blue, and cores taken under the reference plot are shown in red. The soil at the Loam site (panel b) featured a buried profile, which is shown with a dashed gray line

Depth distribution of ¹³C

In general, the ¹³C signature of organic C varied with sampling depth across sites. At the Sandy Loam site, δ¹³C values ranged from −20‰ to −16‰ in the top 30 cm of soil, increased by 3‰–4‰ over 30–90 cm depth, and declined at greater depths (Figure 3a). This pattern appeared under both plant types, but the δ¹³C values were also approximately 2‰–3‰ higher under switchgrass (Figure 3a). At the Loam site, δ¹³C values were also depleted at the surface and comparatively higher at greater depths in a pattern similar to the Sandy Loam site (Figure 3b). The δ¹³C signature was also comparatively higher in cores taken under switchgrass, but this difference attenuated with depth (Figure 3b). At the Clay site, δ¹³C values were highest at the surface and declined with depth (Figure 3c). Patterns under the two plant covers at the Clay site were similar, with slightly higher isotopic values under switchgrass (Figure 3c).

3 FIGURE. C isotopes versus depth. Data are from three replicate cores sampled under switchgrass and paired “reference” annual crops at three sites in Oklahoma characterized by different soil textures (Sandy Loam, panels a and d; Loam, panels b and e; Clay, panels c and f). Data from each replicate core are shown individually. Cores taken under switchgrass are shown in blue, and cores taken under the reference plot are shown in red. The soil at the Loam site featured a buried profile, which is shown with a dashed gray line. The range of ∆14C for the atmosphere over the study period is shown as a gray region on the right of panels (d)–(f). C isotope data could not be collected at all depths at the Sandy Loam site because organic C concentrations were too low; data gaps are interpolated with dashed lines

Depth distribution of ¹⁴C

Radiocarbon values declined with depth at all sites (Figure 3d–f). At the Sandy Loam site, ∆¹⁴C values were near zero at the surface and declined to values near −400‰ at 150 cm. Below 30 cm, ∆¹⁴C values were systematically higher in cores taken under switchgrass (Figure 3d). At the Loam site, ∆¹⁴C values did not decline nearly as steeply as at the Sandy Loam site: at a depth of 150 cm, ∆¹⁴C was approximately 200‰. Between 30 and 90 cm, the ∆¹⁴C values of cores sampled under switchgrass were higher at the Loam site (Figure 3e). In the buried soil at the Loam site, ∆¹⁴C values were higher in cores taken under the reference vegetation (Figure 3e). At the Clay site, ∆¹⁴C values within the top 30 cm were more depleted relative to the atmosphere than at the other two sites (Figure 3f). The ∆¹⁴C values declined steeply with depth at the Clay site, reaching values in the −200‰ to −400‰ range at a depth of 1 m. At this site, the ∆¹⁴C depth profiles were broadly similar under the two vegetation types (Figure 3f).

Total organic C stocks

We obtained ESM estimates of C stocks at each site. The mean C stocks for the top 500 kg/m² of soil (approximately 0–30 cm depth) and the lower 500–1,500 kg/m² of soil (approximately 30–100 cm depth) are reported in Table 3. While we focused on ESM estimates when comparing plots to factor out bulk density differences between plots and sites, we also report total estimates to a depth of 1.2 m—which was the greatest depth at which we were able to collect samples across all sites—and to a depth of 2.4 m, which was attained at the Sandy Loam and Loam sites (Table 3). All soil chemical data and C‐isotope values are reported in Table S1.

3 TABLESoil organic carbon stock estimates. Data are from three replicate cores sampled under switchgrass and paired “reference” annual crops at three sites in Oklahoma characterized by different soil textures. Values are means, with standard deviations in parentheses. The first two columns of data represent stocks estimated on an equivalent soil mass basis; the second two columns represent stocks to a fixed depth. Stocks to 2.4 m are not shown for the Clay site because sampling to this depth was not possible there, possibly due to calcium cementation in the subsoil

Site	Plot	Organic C stock (kg/m2)
0–500 kg/m2	500–1,500 kg/m2	0–1.2 m	0–2.4 m
Clay	Switchgrass	3.9 (0.5)	2.7 (0.6)	7.4 (1.0)	n/a
Reference	3.6 (0.8)	2.6 (0.2)	7.7 (1.1)	n/a
Sandy Loam	Switchgrass	1.7 (0.1)	2.7 (0.1)	5.7 (0.3)	7.3 (0.5)
Reference	1.2 (0.1)	1.3 (0.1)	2.8 (0.3)	3.1 (0.3)
Loam	Switchgrass	4.7 (0.8)	2.8 (0.2)	8.9 (1.3)	13.3 (1.4)
Reference	2.5 (0.1)	2.2 (0.1)	5.9 (0.3)	11.6 (0.5)

We compared the C stocks under switchgrass and reference plots (Figure 4). At the Sandy Loam site, direct comparison of the C stocks suggests that there was slightly more C under switchgrass in the top 500 kg/m² of soil (stock difference = 0.4 kg C/m²; p < .01) and also in the subsoil (stock difference = 1.5 kg C/m²; p < .01). At the Loam site, we observed significantly more C under switchgrass in the top 500 kg/m² of soil (stock difference = 2.2 kg C/m²; p < .01) and in the subsoil (stock difference = 0.6 kg C/m²; p = .01). At the Clay site, the C stock difference in the top 500 kg/m² was comparatively small and not statistically significant (stock difference = 0.2 kg C/m²; p = .4) and the same was true of the subsoil (stock difference = 0.1 kg C/m²; p = .44).

4 FIGURE. Mean difference in C stock between the switchgrass and annual plant cover. Data are from three replicate cores sampled under switchgrass and paired “reference” annual crops at three sites in Oklahoma characterized by different soil textures. Blue bars show estimates for the top 500 kg/m2 of soil, gray bars show estimates for the lower 500–1,500 kg/m2 of soil, and yellow bars show estimates for the lower 500–1,500 kg/m2 of soil based on the stock in the switchgrass plot and the shift in ∆14C between plots (Equations 2–4). Error bars show 95% confidence intervals derived from Monte Carlo uncertainty propagation

Stock differences from ¹⁴C

Using the observed ∆¹⁴C values, the observed C stocks under switchgrass, and Equations (2)–(4), we developed estimates of the difference in subsoil C stocks between the plots independently of the reference plot C stock (Figure 4). Using Equation (3), we estimated that the fraction of additional C introduced after switchgrass planting (f_n) was 0.31 at the Sandy Loam site, 0.21 at the Loam site, and −0.01 (effectively zero) at the Clay site. By multiplying these values by the corresponding C stocks in the switchgrass field, we estimated that the ¹⁴C‐based stock difference at the Sandy Loam site was 0.84 kg C/m² (p < .01)—which was lower than the direct estimate derived from subtracting the observed C stocks. At the Loam site, the ¹⁴C‐based stock difference was 0.6 kg C/m² (p < .01), which overlapped closely with the direct estimate. At the Clay site, the ¹⁴C‐based estimate was near zero and not statistically significant (−0.02 kg C/m²; p = .48).

DISCUSSION

Differences in SOC

At two out of the three sites we sampled, we observed significant differences in SOC between switchgrass and reference plots in both topsoil (0–500 kg/m² soil mass, or approximately 0–30 cm) and subsoil (500–1,500 kg/m² soil mass, approximately 30–100 cm). At these two sites, differences in subsoil C were in the range of 0.6–1.5 kg C/m². This range is comparable to the value observed in subsoils at 42 paired sites where switchgrass was grown across the upper Midwest (1.2 kg C/m²; Liebig et al., 2005). Considering the cumulative soil mass of 1,500 kg/m² (approximately 0–100 cm depth), directly measured C differences normalized by the time since planting are 1.9 Mg C ha⁻¹ year⁻¹ at the Sandy Loam site, 1.4 Mg C ha⁻¹ year⁻¹ at the Loam site, and 0.3 Mg C ha⁻¹ year⁻¹ at the Clay site. These values can be interpreted as “relative changes” in that they estimate the linear rate of divergence between switchgrass and conventionally managed systems. This range of rates is typical of switchgrass systems evaluated to a comparable depth (Frank, Berdahl, Hanson, Liebig, & Johnson, 2004; Qin et al., 2016). Notably, divergence between the two land use types could represent an unknown combination of C sequestration and avoided emissions, depending on the absolute trajectory of C stocks in both fields (Sanderman & Baldock, 2010; Sanford et al., 2012). The discrepancy makes the use of paired plots for C accounting purposes complicated—but at the same time both negative emissions and avoided emissions would be benefits of bioenergy crop production in marginal lands.

Interpreting C isotopes

Both ¹³C and ¹⁴C were sensitive to land use at the three sites, and in general ¹⁴C confirmed that larger C stocks under switchgrass at these sites (or lack thereof) can be attributed to recently fixed C in the subsoil. We did, however, discover some disagreement between the directly measured C stocks and the difference estimated using ¹⁴C: the directly measured difference in subsoil C stocks was largest at the Sandy Loam, but the shift in ¹⁴C values at this site was too small to fully accommodate this difference. The simplest interpretation of this result is that the initial C stocks were greater under the switchgrass field before planting—highlighting the limits of the small sample size (n = 3) plus the spatially pseudo‐replication inherent to the paired sampling design. This interpretation is supported by texture analysis of deeper soil horizons at this site: while soil properties in the upper 90 cm of the soil profile were similar in the reference and switchgrass fields at this site, the reference plot had a higher profile‐averaged sand to silt ratio than switchgrass at depths exceeding 90 cm (mean sand/silt = 14 ± 0.7 vs. 2 ± 0.4 at a depth of 120 cm; Table S1). This indicates that soil physical characteristics did not match perfectly at this site below a certain depth. At the other two sites where the plots were more closely paired, direct and ¹⁴C‐based methods agreed.

Intriguingly, we observed less total C and comparatively depleted ¹⁴C values in the buried soil (paleosol) under switchgrass at the Loam site. The ¹⁴C values in the paleosol were less depleted under the reference plot—and were actually slightly less depleted than the overlying soil (Figure 3e). Given that roots were not observable in the paleosol, we think it is unlikely that patterns in total C and ¹⁴C at the depth are driven by modern plant cover. Instead, we think it is most likely that the soil under the switchgrass and reference plots—while similar now—experienced different histories, resulting in different C stocks and isotopes at depth. The range of ¹⁴C values that we observed in the paleosol (−394‰ to −156‰) suggest that it was buried during the mid‐ or late Holocene (i.e., in the last 5,000 years). It is possible that the paleosol under the switchgrass plot was eroded prior to burial—which would explain its lower C concentrations and ¹⁴C values relative to the reference plot. The material that was subsequently deposited over both paleosols may have been derived from upland soils containing ¹⁴C‐depleted organic matter, which could explain why the upper part of the paleosol is richer in ¹⁴C than the overlying base of the modern soil (Lombardo, Rodrigues, & Veit, 2018). More generally, deep soil sampling in paired plots can reveal inherited soil features that are not identifiable at the surface—particularly in the mid‐continental USA, where paleosols are common under alluvial and eolian deposits (Muhs, 2013).

We generally observed enrichment of ¹³C under switchgrass, particularly at the Sandy Loam site. C₄ plants like switchgrass have tissue δ¹³C values ranging from −16‰ to −11‰, whereas C₃ plant tissue ranges from −30‰ to −20‰ (O'Leary, 1988). The shifts we observed are thus consistent with an increase in the abundance of C₄‐derived C under switchgrass. Interpreting the ¹³C data quantitatively is challenging however, given that these sites have experienced a complex history that has included a mix of C₃ and C₄ crops. The recent C₃ plant contribution may explain ¹³C depth profiles at the Sandy Loam and Loam sites, where δ¹³C values in the top 30 cm of soil were lower than those in the subsoil. However, the relatively higher δ¹³C in the subsoil could also reflect fractionation during decomposition (Menichetti et al., 2015; Werth & Kuzyakov, 2010). Given these complexities, it would be challenging to use ¹³C as an unbiased tracer of switchgrass C in the context of our sites—highlighting the value of ¹⁴C.

Explaining differences between sites

The direct measurements of organic C and the isotopic calculations detected a similar trend: there was more C under switchgrass at the Sandy Loam and Loam sites, and no difference between switchgrass and reference plots at the Clay site. Multiple factors that might explain this pattern given that the three sites have different soil properties and have also experienced different management histories (e.g., tillage, and crop type in the reference fields). Furthermore, the switchgrass stand at the Loam site was 10 years older than the stands at the other two sites. Because these factors are correlated across sites, we have no way to identify which influenced SOC most strongly. Regardless, the large apparent shift at the Sandy Loam site suggests that management effects on C can be substantial even in coarse‐textured soils.

CONCLUSIONS

We found that SOC stocks were significantly larger under switchgrass than in nearby reference plots at two out of three sites in Oklahoma. SOC differences were significant at two sites with coarse‐textured soils, and not detectable at a site with fine‐textured soils. Using ¹⁴C as a tracer of belowground inputs to the subsoil after planting, we were able to confirm that differences in C stocks at the three sites were at least partly attributable to recently fixed C under switchgrass. This demonstrates that ¹⁴C can be a useful tracer for divergence of SOC stocks following shifts in cultivation or land use. Further application of ¹⁴C via repeated measurements and analysis of SOC fractions might help to constrain biogeochemical models that simulate the trajectory of C stock dynamics, improving C accounting following cultivation of perennial bioenergy crops.

ACKNOWLEDGEMENTS

This research is based upon work supported by the U.S. Department of Energy Office of Science, awards SCW1555 and SCW1632 at Lawrence Livermore National Laboratory, award DE‐SC0014079 to UC Berkeley and a subcontract to the Noble Research Institute. Lawrence Livermore National Laboratory's Lab Directed Research and Development program (#19‐ERD‐010) contributed salary support for E.W.S. Work at LLNL was conducted under the auspices of DOE Contract DE‐AC52‐07NA27344. We thank Hugh Aljoe, Kelly Craven, James Rogers, and Shawn Norton for facilitating site access and soil collection at the Nobel Research Institute farms. We thank Yanqi Wu for the opportunity to sample Oklahoma State University's field site near Stillwater, OK. Thanks also to Yuan Wang, Na Ding, Travis Simmons, Josh Barbour, Mellisa McMahon, Jennifer Black, Konstantin Chekhovskiy, Jialiang Kuang, Colin Bates, Ryan Gini, Noah Sokol, and Ricardo Feliciano‐Rivera for their help with harvests and processing of soil samples.

DATA AVAILABILITY STATEMENT

The data that support the findings of this study are available in the Supporting Information of this article.