The Hidden Costs of Land Degradation in US Maize Agriculture

The United States is a world leader in the production of maize and other crops and the agricultural success of the country is directly linked to the intensive use of fertilizers and irrigation. However, even in advanced agricultural systems, soils can become degraded over time due to factors such as soil organic matter (SOM) loss and erosion. Here, we use a series of scenario‐based model analyses to show that about one‐third of current annual US. N fertilizer use in maize agriculture is used to compensate for the long‐term loss of soil fertility through erosion and organic matter loss. This leads to over a half billion dollars per year in extra fertilizer supply costs to US farmers. These results highlight the potential to reduce both the input costs and environmental impacts of agriculture through the restoration of SOM in agricultural soils.

matter over time. These changes can occur through management decisions involving actions such as tillage, residue management, and long-term organic amendments (Diacono & Montemurro, 2011).
The United States is one of the world's largest users of fertilizer and applies more N and P per ha than other high yield agricultural countries including most of the European Union members (albeit less than China and some other rapidly developing economies, FAO, 2020). High rates of synthetic fertilizer offset the removal of crop material from fields but such additions also serve to offset a second category of nutrient losses; specifically, they replace nutrients lost to erosion, SOM decline, and other forms of soil fertility loss. These "compensatory inputs" of nutrients offset the impacts of declining soil fertility but are difficult to quantify because nutrient loss takes place over many decades alongside major changes in agricultural practices and technology. Although difficult to measure, the compensatory additions of fertilizer may have significant economic and environmental impacts.
In this study we used the Environmental Policy Integrated Climate (EPIC) model in a series of 100-year scenarios examining the potential impacts and costs of soil degradation in maize agriculture in the United States. In these simulations, described in more detail in the methods, we focus on the underlying changes in soil properties that accompany land degradation, the resulting declines in soil fertility and the subsequent impacts on yields under different scenarios of fertilizer and irrigation use. To accomplish this, we simulate 100 years of US maize cropping assuming a current (circa 2010) distribution of maize agriculture with broad adoption of conservation tillage practices as discussed below. We carry out four 100-year scenarios of agriculture during which crop growth potential, pest management, climate, and harvest efficiency are all held constant at modern levels. This assumption of current conditions allows us to separate out the effects of genetics and management activities from the underlying effects of degradation. This type of scenario is clearly not intended to represent actual conditions (where all these factors change simultaneously) but allows insight into a key factor (degradation) that would otherwise be hidden under other changes. However, in each of these scenarios, we begin with SOM C, N, and P that is, similar to the conditions present at the start of agricultural cropping and use the combination of scenarios to examine how soil degradation, including organic matter loss, influences the use of fertilizers over time.

Modeling Scenarios
We created four different management scenarios to estimate soil degradation impacts on US maize agriculture. The "no inputs" scenario included 100 years of management without fertilization or irrigation. Although this scenario is not realistic in the US and other highly developed agricultural systems, it is a common occurrence in developing economies where fertilizers are not available or are too expensive (Chianu et al., 2012). This scenario also allows examination of the cumulative impact of soil degradation on yields; a situation that is, difficult to examine under high rates of fertilization. The second scenario is a "fertilization" scenario in which nutrients were added in the model in response to plant demand. The third scenario-"irrigation"-provides water sufficient to meet plant needs (but no nutrient additions). Lastly, the "full replacement" scenario included management with fertilization and irrigation to meet plant demand and is the most similar to contemporary US management practices.

