Boosting climate change mitigation potential of perennial lignocellulosic crops grown on marginal lands

Nitrogen fertilizer (NF) is a major uncertainty surrounding the greenhouse gas (GHG) emissions of lignocellulosic biofuels. NF enhances agronomic yields and soil C inputs via plant litters, but results in soil organic carbon (SOC) decomposition, soil N2O fluxes, and a large fossil energy footprint. Thus, whether NF is beneficial or detrimental to the GHG mitigation of biofuels is unknown. Here, we show the potential GHG mitigation of fertilizing switchgrass (Panicum virgatum) at the NF rate that minimizes net GHG emissions across 7.1 million ha of marginal lands in the Midwest US, with long-term production advantages surpassing emitted GHG by 0.66 Mg CO2e ha−1 yr−1 on the aggregate. Marginal lands limited by poor N fertility could see a much greater benefit, but not SOC-rich lands, limited by low precipitation, or short growing seasons. The objectives of maximizing yield and minimizing GHG overlap only in a few environments, suggesting that maximum yield will reduce the climate benefit of cellulosic biofuels.


Introduction
Renewable lignocellulosic biofuels are one of the few mature technologies with significant potential to reduce fossil fuel consumption and greenhouse gas (GHG) emissions in the transportation sector (Gelfand et al 2020). Supplying sufficient lignocellulose feedstock to satisfy mandates for low-C biofuels will require massive land-use change, possibly involving the conversion of millions of hectares of existing cropland and marginally productive sites into perennial biomass crops (Gelfand et al 2013). Landuse change in such a scale, however, could have significant impacts on food security, biodiversity, and GHG emissions if bioenergy production directly or indirectly displaces food crops or sensitive habitats (Fargione et al 2008, Melillo et al 2009, El Akkari et al 2018. Furthermore, the potential carbon neutrality of lignocellulose-derived biofuel hinges on the ability of the production system to offset upstream and downstream supply-chain GHG emissions with the sequestration of atmospheric CO 2 in soil organic carbon (SOC) pools (Robertson et al 2011), which is a process highly dependent on the pedo-climatic conditions (Agostini et al 2015. Thus, careful selection of sites and management practices used for feedstock cultivation has become the core of strategies seeking the sustainable and efficient deployment of lignocellulose-based bioenergy systems (Somerville et al 2010, Gelfand et al 2013, Field et al 2018, O'Neill and Maravelias 2021. Nitrogen fertilizer (N F ) stands as one of the major uncertainties in the production of perennial bioenergy crops. Application of N F can increase plant growth and agronomic yields (Mclaughlin andKszos 2005, Ruan et al 2016), and the additional biomass may contribute to climate mitigation by increasing availability of biofuel feedstocks and improving efficiency of biofuel supply chains (O'Neill and Maravelias 2021). However, the responses of perennial bioenergy crops to N F additions are quite variable and poorly understood (Fike et al 2017, Monti et al 2019. Enhancing plant productivity with fertilizer additions also generally means increased C inputs to soils via plant litters (Miesel et al 2017), which may result in greater rates of SOC sequestration (Poffenbarger et al 2017). Yet, evidence suggests that in environments with low soil N availability to plants and microbes, N F additions can accelerate the rates of SOC decomposition and reduce soil C stocks, particularly in the mineral-associated fractions (Valdez et al 2017), with large uncertainty of whether N F will result in higher or lower SOC sequestration at a given candidate site. Moreover, N F also promotes microbial nitrification and denitrification processes, resulting in the release of nitrous oxide (N 2 O), a potent GHG with a global warming potential ∼265 times that of CO 2 (IPCC 2014). These emissions are on top of those incurred by the fossil fuel energy used for manufacturing, transportation, and application of the N F itself, which on the aggregate could negate the GHG mitigation benefit from any additional feedstock produced (Crutzen et al 2008).
From a climate mitigation perspective, the question of whether bioenergy crops ought to receive N F -and at what rate-needs to balance the benefits and costs of the practice in terms of GHG emissions. The objective here is to minimize the net GHG contribution of N fertilization to the life-cycle greenhouse warming intensity (GWI, Mg CO 2 e ha −1 yr −1 ) of the producing lignocellulosic biofuel. That can be accounted as the differences in GWI between the fertilized (GWI NF ) and unfertilized (GWI N0 ) systems, as show in equation (1): Solving this expression (methods) reveals that GHG Net can be quantified if we know how much additional GHG is emitted due to N F synthesis and distribution (GHG FertSyn ), field application (GHG FertApp ), the changes in soil N 2 O emissions (∆GHG N2O ), SOC sequestration (∆GHG ∆SOC ), and the negative emissions due to the additional replacement of gasoline by the biofuel (∆GHG BioFuel ): here we examine the potential GHG mitigation benefit of fertilizing perennial bioenergy crops at the N F rate that minimizes GHG Net , at the US Midwest regional scale. Our analysis focuses on the production of upland ecotypes of switchgrass (Panicum virgatum), the leading candidate crop for deployment in the region (Mclaughlin and Kszos 2005), across 7.1 million ha of non-irrigated marginal lands identified in the Midwest. Understanding the magnitude of the potential benefits at each candidate production site, as well as their context-specific drivers, can aid in the design of targeted agronomic recommendations and policy incentives for the effective use of N F in bioenergy feedstock production.

Predictive model development
Given that data on upland switchgrass productivity and ∆SOC is available only for a few experimental locations in the Midwest, we project its N F response to heterogeneous soils, climate and management using a mixed-modeling approach that combines the strengths of process-based crop simulation models (CSMs) and machine learning (ML; figure 1(a)).
To do this, we trained ML algorithms to replicate the ensemble mean of four CSMs using weather, soil, and management data summaries readily computed from available public databases. We selected the gradient boosting machine (GBM) among several algorithms examined because of its superior performance in prediction with this simulated dataset (supplementary figure 1 available online at stacks.iop.org/ERL/17/044004/mmedia). The CSM ensemble consisted of four process-based models (APSIM (Ojeda et al 2017)), DAYCENT (Field et al 2016), EPIC (Jones et al 2018) and SALUS (Martinez-Feria and Basso 2020), all of which have been previously used to simulate upland ecotypes of switchgrass in the Midwest, and represent diverse model structures. We tested each of these CSMs and their ensemble against a dataset with plant and SOC measurements that encompass three decades of experimentation in the region and a range of N F application rates (figures 1(b), (c) and supplementary table 1). Further details about the development of the predictive models for switchgrass yield and ∆SOC are included in the supplemental information.

In-silico nitrogen fertilizer rate experiments on marginal lands
We applied the developed predictive models to estimate geographical variations in switchgrass agronomic yield and SOC sequestration on Midwestern marginal lands across a gradient of N F additions. Growing bioenergy crops on marginal lands has been proposed as a strategy (Gelfand et al 2013) to avoid indirect emissions occur when biofuels production displaces agricultural production and causes additional land-use change that leads to an increase in GHG emissions elsewhere (Melillo et al 2009). Here, we consider marginal lands that could be good candidates for switchgrass establishment as those with soils under recent or current agricultural use, but with very severe limitations to row-crop production. Similar to previous studies (Gelfand et al 2013, Liu andBasso 2017), we identified these prospective sites by selecting fields that were planted with corn at least once during the period between 2008 and 2018 according to the crop frequency layer in the CropScape database (Han et al 2012). Then, we subset these soils to only those identified with land capability class (LCC) IV-VIII in the SSURGO database (Soil Survey Staff) (see supplementary table 2 for detailed   figure 3). The prediction time windows were 20 randomly selected, contiguous weather years from 1980 to 2019. Because we are interested in characterizing responses of established stands, we discarded the first 2 years of the prediction. In total this amounted to more than 63 million individual predictions.

Estimating net emissions associated with nitrogen fertilizer additions
The impact of fertilizer additions on the life-cycle emissions of biofuel production was calculated as the net difference between the total GWI under a given N F rate minus the GWI under unfertilized management (N 0 ) as shown in equation (1). The GWI of switchgrass production can be accounted as the sum of the GHG emitted by fossil fuel used for field operations (GHG Fuel ),seed production (GHG Seed ), agrochemicals other than N F (GHG AgCh ) and N F synthesis (GHG FertSyn ), N F application operations (GHG FertApp ), the biogenic direct and indirect N 2 O emissions (GHG N2O ) and changes in soil C sequestration (GHG ∆SOC ) associated with fertilizer management, and downstream emissions related to conversion of the biomass feedstock into biofuel (GHG Feedstock ): Because GHG Fuel , GHG Seed and GHG AgCh are not a function of N F , these can be ignored in the calculation of GHG Net . Additionally, GHG FertSyn and GHG FertApp equal zero under N 0 . Therefore: GHG FertSyn was assumed to be 4.5 kg CO 2 per kg of N, and GHG FertApp as 26 kg CO 2 e ha −1 . Direct and indirect N 2 O emissions were estimated following the tier 1 method in the IPCC Guidelines for National Greenhouse Gas Inventories (IPCC 2019), which include N 2 O from nitrogen fertilizer, dead biomass residues and SOC loss (supplemental table  3). GHG ∆SOC were the predicted yearly changes in SOC times 44/12 to convert to CO 2 e. Finally, for a given N R , GHG BioFuel is calculated as the product of the predicted agronomic biomass yield (Yield Ag ; Mg dm ha −1 ) and the net GHG offset (GHG offset ; Mg CO 2 e Mg −1 dm), that is, the emissions avoided from gasoline combustion: The offset was estimated as −0.737 Mg CO 2 e Mg −1 of biomass. This value assumes a 1:1 replacement of gasoline by bioethanol on an energy basis and accounts for the supply-chain emissions associated with feedstock transportation (GHG Transp ) and storage (GHG Storage ), transportation and storage, conversion in biorefinery (GHG Bioref ), distribution of biofuel (GHG Dist ), and the CH 4 and N 2 O emissions produced by biofuel combustion (GHG Combust ): where L T and L S (unitless) are the portion of biomass lost during transportation and storage, respectively; Yield BioFuel is the yield from biomass-to-ethanol conversion, NCV BF is the net caloric value of bioethanol and GWI Gasoline is the life-cycle GHG emissions intensity of gasoline. Values for each of the parameters in equation (6)   Projected across all the marginal lands identified in the Midwest and 20 years of weather (figure 2), site-mean agronomic yields without N F are expected to range between 2.1 and 11.9 Mg ha −1 , averaging at 5.6 Mg ha −1 (figure 2(a)). Our predictions, however, show substantial spatial variability due to differences in soil characteristics and climates (figure 2(b)). Highest yields of unfertilized switchgrass are predicted mainly in the upper Midwest, eastern areas of the Northern plains, and scattered pockets of Indiana and Ohio. Unfertilized yields are projected to be relatively low in the drier parts of the central and northern plains (west of the 99th meridian) and western Illinois. Further, most marginal land cultivated with upland switchgrass (∼95%) is expected to accumulate SOC under unfertilized conditions, with ∆SOC projected across all marginal land averaging 0.58 Mg C ha −1 yr −1 (figure 2(d)). Negative ∆SOC under unfertilized conditions is predicted in some areas in North Dakota and South Dakota, and isolated areas in Minnesota, Wisconsin, Ohio and Michigan (figure 2(e)).
Analysis of emulator predictions reveals that N F rate is among the most important determinants of agronomic yields, ranking 3rd most important among all 40 variables included, and accounting for 17% of the total variance (supplementary figure 6). With annual applied N F (up to 200 kg N ha −1 ) upland switchgrass is modeled to respond positively in most marginal sites following diminishing-returns pattern. With full N F , the projected long-term mean yield boost averages 3.1 Mg ha −1 (or 54%), although yield increases are generally predicted to be small beyond 100 Kg N ha −1 rate (12% on average, figure 2(c)). This result is supported by findings from literature syntheses that suggests that 100-180 kg N ha −1 are Conversely, ∆SOC is projected to have a decreasing trend with increasing N F rate. This is consistent with the theory that N F in switchgrass has an overall dampening effect on the sequestration potential of in soil C pools via enhanced SOC decomposition (Valdez et al 2017). Yet, the mean effect is minimal, decreasing from 0.58 Mg C ha −1 yr −1 with no N F to 0.46 Mg C ha −1 yr −1 with full N F (figure 2(f)). In fact, N F rate seems to have much less influence in the emulator for predicting ∆SOC, ranking 9th in variable importance and accounting for less than 2% of the variance (supplementary figure 6). Nevertheless, both predicted responses to N F exhibit substantial spatial variation, explained in part by fundamental differences in soil characteristics and climates. For example, the yield response to N F is expected to be much more marked in sites with high precipitation, coarser texture and low SOC, whereas the negative ∆SOC response is more pronounced with increasing stand age, in sites with low precipitation and in soils with mesic levels of SOC (∼1%-2% C; supplementary figure 7).

Optimized nitrogen fertilizer rate based on minimizing net GHG emissions
Next, we contextualize these projected site-specific N F responses in terms of how they could impact the GHG emissions of the biofuel production system as defined in equation (2). In this calculation, we include the projected changes in yields and SOC sequestration, and the biogenic N 2 O emissions due to N F application using the IPCC tier 1 methodology (IPCC 2019). In addition, we apply a life cycle assessment (LCA) methodology to account for the upstream emissions due to fertilizer synthesis, transportation and application, and changes in the downstream emissions related to the manufacturing of the biofuel, its combustion, and the value replacement of gasoline fuel (methods).
A graphical example of the calculated net GHG balance as a function of N F rates is shown on figure 3. Here we display the mean GHG balances across all marginal land sites. In figure 3(a) we show both GHG FertSyn and GHG N2O increasing monotonically with NF applied, whereas GHG BioFuel (i.e. the emissions offset credit due to the additional biofuel produced) becomes more negative mirroring the yield response (figure 2(c)). Furthermore, GHG ∆SOC is counted as net positive emissions, given the small dampening effect of N F on the SOC sequestration rate predicted by the model (figure 2(f)). Including all GHG sources and credits indicates net GHG emissions mitigation (i.e. negative GHG Net ) for a wide range of N F levels up to 140 kg N ha −1 .
Focusing solely on the average response, however, would masks important variation. Indeed, further examination of a gradient of SOC content levels reveals that yield increases in sites with greater than 2% SOC in the top 25 cm are, on average, predicted to be insufficient to overcome the additional life-cycle emissions incurred by applying fertilizer (figure 3(b)). Thus, from a GHG mitigation perspective, these C-rich soils should not generally receive N F . Yet, these sites represent a small fraction of marginal lands in our dataset (<5%). Optimal N rates are projected to increase monotonically with decreasing topsoil SOC content, up to ∼120 kg N ha −1 in soils with very low topsoil SOC content ( figure 3(b)). . Production outcomes at the projected optimal N fertilizer rate. Maps show (a) the geographical spatial distribution of the predicted optimal NF rate, (b) the expected probability of that the optimum NF rate is greater than zero, and changes in (c) agronomic yields and (d) soil C sequestration rate at the optimum NF, compared to the unfertilized baseline. Pixels (0.01 • spatial resolution) show the area-weighted average of the predictions. Across all marginal lands, the modeled optimum N F rate averages at 74.8 kg N ha −1 , although its distribution is quite wide and spatially variable as shown in figure 4(a). Areas with low optimal rates (<50 kg N ha −1 ) are mainly projected in western states, and large areas in Minnesota, southern Iowa and northern Missouri, and isolated areas in the eastern states. These low optimal rates occur mostly in sites with low precipitation, SOC-rich soils, or a combination of both as is the case of areas of North and South Dakota. Indeed, these conditions in many marginal lands in the Great Plains, and to a lesser extend areas of Missouri and northern Wisconsin, lead to environments where the optimum N F rate is much less likely to be greater than zero ( figure 4(b)). On the other hand, high optimal N F rate (>100 kg N ha −1 ) is projected mainly for the eastern states, particularly in the southern tip of Indiana, Illinois and Missouri ( figure 4(a)).
With optimal N F rate, agronomic yields across all marginal lands are expected to increase by an average of 1.99 Mg ha −1 (figure 4(c)), which aggregated at the regional scale would amount to 14.1 million tons of additional lignocellulosic feedstock. Accounting for transportation and storage biomass losses and assuming a conversion yield of 355 l per Mg of biomass at the biorefinery (Argonne National Laboratory 2020) (supplementary table 3), the boost in productivity would amount to nearly 4.6 million m 3 of additional bioethanol. This represents about 12% of the total statutory cellulosic biofuel production goal of 39 million m 3 for 2020 in the US (National Research Council 2012). The optimal N F is also projected to have a more mixed effect on SOC sequestration, with only a few marginal lands seeing increase SOC sequestration, particularly in the Great Lakes states (Michigan, Minnesota, Wisconsin), areas of Iowa, and the southern tip of Indiana, Illinois and Missouri. Overall, however, the effect for the entire region is projected to be neutral, with a mean change in soil C sequestration rate of 0.006 Mg C ha −1 yr −1 .
Given these responses, optimum N F additions would mitigate emissions by 0.66 Mg CO 2 e ha −1 yr −1 on average across all marginal lands ( figure 5(a)). Areas surrounding the Great Lakes are projected to have the greatest N F -related GHG mitigation potential (up to 2.8 Mg CO 2 e ha −1 yr −1 ), likely a result of concomitant gains in SOC sequestration and a large yield response to N F (figure 4). On the aggregate, optimal NF additions could provide up to  (1) at the site-specific optimal N rate (i.e. the rate that minimizes net GHG emissions), (2) using a blanket rate with 74 kg ha −1 yr −1 , and (3) at the rate that maximizes yield, shown per hectare (a) and total area (b). Maps in (c) show the geographical distribution of the spatially weighted average net GHG resulting from fertilizing using each scenario. 4.63 Tg CO 2 e yr −1 of GHG mitigation compared to unfertilized switchgrass production across all marginal lands examined.
Overall, these results reveal the significant contributions that optimal N F could provide to the GHG mitigation of future lignocellulosic biofuel production in the Midwest. However, these projections should be considered as a hypothetical ceiling given that on-the-ground implementation of site-specific optimal N F management could face important challenges due to the uncertainty in production costs (Fike et al 2017). Therefore, we further examine how optimum N F would compare to a more simplistic blanket recommendation of 75 kg N ha −1 . This value not only represents the average optimum N F rate across all marginal lands, but also approximates the long-term average N removal of the harvested biomass (Monti et al 2019).
The blanket-rate strategy is projected to result in a net GHG mitigation of 0.35 CO 2 e ha −1 yr −1 on average, or 2.48 Tg CO 2 e yr −1 on the aggregate, when compared with the unfertilized scenario (figures 5(a) and (b)). This GHG mitigation is about 54% of that projected to be achievable with the optimal rate. Both of these scenarios compare favorably to a situation where N F is applied at the rate that maximizes longterm yields (∼160-180 kg N ha −1 ; supplementary figure 8 and Monti et al 2019), resulting in a slight increase in the GHG emissions over the unfertilized scenario (0.036 CO 2 e ha −1 yr −1 on average or 0.26 Tg CO 2 e yr −1 aggregated; figures 5(a) and (b)). It should be noted that this latter scenario exhibits the strongest spatial divergence, with areas of large positive net GHG in the Great Plains states and southern Iowa, and large negative net GHG in the Great Lakes states and areas of Illinois and Indiana ( figure 5(b)).

Discussion and conclusion
Our findings seem to challenge the prevalent belief that N F application is generally detrimental to bioenergy GHG balances and that its use should be avoided in the production of lignocellulosic feedstock (Erisman et al 2009, Smith et al 2010. On the contrary, here projections suggest that the long-term yield advantages of relatively low amounts of N F in marginal lands would surpass the GHG emissions stemming from its use. Even a blanket recommendation of 75 kg N ha −1 that would forego nearly half of the potential GHG mitigation benefit, is predicted to be a net positive on the aggregate for Midwestern marginal lands ( figure 5(c)). The empirical evidence, though limited, seems to concur. Indeed, a field study in Michigan (Ruan et al 2016) that examined field-level GHG measurements found that N F additions could support GHG mitigation of switchgrass production up to 140 kg N ha −1 , although GHG Net was minimized with a N F application of 56 kg N ha −1 . Importantly, our projections seem to imply that the benefit could be much greater in marginal lands that are limited by poor N fertility as is the case of much of Michigan, Wisconsin, and isolated areas of the eastern Corn Belt ( figure 4). Meanwhile, areas with SOC-rich soils (Iowa, Minnesota) and/or limited by low precipitation (Nebraska, Kansas) or short growing seasons (North and South Dakota) are predicted to be less likely to respond to N F . Thus, the potential GHG mitigation benefit in these marginal lands is minimal. Lastly, it should be noted that the objectives of maximizing yield and minimizing GHG Net do not seem to overlap in all but a few environments (figure 5(c)). As such, widespread implementation of fertilizing strategies aiming to maximize yield will likely result in the overall reduction of the climate benefit of cellulosic biofuels (figures 5(b) and (c)). (We include a map classifying Midwestern marginal lands in four levels of optimal N rates in the supplementary figure 9 to facilitate discussions with stakeholders).
Of course, it is likely that the optimal N F rates estimated here, which are based on minimizing GHG emissions, are different from those N F rates that are most economically efficient (i.e. where N F marginal revenue equals marginal cost). Indeed, a recent study (Fike et al 2017) that examined production sites in five states concluded that typical N F recommendations (56-112 kg N ha −1 ) added between $37 and $74 ha −1 to production costs and required switchgrass biomass prices in excess of $70 Mg −1 to be economically justifiable under the observed yield responses. Whether the bioenergy industry or ancillary carbon markets can support such pricing is yet to be determined and lies beyond the scope of our analysis.
We advocate caution when interpreting these results given the various sources of uncertainty contained in our projections. Foremost are the uncertainties associated with our characterization of the responses to N F additions across the pedo-climatic conditions of Midwestern marginal lands. Although the process-based CSMs here used are reasonably capable of capturing the variation in biomass production and SOC changes and their dynamic interaction with agricultural management (supplementary figure 4 (Surendran Nair et al 2012, Field et al 2016, Jiang et al 2017, Liu and Basso 2017, Ojeda et al 2017, Jones et al 2018, Jarecki et al 2020, Martinez-Feria and Basso 2020), their complex structure, large sets of parameters, and statistical misspecification (Wallach 2011) generally means that prediction uncertainty is large when CSMs are applied outside the calibrated domain. These uncertainties are partially addressed by using the ensemble mean of independent CSMs, which has been shown to provide more accurate predictions outside of sample populations (Rosenzweig et al 2013, Bassu et al 2014, Martre et al 2015, Basso et al 2018. The emulation step unlikely added significant uncertainty at this scale because of the independence between CSM and ML errors, which would tend to cancel out when aggregating across soil types and geographical locations. Nevertheless, there are a few known weaknesses shared by all CSMs that could bias our projections. Chief among them is the representation of associative N fixation (ANF) by free-living diazotrophs, which is increasingly recognized as a potentially important source of N to switchgrass (40-60 kg ha −1 yr −1 ), especially in N-limiting environments (Roley et al 2018, Bahulikar et al 2021. Yet, in contrast to biological N fixation (BNF) in legumes, ANF is largely episodic and irregularly detectable (Roley et al 2019), and our current understanding of the environmental factors controlling ANF rates is rather incomplete (Bahulikar et al 2021). Thus, current switchgrass CSMs, either do not account for this pathway (EPIC) or use instead BNF routines as a surrogate (APSIM, SALUS). Only the DAYCENT model has an explicit routine that, in the absence of soil phosphorus data, calculates free-living N fixation rates as a linear function of annual precipitation (Necpálová et al 2015). We are not aware of any work verifying the robustness of these routines for capturing ANF environmental variation in switchgrass, or whether it is even possible at this stage given the paucity of data available. The likely effect of underestimating non-N F sources under certain environments would be the overestimation of the yield response, which would lead to an overestimation of the optimal N rate. Although we found no widespread evidence of overestimation of the sitespecific yield response in the experiments (supplementary figure 10), a conservative interpretation of our data may still be warranted until further examination of this mechanism within the models can be performed.
An additional source of potential bias in our optimal N F rate estimates is the use of the offline calculation of direct and indirect soil N 2 O emissions following the IPCC tier 1 methodology. We opted not to use simulated N 2 O because not all the CSMs we selected contain an explicit nitrification-denitrification routine, or have not been comprehensively validated in production of perennial bioenergy grasses. On the other hand, the uncertainties of the IPCC methodologies are more thoroughly understood. Namely, that the response of annual N 2 O fluxes to N F applied tends to be exponential rather than linear (Philibert et al 2012). Thus, the use of a constant emission factor as required by IPCC tier 1, tends to overestimate N 2 O fluxes at lower N F rates, and underestimate them at higher rates. This pattern seems to hold in switchgrass production (Ruan et al 2016). Using a non-constant emissions factor that increases with N F rates to calculate N 2 O fluxes (as shown by Ruan et al 2016) results in overall lower optimum N F rates with a much narrower distribution (supplementary figure 11). Nevertheless, it is not known how well this exponential relationship, which is derived from a few experiments, would represent emissions at the regional scale. Thus, more work is needed to improve the characterization of soil N 2 O emissions when estimating GHG balances of bioenergy.
A final, possibly large source of uncertainty are the assumptions made for calculating the biofuel GHG offset, which essentially scales the projected yield and ∆SOC responses into GHG emissions (supplementary table 3). Here, we used a purely static LCA approach, which assumes constant values derived from the GREET model (methods) notwithstanding the spatial distribution of future infrastructure (e.g. depots, biorefineries, etc) or biomass end use (e.g. pelletization, direct combustion, energy cogeneration). These would result in changes in the efficiency of utilization of the feedstock supply (O'Neill and Maravelias 2021), particularly because lower switchgrass yields would require larger collection radius, and therefore greater GHG emissions of switchgrass transport (Gelfand et al 2013). Yet, the uncertainty associated with these parameters may be substantial. A sensitivity analysis where each of these values are randomly varied within ±20% of the original value reveals a potentially large departure from the mean GHG response (supplementary figure 12), though the computed optima tend to cluster around the mean found in our main analysis (75 kg N ha −1 ) suggesting low bias. Yet, two parameters the yield from biomass-to-ethanol conversion (Yield BioFuel ) and the life-cycle GHG emissions intensity of gasoline (GWI Gasoline ), resulted to have the greatest impact in the calculation of the biofuel offset (supplementary figure 10). Thus, much of the uncertainty can be reduced by finding best values for those parameters.
Despite the uncertainties, our study represents an important step forward, as it is the first to examine prospective responses of switchgrass at such scale and resolution, and the emulator model developed could provide the groundwork for building dynamic tools to guide future decisions in agronomic management, technoeconomic analysis, and policy formulation (Ruan et al 2016). Only through this kind of high-resolution predictive tools, biofuel supply-chain network design can move toward integrating temporally and spatially explicit data to reflect the reality of heterogeneity in biomass productivity and SOC sequestration capacity across the landscape, and to allow consideration of land allocation and crop management decisions in addition to infrastructure and operational decisions (O'Neill and Maravelias 2021). This is becoming more critical given the interest in targeted integration of perennial bioenergy crops within the margins of row-crop fields or other marginal subfield areas (Brandes et al 2018, Basso and Antle 2020), which could not only avoid GHG through the production of biofuels, but also by replacing low yielding, highly polluting croplands (Basso et al 2019). Therefore, further research is needed that comprehensively evaluates the modeling uncertainties here outlined to advance prospects for efficient site-specific N F management in perennial bioenergy crops, and to maximize the GHG mitigation potential of this technology.

Data availability statements
The data that support the findings are within this study in the cited references. R code and the database needed to run the developed GBM emulators has been made public at the following repository: https://doi.org/10.5281/zenodo.5116736.
The data that support the findings of this study are available upon reasonable request from the authors. B Basso  https://orcid.org/0000-0003-2090-4616