Estimation of the P Fertilizer Demand of China Using the LePA Model

Modern phosphate (P) fertilizers are sourced from P rock reserves, a finite and dwindling resource. Globally, China is the largest producer and consumer of P fertilizer and will deplete its domestic reserves within 80 years. It is necessary to avoid excess P input in agriculture through estimating P demand. We used the legacy P assessment model (LePA) to estimate P demand based on soil P management at the county, regional, and country scales according to six P application rate scenarios: (1) rate in 2012 maintained; (2) current rate maintained in low-P counties and P input stopped in high-P counties until critical Olsen-P level (CP) is reached, after which rate equals P-removal; (3) rate decreased to 1–1.5 kg ha⁻¹ year⁻¹ in low-P counties after CP is reached and in high-P counties; (4) rate in each county decreased to 1–8 kg ha⁻¹ year⁻¹ after soil Olsen-P reached CP in low P counties; (5) rate in each county was kept at P-removal rate after reduction; (6) P input was kept at the rate lower than P-offtake rate after reduction. The results showed that the total P fertilizer demand of China was 750 MT P₂O₅, 54% of P fertilizer can be saved from 2013 to 2080 in China, and soil Olsen-P of all counties can satisfy the demand for high crop yields. The greatest potential to decrease P input was in Yangtze Plain and South China, which reached 60%. Our results provide a firm basis to analyze the depletion of P reserves in other countries.

INTRODUCTION

Phosphorus is a major limiting nutrient for crop production in many regions (Koning et al., 2008; Sattari et al., 2014; Jarvie et al., 2015). The application of phosphate (P) fertilizers is a necessary practice for sustaining high crop yields in low P soils. However, P rock reserve, as a finite natural resource, is being rapidly consumed because of high P demand in agriculture and will be used up in the next century (Syers et al., 2008; Cordell et al., 2009; Gilbert, 2009). More recent estimates of global P rock stocks suggest that at current consumption rates for the production of P fertilizer and feed P, P rock reserves will be available in the next 300–350 years [calculated from USGS (2015; 2016)]. China plays a critical role in global P consumption and production (Sattari et al., 2014). Its consumption of P fertilizers increased from 1.0 Tg in 1978 to 5.0 Tg in 2017, and P fertilizer production increased from 1.0 to 7.2 Tg P during the same period (Zhang et al., 2019). This inevitably led to a sharp decline in P rock reserves. With the continuous high P input in recent years, P loss from croplands is responsible for 11.4% of the total P discharge, causing severe eutrophication of water bodies in China (Qiu, 2010; Liu et al., 2016). To reduce the depletion of P rock and alleviate soil P loss, the optimization of P input and improvement of crop P use efficiency are always the first measures of choice.

The building-up and maintenance approach, which is based on soil P fertility and P budget between input and removal, is useful to manage soil P and to avoid excessive P input (Li et al., 2011; 2015). In brief, depending on soil Olsen-P levels: high, low, or optimal, P input is less than, more than, or equal to crop P removal, respectively (Li et al., 2015). The critical level of soil Olsen-P for crop yield is defined as the level above which crop yield does not further increase with additional P fertilizer (Bai et al., 2013). In contrast, the critical level for P loss is the change-point in soil Olsen P and CaCl₂ P above which P loss shows a sharp increase (Heckrath et al., 1995). Maintaining soil Olsen-P between the two critical levels mentioned earlier does not only match the demand for high crops yield but also minimizes P loss to the environment (Johnston and Poulton, 2019).

Only 20% of P fertilizer is reported to be apparently used by crops in the first growing season in China (Zhang et al., 2008), and the remainder accumulates as legacy P in the soil. The amount of legacy P in arable land in China is more than 242 kg P ha⁻¹ (Li et al., 2011). Many researchers suggested that legacy P is not irreversibly fixed by soil minerals and that it could still be absorbed by crops in the following seasons (van Keulen et al., 2000; Mishima et al., 2010; Sattari et al., 2012). Considering legacy P, the inputs of inorganic P fertilizers can be saved, and the risk of P transfer to water bodies can be reduced (Johnston and Poulton, 2019; Zhang et al., 2019; Pavinato et al., 2020). Sattari et al. (2012) have estimated lower P fertilizer demand at continental and global scales using the dynamic phosphorus pool simulator (DPPS) model, in which legacy P is considered as a usable resource for crops (Janssen et al., 1987; Wolf et al., 1987). The 20% of P fertilizer could be saved until 2050 by accounting for legacy P in China (Sattari et al., 2014). The prediction of P demand by the DPPS model is based on accurately simulating long-term historical crop P uptake in a yearly temporal scale at provincial or even larger spatial scales. However, the DPPS model is not effective for soil P fertility simulation and ignores the imbalance between counties in China. The legacy phosphorus assessment (LePA) model, which includes three P pools, can be used to estimate P demand according to various soil P statuses (Yu et al., 2021). The LePA can simulate crop P uptake by season, as well as soil TP and Olsen-P. Moreover, the LePA model considers the influence of soil type on crop P uptake.

In this study, we used the LePA model to estimate the future P demand per county according to the building-up and maintenance approach in the crop production systems. The objectives of this study are as follows: 1) to estimate P fertilizer demand of China from 2012 to 2080 at the county, regional, and country scales under several scenarios, 2) to predict consequential changes of soil Olsen P and legacy P in China during the simulation period, and 3) to explore the optimal P reduction strategy at the county scale.

MATERIALS AND METHODS

Study Site and Experimental Design

In this study, the LePA model was developed according to six typical long-term fertilizer experiments, which were distributed in five P management regions (Li et al., 2015). The six long-term fertilizer experiments were located in: Yangling in Northwest China (34°17′N, 108°00′E, 525 m altitude); Zhengzhou in North China (113°41′N, 35°00′E, 59 m altitude); Gongzhuling in Northeast China (43°30′N, 124°48′E, 220 m altitude); Chongqing in Yangtze Plain (30°26′N, 106°26′E, 266 m altitude); Qiyang in South China (26°45′N, 111°52′E, 120 m altitude); and Jinxian in South China (28°35′N, 116°17′E, 41 m altitude). Experiments began in 1980 at Jinxian and in 1990 at Yangling, Chongqing, Zhengzhou, Gongzhuling, and Qiyang. Soil physical and chemical properties at the beginning of the six experiments are shown in Table 1. The cropping systems are as follows: wheat (Triticum aestivum L.)–maize (Zea mays L.) rotation at Yangling, Qiyang, and Zhengzhou; rice (Oryza sativa L.) –wheat rotation at Chongqing; rice–rice rotation at Jinxian; and continuous maize cultivation at Gongzhuling.

TABLE 1