Configuration of the Agroecosystem Modeling
A process-based agroecosystem model, EPIC (J. R. Williams et al., 1989), was used to estimate crop productivity as affected by erosion and nutrient mining (J. R. Williams et al., 1984). The EPIC model represents soil-crop-atmosphere interrelationships with various management scenarios. Major processes in the model included plant growth, development and production, nutrient cycling, nutrient loss, and plant/soil management practices. The model was used to create several scenarios of agriculture that compare modern practices carried out over 100 years to practices without some or all fertilizer and irrigation inputs.
We used contemporary management, crop yield potential, irrigation, and fertilization to generate the four 100-year scenarios described above. To approximate a realistic climate while removing the effects of recent climate change, we used a 31-year observed US climate record and "looped" this climate record 3.2X for a 100-year simulation. Combined spatial and temporal data on US management practices, including tillage, is unavailable so we assumed conditions that approximate conservation tillage practices as a conservative (e.g., not worse case) approach to the simulation however we also carried out a focused uncertainty analysis evaluating key model inputs and parameters including tillage, the harvest index, CO 2 sensitivity and include these results below and in the supplementary materials (Table S3 to S5).
For all simulations, we initialized the model with soil C, N, and P at equilibrated levels consistent with what they typically would be when agricultural activities are started on a new piece of land (Table 1, Izaurralde et al., 2012). Soil carbon values at the beginning of the time-series were based on the US STATSGO soil carbon data base (Table 1). Because there are no distributed soil carbon values for actual preagriculture conditions in US, we initialize with the carbon content common in native vegetation in our simulation areas. This is as close an approximation of preagricultural conditions as is possible with the information available. We used HWSD to initialize soil properties derived from the soils that are not under management and fertilization, irrigation, and conservation tillage for the beginning of the simulation period (FAO Soil Portal, 2019;FAO et al., 2009). To ensure accurate simulation values, we validated a range of input and output variables against published values. These include SOM content, crop yield and erosion rates as these are key determinants of agronomic model accuracy (P. Smith et al., 1997). By the end of the simulation, soil C, N, and P levels are consistent with contemporary agricultural system values although we did not perform a formal validation given the difference between how these simulations are structured versus the actual history of agriculture in the US.
Each simulation unit in the US was created by the aggregation of climate, soil, land cover, and slope. A simulation unit was defined as any location in the US where maize is currently cultivated at any density (USDA NASS, 2019). Accordingly, the map includes locations where maize cultivation dominates the agricultural landscape (USDA NASS, 2019). EPIC simulations were performed for multiple scenarios in each individual simulation unit. There were a total of 11,206 simulation units in this study. Mean yields presented in this paper represent the average of output in these simulation units and not a weighted average of actual production values across the US. This is the primary reason by the average US maize yield in this estimate is higher than the actual reported average for the country.

Yield and Limitation Calculations
The EPIC model tracks a variety of environmental limitations during a simulation and adjusts crop growth estimates downward based on these limitations (J. R. Williams et al., 1989). The resulting potential growth calculation provides a maximum biomass accumulation which is then reduced based on constraining factors that are most important for growth regulation including water (water stress) and nutrients deficiency (nutrient stress), nonoptimal temperatures (temperature stress), and waterlogging (aeration stress) which are used to reduce the potential growth. Final actual yield is then calculated based on the plant biomass and harvest index related to the grain yield as a ratio of the above-ground biomass (J. R. Williams et al., 1989, Equations 1-5).
Maize production (Mg ha -1 ) where YLD is the amount of crop removed from the field (Mg ha −1 day −1 ), HI is the harvest index, STL is the standing biomass of the crop (Mg ha −1 day −1 ) and the STL accumulated up to the current day of the simulation, STL i−1 is the STL present the day before the current day of the simulation, BIOM' is the actual biomass (Mg ha −1 day −1 ) that can be accumulated during the current day of the simulation considering REG which is the crop growth regulating factor such as water stress, temperature stress, nutrient stress, and aeration stress, BIOM is the potential biomass (Mg ha −1 day −1 ), PAR is the intercepted photosynthetic active radiation (MJ m 2−1 day −1 ), RUE is the radiation use efficiency factor for converting energy to biomass (kg ha −1 )/(MJ m 2−1 ), and WAVP is a crop parameter relating RUE and VPD (vapor pressure deficit [kPa]).
In EPIC, growth constraints are represented as stress factors and reported in days per year where the most limiting stress on a particular day is the primary constraint on growth. These regulating factors include water stress (day yr −1 ), nitrogen stress (day yr −1 ), phosphorus stress (day yr −1 ), and nutrient stress (nitrogen stress and phosphorus stress, day yr −1 , Equations 1-11) described in the supplementary materials).

Nutrient Dynamics
The simulation of nutrients (nitrogen, N and phosphorous, P) cycling and losses is one of the main components of the EPIC model. The organic N cycling and transformations are simulated by coupling the dynamics of organic N and carbon (C) in soil as described by Izaurralde et al. (2006). This approach follows the concept used in the CENTURY model (W. J. Parton et al., 1987Parton et al., , 1993Parton et al., , 1994 where C and N contained in the SOM are split into pools characterized by different turnover times ranging from days to hundreds of years (Gassman et al., 2010;. Nitrogen losses simulated by the EPIC model include leaching, surface runoff, and lateral subsurface flow. The organic N transported by sediment during individual runoff events is simulated using a loading function developed by McElroy et al. (1976) and modified by J. R. Williams (1990) and J. R. Williams and Hann (1978). Denitrification is simulated as function of temperature and soil water content while nitrification is based on the first-order kinetic rate proposed by Reddy et al. (1979). Volatilization is estimated along with nitrification and, for surface-applied ammonia fertilizers, is regulated by temperature and wind speed (J. W. Williams et al., 2012). Other processes such as crop N uptake, fertilization, N fixation, and N input from rainfall are considered and simulated in the EPIC model. Simulation of P includes inorganic P dynamics, mineralization of organic P, immobilization, and losses. The inorganic P simulations follow the approach developed by C. A. Jones et al. (1984) with P transferred between the labile, active, and stable mineral pools. The flow between these pools is regulated by several factors such as temperature, soil water content, P sorption coefficient, and the amount of material in each pool. The mineralization and immobilization of organic P is simulated following the approach developed by C. A. Jones et al. (1984) based on a fresh organic pool and a stable organic pool with turnover governed by C:N and C:P ratios. Simulation of P losses include soluble P in surface runoff and later flow, and P transported by sediment (J. R. Williams, 1990). More detailed information is available in J. R. Williams (1990) and J. W. Williams et al. (2012). All the processes described above interact with other simulated processes and activities like tillage operations, fertilizations, and crop growth.

Erosion and Tillage
We estimated agricultural soil erosion in the US maize agriculture using the modified universal soil loss equation (MUSLE) approach (J. R. Williams, 1975;J. Williams & Berndt, 1977, Equations 6 and 7). This is one of the eight different options available in the EPIC model to estimate the soil erosion and the losses of N and P attached to the soil eroded. In this way, it was possible to analyze the impact of these losses on the US maize production. MUSLE is a physically based model that uses a relatively simple set of inputs including basic soil properties, cover (informed by land management), topography and climate. This spatially explicit estimate is used to calculate erosion N and P loss in the US maize agriculture (Equations 18-25 described in the supplementary materials). Because the EPIC model was designed to assess the effect of soil erosion on soil productivity, reliable estimation of runoff was a key factor since the early development stages. The runoff estimation is based on the SCS curve number approach and in early test ) showed that the model was able to replicate realistically mean surface runoff and sediment measured in three small watersheds located in Texas. Later, the EPIC model was successfully used to simulate surface runoff in many studies (Steiner et al., 1987;Puurveen et al., 1977;Chung et al., 1999;Wang et al., 2012).
Sediment yield (Mg ha -1 yr -1 ) MUSLE is calculated with the equations: where Y is the sediment yield in Mg ha −1 , K is the soil erodibility factor, CE is the crop management factor, PE is the erosion control practice factor, LS is the slope length and steepness factor, ROKF is the coarse fragment factor, Q is the runoff volume in m 3 , q p is the peak runoff rate in m 3 /s. PE is determined by considering the conservation practices and the other equations (i.e., K, CE, LS, and ROKF) are described in the Equations 12-17 in the supplementary materials (J. R. Williams, 1975;Wischmeier & Smith, 1978). All variables (e.g., K, CE, PE, LS, and ROKF) are calculated by Arekhi et al. (2012) and Sadeghi et al. (2014).
Tillage systems in the United States have changed considerably in recent years. For most of the 20th century, most tillage followed conventional practices of plowing and harrowing and incorporation/removal of plant residue. Since the beginning of 1960s the use of conservation tillage has increased (Matson et al., 1997). From 1989 to 2008 in the US, there was a trend of decreasing use of conventional tillage and increased adoption of conservation tillage and no-till farming (S. Smith et al., 2014). During the years (1989)(1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008) conservation tillage was 57%-68% (average 62%) of total land area during the period of 1989-2008 (S. Smith et al., 2014). From 2008 to the present, there was a small reversal of this trend due to herbicide (e.g., glyphosate)-resistant weeds requiring tillage control, as well as the increase in organic farming (Davis et al., 2007;Foresman & Glasgow, 2008;Scott & Vangessel, 2007). Largely due to reduced tillage, soil erosion in the US declined over the 30 years prior to 2008 (Marsh et al., 2006;Xia & Shimmin, 2015). Management practices can play a role in conservation tillage and no-till farming systems and will reduce soil erosion (Mitchell, 2011). However, these specific management interventions are beyond the scope of this study as there is limited information on the temporal and spatial implementation of these interventions.
In this simulation, we take a relatively conservative approach to erosion estimates and use conservation tillage systems for all four simulations on all land throughout the entire period. Conservation tillage is a method of soil cultivation which leaves the previous year's crop residue on the soil surface (Jat et al., 2010). Crop residue cover for conservation tillage typically ranges from 30% to 60% after harvest (Hively et al., 2018) and crop residue cover in our simulation was 46%. Typical conservation tillage depth ranges from 15 to 20 cm (Le et al., 2018); the simulation presented here uses a 20 cm tillage depth. This approach will yield lower erosion estimates than would occur from conventional tillage and so the results presented in this paper are a conservative estimate of the cumulative erosional loss over the full 100-year simulation period.
Soil N loss was calculated based on a soil N balance calculation considering the input of N by fertilization, changes in soil N concentration over time, and N loss through erosion, leaching, denitrification, and other export pathways. Soil N loss is a measure of the cumulative change in soil N due to all inputs and outputs over time. A similar calculation was not possible for P due to the large stocks of inorganic P in soils and the limited information on organic P changes in EPIC.

Data Information
For this study, we used a suite of data for the biophysical process simulations. This information is listed in The Harmonized World Soil Database v1.2 (HWSD, FAO Soil Portal, 2019) was used for soil physical information. We used HWSD so that the framework could eventually be scaled globally. HWSD data include soil properties averaged across two depths (0-30 and 30-100 cm). For grid-based biophysical process simu-lations, soil unit composition of each grid contained soil hydrologic group (e.g., A, B, C, and D), depth (m), bulk density (Mg/m 3 ), sand content (%), silt content (%), soil pH, sum of bases (cmol/kg), calcium carbonate content (%), cation exchange capacity (cmol/kg), and coarse fragment (rock) content (%) ( Table 1).
For the topographic information, a digital elevation model (DEM, 1 km resolution) was obtained from the Global Land One-kilometer Base Elevation (GLOBE) provided by National Oceanic and Atmospheric Administration (NOAA). Slopes (0.01-0.23 m/m) were calculated from the DEM (Table 1).
To create the maize map, land cover data was retrieved from the Cropland Data Layer (CDL, 30 m resolution) in United States Department of Agriculture (USDA)-National Agricultural Statistics Service (NASS, Jang et al., 2019) for the United States (Table 1). The land cover map was overlaid on the prepopulated US grid cells aggregated by climate, soil, and slope for the whole U.S and only active grid cells including at least one 30 m pixel classified as maize within a simulation unit were extracted to create a simulation unit map for EPIC. Average simulation unit grid cell size is 18 km resolution. Due to the different grid resolution size between the simulation unit and land cover map, maize areas in the spatial maps created in this study cover larger areas than the land cover map USDA. For this reason, all simulated yields in this manuscript are reported on an average per unit area basis rather than as a spatially weighted national average.
Crop management information included planting and harvest dates and potential heat unit (PHU). Planting and harvest dates are retrieved from the Center for Sustainability and the Global Environment (SAGE) at the University of Wisconsin (Sacks et al., 2010). PHU was calculated based on planting and harvest dates and long-term AgMERRA climate data (years 1980-2010, Table 1).

Model Testing
The starting point for all four scenarios was set to approximate modern maize yield potential in the United States combined with soils that represent preagricultural conditions with respect to SOM content as described below. As the simulation progressed for 100 years, we compared simulated yield for the full replacement scenario to published data on SOC, yield, and erosion. Mean (standard deviation in square brackets) annual simulated SOC (15.56 [±6.19] Mg ha −1 yr −1 in 0-20 cm depth) at the end of the simulations for the whole US maize agriculture fell within the range of the reported mean annual SOC with some variation between the different simulations ( Figure S5,     was also compared to the reported mean annual yield (7.9 [±1.3] Mg ha −1 yr −1 for 1980-2010, USDA NASS, 2019) and the range of mean annual simulated yield (7.9-8.5 Mg ha −1 yr −1 ) fell within the range of the reported mean annual yield (6.6-9.2 Mg ha −1 yr −1 ). A more comprehensive uncertainty analysis is included in the supplementary materials (Table S3 to S5).
Consistent with Quinton et al., 2010 and because EPIC does not report a dynamic N value in eroded soil, we assume soil N content from erosion is 0.1% and a range of erosion N loss estimates for the scenarios (scenario averages varied from 4.5 to 11 kg N ha −1 yr −1 ) ( Table 2). According to the literature (Batjes, 2014;Quinton et al., 2010;Smil, 1999;Tian et al., 2010), the estimated erosion N loss ranges between 18 and 65 Tg yr −1 based on global erosion range (25-45 Pg yr −1 ) from agricultural lands. A range of erosion N loss is estimated based on published erosion rates and N content. In order to compare the estimated erosion N loss from this study with the literature, we used the percentage range of N in eroded sediment (0.07%-0.14%) from the literature (Batjes, 2014;Quinton et al., 2010;Smil, 1999;Tian et al., 2010). The resulting estimated US erosion N loss ranged between 5.1 and 10.3 kg N ha −1 yr −1 based on the US erosion rate (7.1 Mg ha −1 yr −1 ) for maize areas estimated from this study. The uncertainty associated with the assumption of the N content will scale linearly with the assigned N value and the C:N ratio of typical agricultural soils ranges from 10:1 to 15:1 be found in organic soils (Quinton et al., 2010). At a lower N content (15:1 C:N ratio), the estimated losses would range from 3.4 to 6.9 kg N ha −1 yr −1 . Among the four different management scenarios, "full replacement" scenario was the most representative of current US maize management. The average erosion N loss from "full replacement" scenario in this study was 7.1 kg N ha −1 yr −1 well within the range of estimated US erosion N fluxes.
We estimated erosion P loss using a SOM-P value of 0.1% based on the summaries in Quinton et al. (2010) resulting in a mean estimate of 3.6 kg P ha −1 yr −1 . The estimated erosion P loss in the literature (Quinton et al., 2010;Smil, 2000) derived from the US erosion rate (7.1 Mg ha −1 yr −1 ) for maize areas ranges between 3.4 and 7.3 kg P ha −1 yr −1 . As with N losses in erosion, the total estimated flux will vary with actual P content of simulated soils over time; this is an area for future model improvement.

Scenario-Based Analysis of Degradation Effects on Yields
The time series of national average yields for the four scenarios are shown in Figure 1. In the "full replacement" scenario average yields at the beginning of the scenario were 8.4 Mg ha −1 and decline slightly to about 8.2 Mg ha −1 over 100 years. The relative consistency in yields in the "full replacement" scenario is possible because the model continually adjusts the addition of nutrients and water to meet the difference between plant demand relative to the resources naturally available from soils. The contrasting "no inputs" scenario has no fertilization or irrigation and relies solely on naturally available nutrients and water. This scenario began with average US yields around 6.6 Mg ha −1 yr −1 but yields declined 79% to about 1.4 Mg ha −1 yr −1 after 100 years (Figure 1). There are few modern analogs to this simulation in developed economies. However, in a long-term study of maize production under rainfed and unfertilized conditions in Ontario, Canada average yields declined 81% from 4.3 Mg ha −1 in the early 1960s to 0.8 Mg ha −1 by the early 1990s (Drury & Tan, 1995 (Bakker et al., 2004). Such declines in yields are also broadly consistent with the type of yield decline observed in nonirrigated and low fertilization systems in SubSaharan Africa (Lal, 1988;Sileshi et al., 2010). In these cases, the cumulative loss of nutrients over the time since conversion is responsible for the resulting soil nutrient stocks and nutrient supply available from natural biological or physical nutrient mineralization processes in the agricultural field.
The relative effects of fertilization and irrigation vary over time as nutrient limitation intensifies with continued cropping (absent fertilizer addition) as shown by the two additional scenarios in Figure 1.
At the beginning of the 100-year period, there was a 1.8 Mg ha −1 yr −1 average difference between the "no inputs" versus "full replacement" scenarios. This difference represents inherent site yield differential when fields are at relatively high natural fertility (∆Yield Inherent ). At this early stage, water addition (shown in the "irrigation" scenario) has a larger effect than fertilization on mean national yields because crop nutrients are provided by the mineralization of preexisting SOM. The spatial differences in yields at the beginning of the simulation period are shown in Figure 2a. Spatial differences in yields by the end of the simulation period are shown in Figure 2b and discussed below.
As the simulation proceeds, there is N release from SOM and fertilizer input (including cumulative accumulation and release from fertilizer N and from SOM stocks). The sum of these inputs is shown as the N release rate in Figure 3. At the start of the simulation the initial N available for plant growth in all scenarios is well in excess of annual plant demand but declines very quickly as soil N stocks are depleted. At the beginning of the simulation, the N release rate provides 239%-279% of total plant N demand but was at 30%, 123%, 28%, 130% of plant demand for the no inputs, fertilization, irrigation, and full replacement scenarios at the end of the simulation. The N release rate at the end of simulation (average for the last 5 years) for the "no inputs" scenario is 20.5 kg ha −1 yr −1 (N mineralization of SOM in this scenario) while the amount of N release rate at the end of simulation (average for the last 5 years) for "full replacement" scenario is 90.5 kg ha −1 yr −1 (inclusive of both SOM mineralization and fertilizer release). The difference at the end of simulation between two scenarios is roughly equivalent the annual average fertilizer rate of the "full replacement" scenario.
We also investigated P release from SOM and fertilizer input (including cumulative accumulation and release from fertilizer P and from SOM stocks, Figure 4). The P release rate is the sum of the P derived JANG ET AL. The "no inputs" scenario shows the decline in yields over time without fertilization and irrigation (slope = −0.05 Mg ha −1 yr −1 , r 2 = 0.85, p < 0.001). "Fertilization" includes full nutrient replacement and is relatively stable (slope = 0.00 Mg ha −1 yr −1 , r 2 = 0.04, p < 0.05) whereas the "irrigation" scenario includes full water replacement and a large decline over time (slope = −0.05 Mg ha −1 yr −1 , r 2 = 0.80, p < 0.001).
from the mineralization of SOM plus the amount of P that goes from the active pool to labile pool. At the beginning of the simulation, the initial P available for plant growth in all scenarios is more than enough annual plant demand but declines very quickly as soil P stocks are reduced. Without fertilization, the P release rate tend to be zero after a few years of simulation. At the beginning of the simulation, the P release rate provides 153%-193% of total plant P demand but was at 16% (2.5 kg ha −1 yr −1 ), 77% (12 kg ha −1 yr −1 ), 15% (2.3 kg ha −1 yr −1 ), 96% (14.9 kg ha −1 yr −1 ) of plant demand for the no inputs, fertilization, irrigation, and full replacement scenarios at the end of the simulation.
Over time, the yield differences between the scenarios grow as soil nutrient status declines and there are large simulated shifts between N and P limitation through these simulations that influence the dynamics of the above graphs (Figures 1, 3, and 4). By the end of the scenario, the difference between scenarios was JANG ET AL.
10.1029/2020EF001641 9 of 19 Figure 2. Simulated yield difference between the "full replacement" and "no inputs" scenarios. Panel A shows ∆Yield Inherent (the difference between the "full replacement" and "no inputs" scenarios at the start of the 100-year simulation) and panel B shows ∆Yield Degradation (the difference between the two scenarios as shown in Figure 1). Note that the map coverage is based on locations where maize is grown in any abundance and does not represent spatially weighted yield production in the US.
6.8 Mg ha −1 yr −1 . Subtraction of the initial site yield difference from the total end of simulation yields a measure of the yield difference due to soil degradation (∆Yield Degradation ) equal to 5.4 Mg ha −1 yr −1 for the US. In this scenario, soil degradation is defined to include both erosion losses of nutrients from the system and the decline of SOM following cultivation (IPBES, 2018). Even with the potential uncertainties shown in Table S3 and S4, the broad patterns in nutrient loss and fertility replacement remain consistent across a range of input parameters. In the sections below, we separate these processes and discuss the implications of each process on overall site fertility. At the end of the simulation ∆Yield Degradation can be interpreted as the average amount of US maize yields that were once supported by natural processes that are now supported instead by exogenous nutrient inputs. The spatial distribution of this term is shown in Figure 2b.

Controls on Soil Degradation in US Maize Agriculture
There are several processes that influence soil degradation under long-term agricultural management. These processes include erosional losses of SOM N and P, declines in soil fertility through the leaching of N or P or gaseous loss of N, and lastly through the removal of crop biomass containing N and P. All of these processes can lead to rapid reduction in soil nutrient content and all can be mitigated to at least some extent by changes in management practices. Figures 5 and 6 shows that the ∆Yield Degradation is controlled primarily by the loss of soil fertility in the US maize agriculture. The figures illustrate both the wide range in estimates of ∆Yield Degradation and total NP stress across the US maize growing regions. The overall difference in yield associated with degradation processes is influenced strongly by yield potential. Maize grown in high yielding settings with ideal temperature and moisture conditions will also have the highest demand for nutrients and therefore will experience the greatest degree of nutrient stress in the absence of fertilization (Figure 2b). This result is not surprising but highlights the interaction between continual high yields and nutrient depletion in soils and the intersection between climate, irrigation, and potential nutrient depletion (Lobell, 2014).
To explore how nutrient and water limitations evolve over time in these scenarios, we show simulated nutrient and water stress factors in Figure 7. As expected, combined (NP) nutrient stress is highest without fertilization and at the end of the simulation period as soil nutrient stocks are depleted. N stress is highest at the beginning of the simulation period and declines over time as P stress increases. In EPIC, nutrient stress is simulated as serial stress where growth is most limited by either N or P. The pattern observed here suggests that P stress increases over the course of the simulation and eventually becomes the dominant JANG ET AL.
10.1029/2020EF001641 10 of 19 limitation after ∼50 years in the unfertilized scenarios. Water stress occurs in the two nonirrigated scenarios and is highest in the "fertilization" scenario because of higher simulated productivity and crop demand for We estimated declines in overall site fertility due to soil N losses by the calculation of erosion N loss, crop removal N, leaching, denitrification, and other export pathways (kg N ha −1 yr −1 ). Soil N loss is defined as the net change in soil N stocks after accounting for erosion N losses, crop removal N, and other N losses and ranged from 74 to 177 kg N ha −1 yr −1 ( Table 2). The largest source of N and P loss from the simulated maize fields at high yields is crop biomass. However, the additional losses for nitrogen are significant and dominated by the "other N losses" category which includes leaching and trace gas production. For phosphorus, the simulated erosion losses are similar to the total P removed in crops at low yields. These low yields are unlikely to occur in the US but are common in developing agricultural systems such as those in SubSaharan Africa suggesting that P loss in erosion (and mitigation strategies) could be a key issue in longterm site fertility. P losses are also of concern because of rising demand for P fertilization in agriculture and potential limitations to P supply in the future (Jarvie et al., 2015) as well as off-site impacts on water quality (Pimentel et al., 1995;Rickson et al., 2015). Because EPIC is a field rather that watershed scale model, we cannot simulate depositional processes on downstream fields so in such locations, we would expect to see some replacement of N and P on fields. The relative balance between erosion and depositional processes at a smaller scale is an important issue for future research in this area.

Estimation of the Economic Costs of Soil Degradation in Agriculture Operations
Over time, the amount of fertilizer required to support yields increases as natural nutrient supplies from natural lands decline through nutrient mining, erosion, and other export pathways. Here we refer to this shift from natural to industrial sources as compensatory fertilization and there are several different ways to estimate this value including the change in fertilizer additions, change in total nutrient stocks, and simulated declines in natural SOM mineralization. Average N fertilizer inputs at the start of the "full replacement" simulation are 126 kg N ha −1 yr −1 and increase to an average of 159 kg N ha −1 yr −1 after 100 years yielding an estimate of 33 kg N ha −1 yr −1 for compensatory N additions. For phosphorus, initial P fertilization is 19 kg P ha −1 yr −1 needed at the start of the "full replacement" simulation and increase to an average of 31 kg N ha −1 yr −1 after 100 years yielding an estimate of 12 kg N ha −1 yr −1 for compensatory P additions.  The EPIC model allows for an estimate of nitrogen required to offset the annual N loss in soils due to the decline in soil N content (soil N loss) plus the total export of N in erosion. This estimate varies considerably across the US but averages 28 kg N ha −1 yr −1 for soil N loss and 7.1 kg N ha −1 yr −1 for soil erosion leading to an estimate of about 35 kg N ha −1 yr −1 in compensatory N inputs. Phosphorus supply and compensatory fertilization is more difficult to estimate due to changes in the content of different P fractions in the soil (e.g., some added P fertilizer is stabilized into long residence time mineral fractions, particularly in many highly weathered [older] soils in tropics and subtropics, and becomes unavailable to plant). As an alternative, we calculate the decline in simulated P mineralization from the beginning to the end of the "full replacement" scenario to estimate compensatory P inputs and this totals 15 kg P ha −1 yr −1 . If we assume that similar losses occur in the "full replacement" scenario and add this to the P erosion rates of 3.6 kg P ha −1 yr −1 we estimate a total of 18.6 kg P ha −1 yr −1 in P decline that must be met with compensatory fertilizer inputs. Note this estimate does not include potential P sequestration onto soil surfaces or other forms of P loss.
JANG ET AL.

10.1029/2020EF001641
13 of 19 We estimated the erosion N and P loss offset cost ($ yr −1 ) and soil N loss offset cost ($ yr −1 ) in the US maize agriculture from the estimates of compensatory fertilization (Table S1). Fertilizer price data  are adapted from USDA Economic Research Service (USDA ERS, 2019) and we use the 2009-2014 average price for the cost estimation. N fertilizer costs are based on the average (per unit N) costs of the primary form of nitrogen fertilizer used in the US (i.e., anhydrous ammonia, nitrogen solutions, urea nitrogen, ammonium nitrate, and sulphate of ammonium) and estimated cumulative N fertilizer costs for the US maize agriculture per year which ranges from 531 to 3,351 million dollars. Erosion N loss (7.1 kg N ha −1 yr −1 for "full replacement" scenario) accounts for 4.5% of total N fertilizer applied (159 kg N ha −1 yr −1 ) and the erosion N loss offset costs to compensate soil fertility due to N losses by soil erosion range from 24 to 151 million dollars. Based on estimated soil N loss (28-33 kg N ha −1 yr −1 ) which is 18%-21% of total N fertilizer applied (159 kg N ha -1 yr -1 ), soil N loss offset costs are calculated and indicate a range of 144-1,154 million dollars (Table S1).
We estimated a fraction of phosphorus from the two primary sources of phosphate (i.e., superphosphate and diammonium phosphate) to calculate a range of P fertilizer costs for US maize agriculture (48-105 million dollars). The erosion P loss (3.6 kg P ha −1 yr −1 for "full replacement" scenario) offset cost accounts for 11.6 % of total P fertilizer applied (31 kg P ha −1 yr −1 ) and the erosion P loss offset cost ranges from 6 to 12 million dollars per year. We also calculated the decline in simulated P mineralization from the beginning to the end of "full replacement" scenario to estimate compensatory P inputs (15 kg P ha −1 yr −1 ). Considering P erosion rates for "full replacement" scenario (3.6 kg P ha −1 yr −1 ) a total of P decline was calculated as 18.6 kg P ha −1 yr −1 (Table S2). Although there will be large losses in organic P during agriculture (similar to soil N loss), there is insufficient mechanistic representation of the soil P pools and transformation in EPIC to accurately model or otherwise estimate these losses at the current time. Timeseries plots for serial (e.g., most limiting) nutrient and water stress (day yr −1 ) in the four simulation scenarios. Stress occurs when there is insufficient nutrient or water supply to meet daily growth potential. Panel A shows the number of days with combined NP stress, Panel B shows the number of days with N stress, Panel C shows the number of days with P stress, and Panel D shows the number of days with water stress.

(a) (b) (c) (d)
The patterns illustrated in this model simulation reflect broadly understood patterns in agriculture and reproduce the dynamics of SOM in intensely cultivated soils around the world (FAO, 2017). However, to date there has been limited evidence for the magnitude of nutrient loss (and associated compensatory fertilization) or estimation of the costs of the losses. For the US as a whole, approximately one-third of fertilizer application by the end of the 100-year simulation period is used to offset lost fertility due to degradation. At current fertilizer costs, these compensatory inputs amount to a total economic (for fertilizer purchases) ranging between 148 and 953 million US dollars depending on the form of fertilizer used (Table 3). These estimates are for fertilizer inputs only and are likely to be a small fraction of the larger environmental and economic impacts of increased fertilizer use which may reach well into the billions of dollars (Pimentel et al., 1995;Rickson et al., 2015).

Conclusions
The estimates of increased fertilizer input costs due to the loss of soil fertility over time presented here are conservative (low) because they are based on simulations that mirror contemporary conservation tillage instead of the conventional tillage systems more common in the US over most of the past 100 years. Conventional tillage typically leads to a 1.3-4.0 fold increase in erosion relative to conservation tillage (Seitz et al., 2019;Uri et al., 1999). Despite these conservative assumptions, the cost of compensatory fertilizer addition alone represents a significant operating expense in these systems. On the other hand, some aspects of the complex agricultural activity are not considered in this study. In fact, since all the models are an approximation of the reality (Box & Draper, 1987), what is simulated here with the EPIC model is a simplified version of the real agricultural systems and aspects related to the farmers' decision making process. In reality, farmers adapt to variations of factors that affect the crop production process (Gibbons & Ramsden, 2008;Keshavarz & Karami, 2014;Roesch-McNally et al., 2017).
As a response to soil quality decline, we can imagine different adaptation strategies such as the cultivation of different varieties, different crops, or the improvement of the crop residues management. This aspect can slow down the soil quality degradation over time reducing the cost of the compensatory fertilizer addition while introducing some uncertainty to the results of this study. Nonetheless, the decline in global soil quality and associated impacts on food security in developing economies has been well established as a critical issue within the global food system and one that includes a range of additional degradations mechanisms such as compaction or salinization not evaluated here (IPCC, 2019). Conversely, the restoration of soil quality through improved management in degraded systems can sequester soil carbon, improve water holding capacity, and decrease nutrient losses to the atmosphere and freshwater systems (Goyal et al., 1999;Lal, 2015). Supplementary organic matter inputs in organic farming systems have been shown to increase soil nutrient content over time (Clark at al., 1998) and the combination of inorganic and organic fertilization in conventional systems also can improve soil fertility across a diverse range of agricultural ecosystems (Guo et al., 2016;Mando et al., 2005;Šimon & Czakó, 2014

Table 3 Estimation of the Annual Economic and Environmental Costs of Soil Degradation in Maize Production in 2014 US Dollars
The work described here highlights the opportunity even in advanced high yield maize agriculture to reduce costs through practices that reverse land degradation and restore soil fertility (Crews & Rumsey, 2017).
In addition, there are other degradation processes including soil compaction that are currently difficult to simulate in models but which can result in increased erosion and decreased crop production (Colombi & Keller, 2019;Hargreaves et al., 2019;Langmaack et al., 2002). In general, these results highlight that practices that increase SOM and natural nutrient cycling over time could result in significant reductions in the application of fertilizers that have substantial environmental and economic impacts on freshwater quality and greenhouse gas production (Sobota et al., 2015;Von Blottnitz et al., 2006).

Data Availability Statement
This study uses a simulation model, developed with the open-source/process-based model EPIC available for download at https://epicapex.tamu.edu/model-executables/. The model input data presented in this study may be seen in Table 1 which is obtained from publicly available online sources. The climate data can be obtained from AgMIP Climate Forcing Data sets (https://data.giss.nasa.gov/impacts/agmipcf/agmerra/). The soil data is publicly available at http://www.fao.org/soils-portal/soil-survey/soil-maps-and-databases/ harmonized-world-soil-database-v12/en/. The land cover data can be obtained from CropScape -Cropland Data Layer (https://nassgeodata.gmu.edu/CropScape/) and the slope data is available publicly available at https://www.ngdc.noaa.gov/mgg/topo/globe.html. The crop calendar (management) data can be retrieved from Crop Calendar Data set (http://nelson.wisc.edu/sage/data-and-models/crop-calendar-dataset/index. php) and the soil carbon data can be obtained from https://daac.ornl.gov/cgi-bin/dsviewer.pl?ds_id=1238. The output data generated by the model is analyzed by modified Python scripts which are available upon request from the corresponding authors.