TABLE 1. Initial soil properties in long-term experiment sites.

The treatments arranged in a randomized block design in the six experimental sites included: 1) N (only nitrogen fertilizer application), 2) NP (N and P fertilizer application), 3) NK (N and potassium fertilizer application), 4) NPK (N, P, and K fertilizer application), 5) NPKS (N, P, and K fertilizer application combined with straw return), 6) NPKM [NPK combined with manure (100%)], 7) 1.5NPKM (NPK combined 150% manure), and 8) NPKM2 (NPK combined with 200% manure). To avoid the effects of N deficiency on crop yield, only data of treatments with N input were investigated to determine the correlation between soil Olsen-P and crop yield. Annual applications of fertilizers, manures, and straw in each experimental site are listed in Table 2. The N input from both chemical fertilizers and manure was at a ratio of 7:3 in manure treatments. All P and K fertilizers were applied as basal dressings.

TABLE 2

TABLE 2. Rate of fertilizer application each year for all experimental sites.

Model Description and Data Used

The LePA—a simple three-pool P-model including available (P_A), unavailable (P_U), and legacy P pools (P_LE)—was used (Yu et al., 2021). The soil total P amount (TPA) consists of the three P pools. Phosphorus can be absorbed by crops from P_A in the growing season, whereas phosphorus in P_U cannot be absorbed. P_LE is composed of the remaining P in P_A after being absorbed by crops and P in P_U. Also, phosphorus in P_LE can be absorbed and utilized by the next crop. This model thus considers the contribution of legacy soil P to subsequent crop production, which cannot be absorbed by seasonal crops. The LePA requires some parameters to perform simulations, and these parameters are given by the empirical equations between soil total P and Olsen-P, soil Olsen-P and crop yield, and soil TPA (total P amount) and P_A (available pool) of the long-term fertilizer experiments. These equations are shown in Supplementary Table S1. The model can calculate the P uptake by crops, pool sizes, and soil Olsen-P with a seasonal time step. It also can be applied for calculating the future P fertilizer requirement for a target yield according to various soil P statuses. The model was used to predict future P fertilizer demand and the changes in soil P fertility in two field experiments of different soil types (Yu et al., 2021). The model description and its application to estimate future P requirements in detail are given in Yu et al. (2021).

Data on crop yield, soil TP, and Olsen-P were collected from the published literature shown in Supplementary Tables S2–S8. The data reported as graphs in the original publications were extracted using GetData Graph Digitizer (version 2.25, developer: S. Fedorov). Harvested P was calculated as yield × P concentration of the harvested product. The rate of P supply from weathering was 1.6 Tg year⁻¹ (Liu et al., 2008) and equal to 1 kg ha⁻¹ year⁻¹. Atmospheric P deposition was 0.25 kg ha⁻¹ year⁻¹ (Liu et al., 2008). The P loss in Northwest, North, Northeast, Yangtze Plain, and South China were 0.36, 1.30, 0.22, 1.05, and 1.34 kg P ha⁻¹ year⁻¹, respectively, from calculation data of the cultivated land area and P loss in each province (MNR, 2012; Zhang et al., 2019).

Evaluation of Legacy Phosphorus Assessment Performance

The performance of the LePA model was evaluated by comparing the simulated and observed crop P uptake, soil TP, and Olsen-P for the NK, NP, NPK, HNPK, and NPKS treatments. The LePA model is an empirical model based on the correlations between soil TP, Olsen-P, and crop yield of the long-term fertilizer experiments, and thus, no parameters need to be calibrated. Model performance in validation was evaluated using the coefficient of determination (R²) and root mean squared error (RMSE) between simulated and observed values. The higher the R², the lower the RMSE, and the better is the model performance. The slope and intercept are also used to confirm the performance of the modeling. Also, the closer the slope is to 1.0, the lower the intercept, and the better is the model performance. The comparison between model simulations and observations from pieces of literature was performed with Microsoft Excel (version 2016, Microsoft Corporation, Redmond, Washington, United States).

$\begin{array}{l}{R}^2=\frac{{\displaystyle \sum {\left(X_i-{\overline{Y}}_i\right)}^2}}{{\displaystyle \sum {\left(Y_i-{\overline{Y}}_i\right)}^2}}\\ RMSE=\sqrt{\frac{{\displaystyle \sum _{i=1}^n{\left(X_i-Y_i\right)}^2}}{n}}\end{array}$

where the sample size is denoted by n, Yi representing the ith observed value, and Xi corresponds to the ith model simulated value.

Application Methodology

The year 2012 was defined as a basal year. The soil Olsen-P per county in 2012 is shown in Figure 1A (https://kxsf.soilbd.com/Index.aspx). The critical levels of soil Olsen-P for crop production in the five P management regions were 20 mg kg⁻¹ in Northwest China, 15 mg kg⁻¹ in North, Northeast China, and Yangtze Plain, and 35 mg kg⁻¹ in South China (Supplementary Figure S2; Li et al., 2015). Maintaining soil Olsen-P near the critical value is considered to be in a good soil P status. The critical level of soil Olsen-P for P leaching potential and the highest annual crops P uptake assessed on average according to the six long-term fertilizer experiments in China and ecological zones are shown in Supplementary Table S10. We ran the LePA model according to the building-up and maintenance approach at the county scale to estimate P fertilizer demand and change of soil P fertility of each county during 2012–2080, assuming no changes in the cultivated area and main cropping systems during the simulation period. The distribution of the considered dominant cropping systems in China is shown in Figure 1B, including single-season maize (M), wheat-maize (W-M), wheat-rice (W-R), and rice-rice (R-R). This simulation covers 2,283 counties (including districts in the urban area) except those in Tibet and Qinghai provinces, which are not main food cropping areas. The amounts of P input in farmer’s practice were 53 kg P ha⁻¹ year⁻¹ for wheat, 41 kg P ha⁻¹ year⁻¹ for maize, and 30 kg P ha⁻¹ year⁻¹ for rice, which are averages of P inputs during 2008–2018 according to a farmer survey. The initial value of soil TPA (TPA₀) in each county in 2012 can be calculated by soil Olsen-P and several corresponding correlations between soil TP, Olsen-P, and crop yield.

FIGURE 1

FIGURE 1. Soil Olsen-P in each county in 2012 (A) and distribution of considered dominant cropping systems in China (B). M, W-M, W-R, and R-R denote single maize, wheat-maize, wheat-rice, and rice-rice cropping systems, respectively.

Scenarios for the Period 2012–2080

In scenario 1, the P application rate per county is the same as in 2012 from 2013 to 2080. Besides, we set up the other five scenarios aimed to reduce P input without any food security risk. In scenario 2, the current rate of P application according to farmer’s practice is kept in low P counties until the soil Olsen-P rises to the critical level; in high P, P fertilizer application stops until soil Olsen-P drops to the critical value; and then, P input is equal to crop P removal to maintain soil Olsen-P always at the critical level.

We assumed that the current risk tolerance of the P fertilizer industry is a 5–10% decrease of P fertilizer consumption each year in scenarios 3–6, above which P fertilizer producers cannot afford the dramatic shrinking market. In scenario 3, P application decreases at a rate of 1–1.5 kg ha⁻¹ year⁻¹ from 2012 to 2080 in each high P county, whereas soil Olsen P in low P counties is first built up to the critical levels at the current rate of P fertilizer application, and then, it is treated as a high P county. In scenario 4, P input in each high P county is reduced to a rate that is affordable by the P-fertilizer industry (1–8 kg ha⁻¹ year⁻¹), and P input is maintained at the rate in 2012 in low P counties until soil Olsen P reaches to a critical level and then reduced as the high P county. In this scenario, an increase of P application occurs at a rate of 0.1–20 kg ha⁻¹ year⁻¹ per county. In scenario 5, the only difference with scenario 4 is that the rate of P application in each county is kept at the P-offtake or P-removal rate after reduction. In scenario 6, P input was kept at a rate lower than the P-offtake rate after a reduction so that soil Olsen-P in each county reaches the critical level by 2080.

RESULTS

Model Validation

Crop P uptakes simulated by LePA in the NK, NP, NPK, HNPK, and NPKS treatments at Jinxian, Chongqing, Zhengzhou, and Gongzhuling sites are shown in Supplementary Figures S4–S7. The R² of the correlations between simulated and observed P uptake of early rice in Jinxian (n = 84, p < 0.05) was 0.03 with an RMSE of 2.01, and the R² and RMSE for late rice (n = 84, p < 0.01) were 0.09 and 2.19, respectively (Supplementary Table S9). In Chongqing, the R² for wheat P uptake (n = 76, p < 0.01) was 0.55 with an RMSE of 2.43 and was much higher than that of rice (R² = 0.06; RMSE = 2.74). A large variation in crop yield caused the poor fitting of the LePA model. The R² and RMSE for wheat P uptake in Zhengzhou (n = 56, p < 0.01) were 0.70 and 4.93, respectively, and they were 0.25 and 5.34, respectively, for maize P uptake (n = 56, p < 0.01). In Gongzhuling, the R² for maize P uptake was 0.08 with an RMSE of 4.20 (n = 100, p < 0.01). The simulation results of crop P uptake, soil TP, and Olsen-P by LePA at Yangling and Qiyang sites are shown in Yu et al. (2021).

The comparison of simulated and observed soil TP and Olsen-P in Jianxian, Chongqing, Zhengzhou, and Gongzhuling are shown in Supplementary Figures S8, S9. The R² values ranged from 0.78 to 0.99 for total P at the above four sites, and the regression lines are close to the 1:1 lines (Supplementary Figure S8). The R² was 0.78, and the RMSE was 8.58 for soil Olsen-P in Jinxian, and they were 0.69 and 4.87 and 0.58 and 12.68 in Chongqing and Zhengzhou, respectively (Supplementary Figure S9). For Gongzhuling, the R² for soil Olsen-P was only 0.08 with an RMSE of 11.26 (n = 57, p < 0.05).

P Demand for Crop Production in China During 2012–2080

In scenario 1, the annual P input of each county is the same as in 2012 from 2013 to 2080 (Supplementary Excel S1; Supplementary Movie S1A). The annual P application rates were thus maintained at 84 kg ha⁻¹ in Northwest China, 94 kg ha⁻¹ in North China, 41 kg ha⁻¹ in Northeast China, 83 kg ha⁻¹ in Yangtze Plain, 74 kg ha⁻¹ in South China, and 78 kg ha⁻¹ in whole China during this period (Figures 2A–F). In scenario 2, the rate of P application is kept or halted in low or high counties in the first few years, and then, P input is equal to crop P removal (Supplementary Excel S1; Supplementary Movie S1B). The change in P application rates according to scenario 2 was thus different for each region over time. The average application rate was reduced to 54 kg ha⁻¹ in China. The lowest P application rate was 4 kg ha⁻¹ in Northeast China, and the highest was 69 kg ha⁻¹ in South China. The annual P application rate gradually decreased to crop P-offtake from 2013 to 2080 in Northwest China, Yangtze Plain, and South China. In contrast, annual P application rates increased each year to meet crop P-offtake in North China and Northeast China. The annual P input per county gradually decreased from 2012 to 2080 in scenario 3 (Supplementary Excel S1; Supplementary Movie S1C). Annual P demand decreased from 84 to 1 kg ha⁻¹ year⁻¹ in Northwest China, from 94 to 6 kg ha⁻¹ year⁻¹ in North China, from 41 to 13 kg ha⁻¹ year⁻¹ in Northeast China, from 83 to 1 kg ha⁻¹ year⁻¹ in Yangtze Plain, and from 74 to 1 kg ha⁻¹ year⁻¹ in South China. On average, P application rates decreased from 78 to 3 kg ha⁻¹ year⁻¹at a rate of 1 kg ha⁻¹ year⁻¹ in China from 2012 to 2080 in scenario 3. The P demand at the county scale initially decreased and then increased to meet crop P-offtake and maintained this level in scenario 4 (Supplementary Excel S1; Supplementary Movie S1D). In this scenario, the rate of the decrease in annual P demand was 2 kg ha⁻¹ year⁻¹ in China, which is higher than that in scenario 3. In the whole of China, the simulated annual P demand increased from 17 kg ha⁻¹ year⁻¹ in 2037 to 35 kg ha⁻¹ year⁻¹ by 2063 (Figure 2F). Similarly, the modeled annual P fertilizer demand increased to equal the rate of uptake-offtake during a period of approximately 15 years in the five P management regions. The rate of decrease in the P application per county simulated in scenarios 5 and 6 was the same as in scenario 4 (Supplementary Excel S1; Supplementary Movies S1E,F). Scenarios 5 and 6 project that the annual P demand of the five regions and China decreased to crop P-offtake or lower after the first 20 years and then maintained until 2080.

FIGURE 2

FIGURE 2. Trends of annual P application and P uptake in cropland for period 2012–2080 according to six scenarios in (A) Northwest China, (B) North China, (C) Northeast China, (D) Yangtze Plain, (E) South China, and China (F). Rate in 2012 maintained (PF1); current rate maintained in low-P counties and P input stopped in high-P counties until critical Olsen-P level (CP) is reached, after which rate equals P-removal (PF2); rate decreased to 1–1.5 kg ha⁻¹ year⁻¹ in low-P counties after CP is reached and in high-P counties (PF3); rate in each county reduced to one affordable by P fertilizer industry (1–8 kg ha⁻¹ year⁻¹) after soil Olsen-P reached CP in low P counties (PF4); rate in each county was kept at P-removal rate after reduction (PF5); P input was kept at rate lower than P-offtake rate after reduction so that soil Olsen-P in each county reaches critical level by 2080 (PF6).

The P fertilizer demand of each county under the analysis of six scenarios between 2013 and 2080 is shown in Supplementary Excel S1. Total P demand per county was 1,000–3,900 kg P ha⁻¹ in 68 years in scenarios 2–6, whereas it was 2,788–6,392 kg P ha⁻¹ in scenario 1. On a regional scale, up to 51% of P demand decreased in Northwest China, 47% in North China, 33% in Northeast China, 63% in Yangtze Plain, and 60% in South China in the scenarios 2–6 compared with scenario 1 (Figure 3). By 2080, the P fertilizer demand of China is 1,640 MT P₂O₅ in scenario 1. After adopting the building-up and maintenance approach, the P fertilizer demand of China decreased to 750 MT P₂O₅ in scenario 2, 880 MT P₂O₅ in scenario 3, 750 MT P₂O₅ in scenario 4, 840 MT P₂O₅ in scenario 5, and 760 MT P₂O₅ in scenario 6 without yield loss. Also, it led to a 46–54% decrease in the demand for P fertilizer on a country scale.

FIGURE 3

FIGURE 3. Total amount of P demand for five regions and China in 2013–2080 according to six scenarios [Total P demand (P₂O₅, 100 MT) = total P demand for 69 years (P₂O₅, kg ha⁻¹) × cultivated land area in each region (ha) MNR, 2012]. Rate in 2012 maintained (scenario 1); current rate maintained in low-P counties and P input stopped in high-P counties until critical Olsen-P level (CP) is reached, after which rate equals P-removal (scenario 2); rate decreased to 1–1.5 kg ha⁻¹ year⁻¹ in low-P counties after CP is reached and in high-P counties (scenario 3); rate in each county reduced to one affordable by P fertilizer industry (1–8 kg ha⁻¹ year⁻¹) after soil Olsen-P reached CP in low P counties (scenario 4); rate in each county was kept at P-removal rate after reduction (scenario 5); P input was kept at rate lower than P-offtake rate after reduction so that soil Olsen-P in each county reaches critical level by 2080 (scenario 6).

Changes of Soil Olsen-P and Legacy P Between 2012 and 2080

The soil Olsen-P of each county was gradually increased with the traditional P application rates from 2012 to 2080 in scenario 1, and it reached 18–400 mg kg⁻¹ by 2080 (Figure 4A; Supplementary Movie S2A). In scenario 2, soil Olsen-P per county reached the critical level for crop yield by 2036 and was then maintained at this level until 2080 except the soil Olsen-P level in eight counties in South China, which were less than 10 mg kg⁻¹ by 2080 (Figure 4B; Supplementary Movie S2B). That was because the soil pH was below 5.0 in these counties and no more P fertilizer was applied during this simulation period. There was a continuous increase in soil Olsen-P in all counties with reducing P input at a low rate in the first 40 years according to scenario 3 (Figure 4C; Supplementary Movie S2C). Although the soil Olsen-P level of each county began to decline after 2052, the soil Olsen-P of counties in Yunnan, Sichuan, and Guizhou provinces were over 100 mg kg⁻¹ by 2080. Additionally, the soil Olsen-P level exceeded the critical level for P leaching potential (76.1 mg kg⁻¹ for Yangtze Plain and 78.2 mg kg⁻¹ for South China) in many counties in Yangtze Plain and South China. In scenario 4, soil Olsen-P per county increased firstly and then decreased with P application and finally remained near the critical level from 2060 to 2080 (Figure 4D; Supplementary Movie S2D). Furthermore, soil Olsen-P did not exceed 40 mg kg⁻¹ in 2080 in all counties. In scenario 5, soil Olsen-P in each county increased in the first 30 years, and then, it was mostly maintained at a 30–78 mg kg⁻¹, and it was much higher than the critical level for crop yield until 2080 (Figure 4E; Supplementary Movie S2E). Moreover, the soil Olsen-P level in some counties in Northwest China was above the P leaching level (36.9 mg kg⁻¹) from 2012 to 2080. In scenario 6, the soil Olsen-P of each county increased in the first 20 years and then decreased to the critical levels in 2080 (Figure 4F; Supplementary Movie S2F).

FIGURE 4

FIGURE 4. Soil Olsen-P at county level in 2080 under scenarios 1–6 (A–F). Rate in 2012 maintained (scenario 1); current rate maintained in low-P counties and P input stopped in high-P counties until critical Olsen-P level (CP) is reached, after which rate equals P-removal (scenario 2); rate decreased to 1–1.5 kg ha⁻¹ year⁻¹ in low-P counties after CP is reached and in high-P counties (scenario 3); rate in each county reduced to one affordable by P fertilizer industry (1–8 kg ha⁻¹ year⁻¹) after soil Olsen-P reached CP in low P counties (scenario 4); rate in each county was kept at P-removal rate after reduction (scenario 5); P input was kept at rate lower than P-offtake rate after reduction so that soil Olsen-P in each county reaches critical level by 2080 (scenario 6). Supplementary Movie S2 showing changes in soil Olsen-P in 2012–2080 is available in Supplementary Information. Panel 4 was created with ArcMap 10.2.

The changes in soil legacy P in each county from 2012 to 2080 according to the scenarios are shown in Figure 5 and Supplementary Movie S3. In 2012, the average amounts of soil legacy P calculated by LePA were 1,825 kg P ha⁻¹ in Northwest China, 1,785 kg P ha⁻¹ in North China, 1,829 kg P ha⁻¹ in Northeast China, 1,610 kg P ha⁻¹ in Yangtze Plain, 1,590 kg P ha⁻¹ in South China, and 1,699 kg P ha⁻¹ in China. Similar to the changes in soil Olsen-P, the soil legacy P in all counties was gradually increased during the simulation period and reached 1,532–8,054 kg P ha⁻¹ in 2080 in scenario 1 (Figure 5A; Supplementary Movie S3A). The soil legacy P of each county reached an appropriate level (corresponding to the critical level of soil Olsen-P) in 2036 in scenario 2 and in 2060 in scenario 4 and then was maintained at this level for several years until 2080 (1,970 kg P ha⁻¹ in Northwest China, 1,710 kg P ha⁻¹ in North China, 1,390 kg P ha⁻¹ in Northeast China, 1,620 kg P ha⁻¹ in Yangtze Plain, and 1,900 kg P ha⁻¹ in South China) (Figures 5B,D; Supplementary Movies S3B,D). That indicated the soil legacy P that can be used for crop production was 2,860 kg P ha⁻¹ in Northwest, 2,975 kg P ha⁻¹ in North, 965 kg P ha⁻¹ in Northeast, 3,535 kg P ha⁻¹ in Yangtze Plain, 3,030 kg P ha⁻¹ in South, and 2,865 kg P ha⁻¹ in the whole China. The soil legacy P increased first and then decreased to the appropriate level in 2080 in scenario 6 (Figure 5F; Supplementary Movie S3F). In scenario 3, the legacy soil P per county was gradually increasing from 2012 to 2052 and was more than 3,000 kg P ha⁻¹ in counties in South China (Figure 5C; Supplementary Movie S3C). Additionally, the legacy P at the county level reached a maximum of 4,227 kg P ha⁻¹ in 2080. The legacy P of each county increased to 1,240–3,860 kg P ha⁻¹ after 30 years, and then, it was maintained at this level until 2080 in scenario 5 (Figure 5E; Supplementary Movie S3E). In 2080, the legacy P at the county level in this scenario was the highest in scenarios 2–6 except in counties in Yangtze Plain and South China.

FIGURE 5

FIGURE 5. Soil legacy P at county level in 2080 under scenarios 1–6 (A–F). Rate in 2012 maintained (scenario 1); current rate maintained in low-P counties and P input stopped in high-P counties until critical Olsen-P level (CP) is reached, after which rate equals P-removal (scenario 2); rate decreased to 1–1.5 kg ha⁻¹ year⁻¹ in low-P counties after CP is reached and in high-P counties (scenario 3); rate in each county reduced to one affordable by P fertilizer industry (1–8 kg ha⁻¹ year⁻¹) after soil Olsen-P reached CP in low P counties (scenario 4); rate in each county was kept at P-removal rate after reduction (scenario 5); P input was kept at rate lower than P-offtake rate after reduction so that soil Olsen-P in each county reaches critical level by 2080 (Scenario 6). Supplementary Movie S3 showing changes in soil legacy P in 2012–2080 is available in Supplementary Information. Panel 5 was created with ArcMap 10.2.

Analysis of Optimal P Reduction Strategy at County Scale

Annual P demand per county according to the optimal P reduction scenario is shown in Supplementary Excel S1 and Supplementary Movie S4A. This scenario is a combination of scenarios 2 and 4 according to different soil P statuses in counties in 2012. For P-deficient counties, P application rates were maintained in the first few years, and then, P input was equaled to crop P removal; for P-optimal counties, P input was equaled to P uptake from 2013 to 2080, and for P-rich counties, the P input was initially reduced and then increased to meet crop P removal and maintained at this level. In this scenario, around 65% of the counties seem to keep the P application rates in 2012 for the first few years until the soil Olsen-P reached the critical level for crop production, and 4% of the counties seem to keep the P input equal to crop P removal from 2013 to 2080. There are also more than 700 counties with high P status that need to slowly decrease P input at a rate of 2 kg P ha⁻¹ year⁻¹ so that the changes in annual P demand in China have less impact on the P fertilizer industry. On average, P application rates decreased from 78 to 26 kg ha⁻¹ year⁻¹ at a rate of 2 kg ha⁻¹ year⁻¹ in China from 2012 to 2036 (Figure 6A) and increased to 34 kg ha⁻¹ year⁻¹ in 2063. Then, the annual P input of China remained at 34 kg ha⁻¹ year⁻¹ until 2080. The total P input of each county ranged from 1,000 to 3,570 kg P ha⁻¹ from 2013 to 2080. The total P fertilizer input of China was 750 MT P₂O₅, which was the same as in scenarios 2 and 4. In the optimal scenario, soil Olsen-P in P-optimal counties was always maintained at the critical level from 2012 to 2080 (Figure 6B; Supplementary Movie S4B). Soil Olsen-P in P-deficient counties was first gradually increased to the critical level and then kept at this level until 2080, whereas soil Olsen-P in P-rich counties increased firstly and then decreased to the critical level with P input and finally remained near the critical level. The changes of soil legacy P per county were similar to that of soil Olsen-P during this simulation period (Figure 6C; Supplementary Movie S4C).

FIGURE 6

FIGURE 6. Trends of annual P application and P uptake in cropland in China for period 2012–2080 (A), soil Olsen-P (B), and legacy P (C) per county in 2080 according to optimal scenario. In this scenario, rate maintained in first few years, and then P input was equaled to crop P removal in low P counties; P input was equaled to P uptake from 2013 to 2080 in P-optimal counties; rate was initially reduced and then increased to meet crop P removal and maintained at this level in high P counties. Supplementary Movie S4 showing P fertilizer demand per county, changes in soil Olsen-P, and legacy P in 2012–2080 is available in Supplementary Information.

DISCUSSION

Contribution of Legacy P to Future P Demand

Assuming the empirical equations and cropping systems used in the analysis in each county are the same as the near long-term fertilizer experiments, this study ignores the specificity of regions and cropping systems to estimate P demand from county to regional and national scales. Our results show that accounting for the critical role of legacy soil P leads to lower projections of P demand in scenarios 2–6 without yield loss for the period 2012–2080. Many previous studies have shown that crop yields in soils with a continuous P surplus for many years can be maintained or even increased a bit when P fertilizer application is reduced (Kamprath, 1967; Mishima et al., 2010; Verloop et al., 2010). As described in the conceptual model, only low rates of P fertilizer need to be applied to obtain good crop yields when soil available P (Olsen-P) is sufficient (Li et al., 2011). In this study, in counties where soil Olsen-P is lower than the critical level in 2012, the P input should be maintained the same as in 2012 for several years until the soil Olsen-P reaches the critical level to build up a good soil P status. This result is consistent with conclusions of Steen (1998), who suggested that 30–50% more P fertilizer than crop requirements must be applied for several decades to restore soil P in depleted soils. That indicates that a critical concentration of readily available P must be maintained to obtain a good crop yield. In addition to scenario 2, the demand for phosphate fertilizers estimated in scenarios 3–6 will increase or decrease at a rate of 1–2 kg P ha⁻¹ each year. This has a limited impact on the development of the P fertilizer industry in China. Moreover, the cumulative P fertilizer consumption of China in scenario 1 between 2013 and 2080 is 1,640 MT P₂O₅, and P rock demand is therefore 5.5 BT (P₂O₅: 30%). It exceeds the P rock reserves of China, which are only 3.3 BT (NBS, 2016). Actually, the P fertilizer demand can be reduced to 750–880 MT P₂O₅ through P management in scenarios 2–6, which enables saving up to 54% of P reserves. Sattari et al. (2014) also estimated future P demand in China’s arable land by DPPS model and showed that only approximately 20% of the P fertilizer input can be saved if the residual P is accounted for from 2010 to 2050 in China relative to the Rio+20 Trend scenario. That may be because the DPPS model ignores the difference in soil P status between the county levels when predicting the future demand for P fertilizer and the predicted years of the two studies are also different. For the five P management regions, P surpluses increased dramatically in North, South, and Yangtze plain between 1970 and 2010, causing a huge P reduction potential in these regions (Sattari et al., 2014). This result was consistent with our findings that P demand can be saved by nearly 60% in the Yangtze Plain and South China. The highest crop P uptake in Yangtze Plain and South China was lower than that of the other three regions, which were only 33 and 28 kg P ha⁻¹ year⁻¹, respectively (Supplementary Table S10). That caused the soil legacy P to increase dramatically with P application and lower demand for phosphate fertilizer to maintain the final steady-state (P inputs = P outputs) based on the building-up and maintenance approach.

Changes in Soil P Fertility Between 2012 and 2080

The changes in soil Olsen-P and legacy P simulated by LePA are different in the five P reduction scenarios 2–6 from 2012 to 2080. The decrease in soil Olsen-P and legacy P in the last 30 years under scenarios 3 and 6 indicated that soil P stocks would be depleted with continuously less P input than removal. This was consistent with the results of Sattari et al. (2014), who found the simulated trend of increasing P removal and declining P application rates could not continue over prolonged periods, as soil P stocks would be depleted, leading to soil degradation. Ideally, P removal and run-off losses need to be replenished by P inputs (Sattari et al., 2012). In this study, scenarios 2, 4, and 5 finally led to this ideal state where P input matched P removal and run-off losses. Many studies have shown that P-use efficiency can be improved without a yield loss when soil Olsen P is near the critical level for crop yield (Syers et al., 2008; Johnston and Poulton, 2019). However, the soil Olsen-P of some counties in Yangtze Plain, Northwest, and South China in scenarios 3 and 5 is much higher than the critical level from 2012 to 2080. As confirmed by some previous studies, excessive Olsen-P values above the critical P leaching level greatly increase the risk of P losses, leading to serious-environment problems (Sharpley et al., 2000; McDowell and Sharpley 2004; Aulakh et al., 2007). By using the globally differentiated eutrophication potential (eutrophication potential = P discharge × eutrophication potential factor) characterization method, previous research has also found out that some areas in China are facing huge eutrophication risks, especially in Yangtze Plain (Azevedo et al., 2013; Hellweg and Milà i Canals, 2014; Liu et al., 2016). Through the analysis of the P demand and the changes in soil P fertility according to scenarios 2–6, we can find that scenario 4 is the optimal P reduction strategy on a national scale.

Suggestions for Soil P Management at County Level

The model results show P input can be greatly saved by soil P management for optimizing P supply to meet the crop demand of each county. In this study, the building-up and maintenance approach based on soil P fertility and P budget between input and removal is used to manage soil P (Li et al., 2011; 2015). The P input is, respectively, less than, more than, or equal to crop P removal when the soil Olsen-P level is high, low, or optimal. Among scenarios 2–6 of different P reduction methods, considering the negative effect of scenario 2 on the phosphate fertilizer industry in China, scenario 4 is considered as the optimal P reduction strategy on a national scale. In China, the soil Olsen-P of cropland is varied among counties. To better manage soil P at the county level, we explored the optimal scenario for P reduction strategy at the county level in China according to different soil P statuses in counties in 2012. It is a combination of scenarios 2 and 4. When the soil Olsen-P in counties is low, P application rates were maintained in the first few years until soil Olsen-P reached the critical level, and then, P input was equaled to crop P removal, whereas P input was equaled to P uptake during the simulation period in optimal-P counties, and the P input was initially reduced and then increased to meet crop P removal and finally maintained at this level in high-P counties. In this optimal scenario, the annual reduction of phosphate fertilizer in China will not affect the phosphate fertilizer industry. The model results also provide how much P is needed in each county to achieve sufficient soil P fertility for crop production until 2080. Local governments could easily use it as targets at the county level. The LePA model can predict future P demand based on the soil P status of each county, which makes up for the deficiencies of previous related models such as the DPPS model (Wolf et al., 1987; Sattari et al., 2012; Sattari et al., 2014). This can effectively aid to prevent certain areas from facing environmental problems, whereas soil P fertility still needs to be improved in some areas when the total demand for P fertilizer is certain (Cao, 2006).

Although we explored the optimal P reduction strategy at the county level based on soil P management, the approach to extension is still a challenge for soil P management because smallholder farming dominates the agricultural landscape, and more than 90% of farmers hold only 0.1 ha of land in China (Zhang et al., 2013) The transfer rate of agricultural, scientific, and technological achievements is thus only 30–40% in China, whereas it is 70–90% in developed countries (Zhang et al., 2019). Some strategies need to be provided to encourage farmers to adopt soil P management approaches. Zhang et al. (2016) found that the Science and Technology Backyard platform was an efficient method for enabling smallholders to achieve yield and economic gains sustainably. Also, the conversion rate of technology increased from 31 to 53% via the Science and Technology Backyard platform. Moreover, the Chinese government has recently made some policies to develop large-scale management schemes (Long and Qu, 2018). A relatively large-sized farm is more conducive to the extension and implementation of soil P management approaches for reducing production costs. In addition, China lacks an environmental law focused on nutrient management in agricultural production, and profit is thus the major driver for farmers to adopt soil P management approaches (Zhang et al., 2019). Soil P management for agriculture in China requires both the scientists and the government to formulate relevant strategies jointly.

CONCLUSION

Our results suggest that legacy soil P can contribute to crop production. The LePA model shows that the total P demand of each county ranged from 1,000 to 3,570 kg P ha⁻¹ from 2013 to 2080. Also, the total P fertilizer demand of China was 750 MT P₂O₅, and 54% of phosphate fertilizer could be saved from 2013 to 2080. This greatly reduced the consumption of P rock reserves. Through the building-up and maintenance approach, the imbalance of soil P fertility between counties in each region is gradually eliminated. The demand for P fertilizer in China reduced from 78 to 35 kg P ha⁻¹ year⁻¹ based on soil P management, which indicated 55% of P fertilizer could be replaced by soil legacy P. The optimal P reduction strategy at the county level in China is a combination of scenarios 2 and 4, which achieves the greatest reduction potential of P fertilizers application with little effect on the phosphate fertilizer industry. For improved soil P management for agriculture in China, the effective extension requires both the scientists and the government to formulate relevant strategies jointly. Further work is needed to consider the contribution of root management to crop production, such as manipulating cropping systems, rhizosphere regulation, and microbial engineering. Our results provide a firm basis to analyze the depletion of P reserves in other countries.

Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material; further inquiries can be directed to the corresponding author.

Author Contributions

WY used the LePA model to analyze the data, estimate P demand, and write the manuscript. HL designed the study and revised the manuscript. PN revised the manuscript. GL, TM, JZ, and JS offered the suggestion.

Funding

This study is financially supported by the National Natural Science Foundation of China (31960627), the National Key Research and Development Program of China (2017 YFD0200200/2017YFD0200202), and the Deutsche Forschungsgemeinschaft (German Research Foundation)-328017493/GRK 2366.

Conflict of Interest

Author GL was employed by the company Israel Chemicals Ltd (ICL).

The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

The authors would like to thank Hanbing Cao for providing averages of P inputs of wheat, maize, and rice during 2008–2018 according to a farmer survey and Yuntianhua’s Yong Luo for advising on an appropriate annual reduction rate of P fertilizer consumption for the current risk tolerance of the P fertilizer industry.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fenvs.2021.759984/full#supplementary-material

References

Aulakh, M. S., Garg, A. K., and Kabba, B. S. (2007). Phosphorus Accumulation, Leaching and Residual Effects on Crop Yields from Long-Term Applications in the Subtropics. Soil Use Manage. 23, 417–427. doi:10.1111/j.1475-2743.2007.00124.x

CrossRef Full Text | Google Scholar

Azevedo, L. B., Henderson, A. D., van Zelm, R., Jolliet, O., and Huijbregts, M. A. J. (2013). Assessing the Importance of Spatial Variability versus Model Choices in Life Cycle Impact Assessment: the Case of Freshwater Eutrophication in Europe. Environ. Sci. Technol. 47 (23), 13565–13570. doi:10.1021/es403422a

PubMed Abstract | CrossRef Full Text | Google Scholar

Bai, Z., Li, H., Yang, X., Zhou, B., Shi, X., Wang, B., et al. (2013). The Critical Soil P Levels for Crop Yield, Soil Fertility and Environmental Safety in Different Soil Types. Plant Soil 372 (1-2), 27–37. doi:10.1007/s11104-013-1696-y

CrossRef Full Text | Google Scholar

Cao, N. (2006). Phosphorus Management in Chinese Agroecosystem Based on Forecasting the Change of Soil Phosphorus Fertility. Shanxi: Northwest A & F University in Chinese.

Google Scholar

Cordell, D., Drangert, J.-O., and White, S. (2009). The story of Phosphorus: Global Food Security and Food for Thought. Glob. Environ. Change 19 (2), 292–305. doi:10.1016/j.gloenvcha.2008.10.009

CrossRef Full Text | Google Scholar

Gilbert, N. (2009). Environment: The Disappearing Nutrient. Nature 461, 716–718. doi:10.1038/461716a

PubMed Abstract | CrossRef Full Text | Google Scholar

Heckrath, G., Brookes, P. C., Poulton, P. R., and Goulding, K. W. T. (1995). Phosphorus Leaching from Soils Containing Different Phosphorus Concentrations in the Broadbalk experiment. J. Environ. Qual. 24, 904–910. doi:10.2134/jeq1995.00472425002400050018x

CrossRef Full Text | Google Scholar

Hellweg, S., and Milà i Canals, L. (2014). Emerging Approaches, Challenges and Opportunities in Life Cycle Assessment. Science 344 (6188), 1109–1113. doi:10.1126/science.1248361

PubMed Abstract | CrossRef Full Text | Google Scholar

Janssen, B. H., Lathwell, D. J., and Wolf, J. (1987). Modeling Long‐Term Crop Response to Fertilizer Phosphorus. II. Comparison with Field Results 1. Agron.j. 79, 452–458. doi:10.2134/agronj1987.00021962007900030009x

CrossRef Full Text | Google Scholar

Jarvie, H. P., Sharpley, A. N., Flaten, D., Kleinman, P. J. A., Jenkins, A., and Simmons, T. (2015). The Pivotal Role of Phosphorus in a Resilient Water-Energy-Food Security Nexus. J. Environ. Qual. 44, 1049–1062. doi:10.2134/jeq2015.01.0030

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnston, A. E., and Poulton, P. R. (2019). Phosphorus in Agriculture: A Review of Results from 175 Years of Research at Rothamsted, UK. J. Environ. Qual. 48 (5), 1133–1144. doi:10.2134/jeq2019.02.0078

PubMed Abstract | CrossRef Full Text | Google Scholar

Kamprath, E. J. (1967). Residual Effect of Large Applications of Phosphorus on High Phosphorus Fixing Soils 1. Agron.j. 59, 25–27. doi:10.2134/agronj1967.00021962005900010007x

CrossRef Full Text | Google Scholar

Koning, N. B. J., Van Ittersum, M. K., Becx, G. A., Van Boekel, M. A. J. S., Brandenburg, W. A., Van Den Broek, J. A., et al. (2008). Long-term Global Availability of Food: Continued Abundance or New Scarcity? NJAS: Wageningen J. Life Sci. 55, 229–292. doi:10.1016/s1573-5214(08)80001-2

CrossRef Full Text | Google Scholar

Li, H., Huang, G., Meng, Q., Ma, L., Yuan, L., Wang, F., et al. (2011). Integrated Soil and Plant Phosphorus Management for Crop and Environment in China. A Review. Plant and Soil 349 (1-2), 157–167. doi:10.1007/s11104-011-0909-5

CrossRef Full Text | Google Scholar

Li, H., Liu, J., Li, G., Shen, J., Bergström, L., and Zhang, F. (2015). Past, Present, and Future Use of Phosphorus in Chinese Agriculture and its Influence on Phosphorus Losses. Ambio 44 (Suppl. 2), 274–285. doi:10.1007/s13280-015-0633-0

CrossRef Full Text | Google Scholar

Liu, X., Sheng, H., Jiang, S., Yuan, Z., Zhang, C., and Elser, J. J. (2016). Intensification of Phosphorus Cycling in China since the 1600s. Proc. Natl. Acad. Sci. USA 113 (10), 2609–2614. doi:10.1073/pnas.1519554113

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, Y., Villalba, G., Ayres, R. U., and Schroder, H. (2008). Global Phosphorus Flows and Environmental Impacts from a Consumption Perspective. J. Ind. Ecol. 12 (2), 229–247. doi:10.1111/j.1530-9290.2008.00025.x

CrossRef Full Text | Google Scholar

Long, H., and Qu, Y. (2018). Land Use Transitions and Land Management: A Mutual Feedback Perspective. Land Use Policy 74, 111–120. doi:10.1016/j.landusepol.2017.03.021

CrossRef Full Text | Google Scholar

McDowell, R. W., and Sharpley, A. N. (2004). Variation of Phosphorus Leached from Pennsylvanian Soils Amended with Manures, Composts or Inorganic Fertilizer. Agric. Ecosyst. Environ. 102, 17–27. doi:10.1016/j.agee.2003.07.002

CrossRef Full Text | Google Scholar

Mishima, S., Endo, A., and Kohyama, K. (2010). Recent Trends in Phosphate Balance Nationally and by Region in Japan. Nutr. Cycl Agroecosyst 86 (1), 69–77. doi:10.1007/s10705-009-9274-7

CrossRef Full Text | Google Scholar

MNR (2012). Ministry of Natural Resources of the People's Republic of China. Cultivated Land Area Database [Online]. Available at: http://www.mnr.gov.cn/sj/sjfw/td_19262/ (Accessed January 15, 2021).

Google Scholar

NBS (2016). National Bureau of Statistics. Production and Consumption of P Fertilizer in China [Online]. Available at: http://www.stats.gov.cn/ (Accessed January 15, 2021).

Google Scholar

Pavinato, P. S., Cherubin, M. R., Soltangheisi, A., Rocha, G. C., Chadwick, D. R., and Jones, D. L. (2020). Revealing Soil Legacy Phosphorus to Promote Sustainable Agriculture in Brazil. Sci. Rep. 10, 15615. doi:10.1038/s41598-020-72302-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Qiu, J. (2010). Phosphate Fertilizer Warning for China. Nature. doi:10.1038/news.2010.498

CrossRef Full Text | Google Scholar

Sattari, S. Z., Bouwman, A. F., Giller, K. E., and van Ittersum, M. K. (2012). Residual Soil Phosphorus as the Missing Piece in the Global Phosphorus Crisis Puzzle. Proc. Natl. Acad. Sci. 109 (16), 6348–6353. doi:10.1073/pnas.1113675109

PubMed Abstract | CrossRef Full Text | Google Scholar

Sattari, S. Z., van Ittersum, M. K., Bouwman, A. F., Smit, A. L., and Janssen, B. H. (2014a). Crop Yield Response to Soil Fertility and N, P, K Inputs in Different Environments: Testing and Improving the QUEFTS Model. Field Crops Res. 157, 35–46. doi:10.1016/j.fcr.2013.12.005

CrossRef Full Text | Google Scholar

Sattari, S. Z., van Ittersum, M. K., Giller, K. E., Zhang, F., and Bouwman, A. F. (2014b). Key Role of China and its Agriculture in Global Sustainable Phosphorus Management. Environ. Res. Lett. 9 (5). 054003. doi:10.1088/1748-9326/9/5/054003

CrossRef Full Text | Google Scholar

Sharpley, A., Foy, B., and Withers, P. (2000). Practical and Innovative Measures for the Control of Agricultural Phosphorus Losses to Water: an Overview. J. Environ. Qual. 29, 1–9. doi:10.2134/jeq2000.00472425002900010001x

CrossRef Full Text | Google Scholar

Steen, I. (1998). Phosphorus Availability in the 21st century: Management of a Non-renewable Resource. Phosphorus and Potassium 217, 25–31.

Google Scholar

Syers, J. K., Johnston, A. E., and Curtin, D. (2008). Efficiency of Soil and Fertilizer Phosphorus Use. Rome: FAO.

Google Scholar

USGS (2015). United States Geological Survey. Mineral Commodity Summaries. Available at: http://minerals.usgs.gov./minerals/pubs/commodity/phosphate_rock/mcs-2015-phosp.pdf (Accessed June 15, 2021).

Google Scholar

USGS (2016). United States Geological Survey. Mineral Commodity Summaries. Available at: http://minerals.usgs.gov./minerals/pubs/commodity/phosphate_rock/mcs-2016-phosp.pdf (Accessed June 15, 2021).

Google Scholar

van Keulen, H., Aarts, H. F. M., Habekotté, B., van der Meer, H. G., and Spiertz, J. H. J. (2000). Soil-plant-animal Relations in Nutrient Cycling: the Case of Dairy Farming System 'De Marke'. Eur. J. Agron. 13, 245–261. doi:10.1016/s1161-0301(00)00077-0

CrossRef Full Text | Google Scholar

Verloop, J., Oenema, J., Burgers, S. L. G., Aarts, H. F. M., and van Keulen, H. (2010). P-equilibrium Fertilization in an Intensive Dairy Farming System: Effects on Soil-P Status, Crop Yield and P Leaching. Nutr. Cycl Agroecosyst 87 (3), 369–382. doi:10.1007/s10705-010-9344-x

CrossRef Full Text | Google Scholar

Wolf, J., Wit, C. T., Janssen, B. H., and Lathwell, D. J. (1987). Modeling Long‐Term Crop Response to Fertilizer Phosphorus. I. The Model 1. Agron.j. 79, 445–451. doi:10.2134/agronj1987.00021962007900030008x

CrossRef Full Text | Google Scholar

Yu, W. J., Li, G. H., Hartmann, T., Xu, M. G., Yang, X. Y., Li, H. G., et al. (2021). Development of a Novel Model of Soil Legacy P Assessment for Calcareous and Acidic Soils. Front. Environ. Sci. 8, 1–12. doi:10.3389/fenvs.2020.621833

CrossRef Full Text | Google Scholar

Zhang, F., Chen, X., and Vitousek, P. (2013). An experiment for the World. Nature 497, 33–35. doi:10.1038/497033a

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, W., Cao, G., Li, X., Zhang, H., Wang, C., Liu, Q., et al. (2016). Closing Yield Gaps in China by Empowering Smallholder Farmers. Nature 537, 671–674. doi:10.1038/nature19368

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, W., Ma, W., Ji, Y., Fan, M., Oenema, O., and Zhang, F. (2008). Efficiency, Economics, and Environmental Implications of Phosphorus Resource Use and the Fertilizer Industry in China. Nutr. Cycl Agroecosyst 80 (2), 131–144. doi:10.1007/s10705-007-9126-2

CrossRef Full Text | Google Scholar

Zhang, W., Tang, X., Feng, X., Wang, E., Li, H., Shen, J., et al. (2019). Management Strategies to Optimize Soil Phosphorus Utilization and Alleviate Environmental Risk in China. J. Environ. Qual. 48 (5), 1167–1175. doi:10.2134/jeq2019.02.0054

PubMed Abstract | CrossRef Full Text | Google Scholar