Optimization of management strategies for reducing nitrogen loading in China

Journal Pre-proofs Optimization of management strategies for reducing nitrogen loading in China Shanshan Hua, Xin He, Chunmiao Zheng PII: DOI: Reference:

S0048-9697(19)34611-X https://doi.org/10.1016/j.scitotenv.2019.134620 STOTEN 134620

To appear in:

Science of the Total Environment

Received Date: Revised Date: Accepted Date:

24 July 2019 19 September 2019 21 September 2019

Please cite this article as: S. Hua, X. He, C. Zheng, Optimization of management strategies for reducing nitrogen loading in China, Science of the Total Environment (2019), doi: https://doi.org/10.1016/j.scitotenv.2019.134620

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2019 Elsevier B.V. All rights reserved.

Oimization of management strategies for reducing nitrogen loading in China Shanshan Hua a, b, Xin He c, Chunmiao Zheng b, a School b



of Environment and Energy, Peking University, Shenzhen 518055, China

State Environmental Protection Key Laboratory of Integrated Surface Water-

Groundwater Pollution Control, School of Environmental Science and Engineering, Southern University of Science and Technology, Shenzhen 518055, China c Department

of Water Resources, China Institute of Water Resources and Hydropower

Research, Beijing 100038, China

Corresponding author, Chunmiao Zheng, School of Environmental Science and Engineering, Southern University of Science and Technology, Shenzhen 518055, China. Tel.: +86 0755-88010020; E-mail address: [email protected] 1

Abstract: Terrestrial Nitrogen (N) loading in the environment has seen a steady increase over the past several decades as a result of more intensive anthropogenic activities. Quantifying N loading for an extended period is important for effective N management. In this study, a statistical model is constructed to describe the relationship between N loading and anthropogenic activities at watershed scale for 211 watersheds covering the entire land area of China. Subsequently, a portfolio optimization model is used to optimize the future management efforts of the long-term N loading. Our results show that N loading in China due to anthropogenic activities has increased significantly over the past 60 years (1949-2010), with the rate of increase at approximately 1 Tg N/year. When designing future N loading management strategies, the next 30 years is divided into three temporal stages and assume that the total amount of expenditure is fixed. The results of portfolio optimization analysis show that the best allocations of management efforts (e.g. capital investments, making new policies, improving technology, or alike) among three temporal stages are 28.55% (2021-2030), 71.45% (2031-2040) and 0 (2041-2050). Furthermore, it is suggested that the future population growth scenario has the largest influence on the results of the portfolio optimization analysis. Keywords: Nitrogen loading; Watershed management; Water pollution; Portfolio optimization; China

2

1. Introduction Anthropogenic activities have greatly altered the nutrient cycle at catchment scale with population and agricultural land use, and therefore accelerated the rate of nutrient fixation across landscapes as well as the discharge of nutrients to the aquatic environment (Smith, 1998; Vitousek et al., 1997; Zhou et al., 2000; Santos et al., 2017). Due to anthropogenic activities, the rate of biologically available nutrients entering the terrestrial environment have more than doubled compared to the level in the preindustrial era (Galloway et al., 2004). In the past several decades, the world has witnessed degradation of aquatic ecosystems and impaired water quality for drinking, agriculture, industry, recreation etc. (Carpenter et al., 1998, Peterson et al., 2001, Zhang et al. 2015). One of the key nutrients that strongly influences the degradation of aquatic and terrestrial ecosystems is Nitrogen (N). The environmental degradation caused by excessive N discharge has become a serious problem in China (Mcdonald et al., 2016). Surface water quality in most Chinese watersheds is substandard due to human wastewater discharge, and agricultural run-offs (e.g. fertilizers, manures, and pesticides) (Yang et al., 2016). According to Kahrl et al., (2010), China is the world’s largest consumer of N fertilizers, accounting for nearly one-third of the global consumption, where the massive input of N has led to serious environmental problems (Galloway et al., 2004; Galloway, 2008). Therefore, in order to restore the watershed quality in China, it is important to control the N discharge to water bodies. 3

Several management efforts have been taken in the past decades in China to minimize N loading in the environment. For instance, the capacity of sewage treatment plants has increased significantly in recent years, where the urban sewage treatment rate has increased from less than 7% in 1985 to 89.3% in 2013 (CSY, 2016). However, those efforts do not match the rate of increase in the N loading, and thus do not successfully prevent deteriorating watershed quality. In order to deal with this situation, more effective measures must be introduced and added to the current N loading management practices. Due to the fact that the majority of the N loading is caused by anthropogenic activities, a deeper understanding of the link between the two factors is urgently needed. Various types of models have been used to build the relationship between N loading and anthropogenic activities in numerous watersheds (He et al., 2011; Mockler et al., 2017). These models range from simple statistical models, such as linear regression models and Mass Balance Models (MBM) (Mouri et al., 2010), to complex process-based models, such as the SPARROW model (Li et al., 2015), Integrated Nitrogen in Catchments models (INCA) (Pathak et al., 2018), the Storm Water Management Model (SWMM) (Bo et al., 2013), and Visualizing Ecosystems for Land Management Assessments (VELMA) (Abdelnour et al., 2013). Several studies have estimated the N loading in watersheds using these models and recommended management policies accordingly (Wu et al., 2009; Zhang et al., 2009; Gu et al., 2015; Yu et al. 2019). However, many of these models are primarily based on describing 4

physical or chemical processes that lead to watershed pollution (Beusen et al., 2015), and thus would need detailed spatial information regarding these processes. However, such information is still insufficient in most watersheds in China, especially for longer time spans (He et al., 2011). Therefore, a relatively simple model to estimate the N loading for a long time series is needed. Current watershed management goals that aim to reduce N loading can be achieved through cultural, policy and/or technological measures (Bernhardt et al., 2008). The historical data on the N loading can guide these efforts, but historical trends may not persist into the future. A key issue to manage future N loading is to incorporate uncertainty (He et al., 2011). A portfolio optimization model is a management tool used to maximize expected return or minimize overall risks for the summation of assets that give consideration to uncertainty (Sharpe, 1970). Application of such a model has moved from the financial sector to environmental management, like the analysis of optimal species and genetic diversity (Figge, 2004; Koellner and Schmitz, 2006). Some studies used the model to select seed sources (Crowe, and Parker, 2008) and protected lands (Ando and Mallory, 2012, Liang et al., 2017). However, these studies have not combined long-term N loading with the portfolio optimization model. In this study, portfolio optimization modeling is attempted to better allocate management efforts (e.g. capital investments, making new policies, improving technology, etc.) for mitigating the N loading to the environment. The objectives of this study are: (1) to construct a statistical model that connects 5

the anthropogenic activities and the N loading for all watersheds in China, (2) to estimate past and future N loading in China at the national scale; and (3) to estimate the optimal allocation of management efforts to deal with excessive N loading in the next 30 years.

2. Methods and data sources Population growth and human induced land use change was used to represent anthropogenic activities in our study. A mass balance model first was used to calculate N loading in 211 catchments in China. Then, a statistical model was constructed by relating the above-mentioned anthropogenic factors and the calculated N loading at each catchment. Lastly, a portfolio analysis was performed to determine the weights for optimal allocation of management efforts in the next 30 years. 2.1 Calculation of nitrogen loading for selected catchments A mass balance model was used to estimate the N loading in China at watersheds (data available at http://www.resdc.cn/Default.aspx) scale (Fig 1) in 2005 and 2010, respectively. The equation is shown below: 𝑁𝑙𝑜𝑎𝑑𝑖𝑛𝑔 = 𝑁𝑖𝑛𝑝𝑢𝑡 ― 𝑁𝑜𝑢𝑡𝑝𝑢𝑡 = 𝑁𝑓𝑒𝑟𝑡 + 𝑁𝑓𝑖𝑥𝑎 + 𝑁𝑟𝑒𝑠𝑖 + 𝑁𝑎𝑛𝑖𝑚 + 𝑁𝑑𝑒𝑝𝑜 ― 𝑁𝑐𝑟𝑜𝑝 ― 𝑁𝑟𝑢𝑛𝑜

(1)

The N input to the catchments are from fertilization, atmospheric deposition, fixation and their subsequent redistribution spatially across the landmass through food and feed trade (Anderson et al., 1997; Diamond, 1995; Green, 2004). The output of the

6

MBM model is the N due to denitrification losses on landmass, in aquatic systems and for potential sequestration. Therefore, N loading in this study is interpreted as the N leftover in the catchment after internal movements and fades (growth, harvesting and transport of crops, feed and forage, production and transport of animal products, denitrification, ammonia volatilization, emission from livestock and soils, removal by sewage treatment systems, and biomass burning (Galloway et al. 2004)) in each catchment and includes all forms of N. It is not distinguished among hydrological compartments, and these may include surface water, groundwater, pore water, etc. The details of the parameters in the equation are described below. Industrial nitrogen fertilizer application (Nfert): net fertilizer inputs were computed accounting for fertilizer application and volatilization losses from field-application. Fertilizer loads were considered for industrially manufactured fertilizers only, separately dealing with animal manure fertilizers. N fertilizer consumption totals and volatilization losses for 2010 and 2005 were taken from the FAOSTAT Statistical Databases (FAO 2012). These values were evenly distributed among cropland. Nitrogen fixation (Nfixa): including nitrogen fixation by undisturbed vegetation and cultivated land. Nitrogen fixation by undisturbed vegetation was developed from potential rates defined for undisturbed vegetation classes (Cleveland et al. 1999). The crop fixation was based on country-level estimates of crop area (FAO 2010) and fixation rates for the predominant N fixing crops (seed legumes, open fields, rice, and sugar cane) (Smil 1999). Soil emissions of nitrogen in the form of ammonia (NH3) and 7

nitrogen losses from biomass burning were deducted from the fixation estimates. The urban areas were assigned a fixation value of 0. Nitrogenous compounds generated by residential living (Nresi): human waste loads were calculated using total nitrogen intake estimates and applied to a spatially distributed global data set of urban and rural populations. In this portion, human loads of nitrogen from rural and urban populations with no access to sewer systems were considered directly to the land. For populations with access to sewage treatment (OECD 2015), the decrement of nitrogen was calculated and deducted from the total nitrogenous compounds generated by the residences in the corresponding area. Animal intake and excretion of nitrogen (Nanim): animal waste load was derived by applying species specific N emission units based on animal intake minus ammonia volatilization. This yielded the effective nitrogen loads to the land per livestock head: dairy cattle, non-dairy cattle, swine, sheep and goats, horses, caribou, camels, and water buffalo (NRC 1985; Smil 1999). The animal waste load was determined from the province-level production statistics from the China Statistical Yearbook for 2010. These values were evenly distributed among agriculture land. Atmospheric deposition (Ndepo): Atmospheric inputs for nitrogen were based on modeled estimates of total (wet and dry) inorganic (NOy plus NHx) deposition estimates from Dentener and Crutzen (1994). Nitrogen in crops (Ncrop): the nitrogen content in crops was determined from the country-level production statistics from the FAOSTAT database for 2010. The percent 8

of nitrogen in major crop classes was taken from estimates of nitrogen in harvested crops. This was done for seven harvested types globally (cereals, legumes, sugar crops, roots and tubers, vegetables and fruits, forages, and other crops). A part of the nitrogen in crops was consumed by the residential area and animals. Nitrogen loading to the aquatic system by runoff (Nruno): the spatially distributed nitrogen loadings to the landmass that are potentially mobilizable into the aquatic system by runoff. It was calculated by watershed runoff, precipitation and nitrogen delivery coefficients on the landmass. Nitrogen delivery coefficients on the landmass (coE) were calculated by residence time for nitrogen on the landmass (τ) which was estimated as the sum of the mean annual volume of water stored as soil moisture and shallow groundwater, (from 1-year estimates) divided by mean annual basin runoff (Fekete et al., 2002). 𝐶0

𝐶𝑇 = coE = 𝑒

( ― 𝑘𝑇𝜏)

which kT (year-1) is the Rate Constant of the First Order Reaction (KT=1 year-1 in this paper). C0 is influent concentration; Ct is effluent concentration; 𝑁𝑟𝑢𝑛𝑜 = 𝑁𝑁𝐸𝑇 ∗ coE ∗ runoff/precipitation

2.2 Statistical model between anthropogenic activities and N loading The goal of this paper is to calculate the N loading in China over a long period of time. In principle, the best solution is to use the MBM to calculate the N loading for the entire time duration (1949-2050). However, N loading estimated using MBM requires 9

detailed maps of numerous types of data, which are not available for most layers, especially before 2000. To overcome this problem, a regression model is used to compensate for the years when the MBM model cannot be used due to lack of data. Multiple linear regression modeling has been widely used to predict sediment and nutrient loading over time based on anthropogenic activity. (Beusen et al., 2005; McDonald et al., 2016). Accordingly, the second step in our analysis is to construct a multiple linear regression based on the N loading which calculated by MBM and anthropogenic factors. This regression model is used to estimate N loading of China since 1949 based on watershed units, this N loading is called Watershed N-loading (WNL) and the model as the "statistical model". The details of this statistical model are described below. Four anthropogenic factors are considered in the model: (1) crop area× fertilizer N application rate, (2) agricultural area × manure N application rate, (3) population, and (4) total catchment area. The difference between agricultural and crop land is that only manure is applied on agricultural land, while both fertilizers and manure are applied on crop land. The equation of the statistical model is shown below:

Watershed N-loading (WNL) = a*Xcrop area× + b*Xagricultural area×

manure N application rate

fertilizer N application rate

+ c*Xpopulation

+ d*Xtotal catchment area+e

(2)

Where a, b, c and d are the regression coefficients for each parameter; e is the intercept of the equation. The multiple linear regression was performed to find out parameter a

10

to e using all results from the 211 catchments in 2010. To improve normality, all variables were standardized by log10-transformed. To validate our statistical model, the estimate of N loading from the model in 2005 was compared to that predicted by the Mass Balance Model (McDonald et al., 2016). Afterwards, the "validated" statistical model was used to predict N-loading for all of China as a whole for the past years (from 1949 to 2010) and the next 30 years (from 2020 to 2050). This approach assumes that the relationship between the N loading and the variables in the statistical model have been constant over time. When the N loading is predicted in the future years, only the population and N application rate change with time, while other impacting factors, mainly agricultural area size, remain constant. The trend for how populations grows in the future was projected by World Population Projections 2017 revision: S1 “No change variant” scenarios, both fertility and mortality are kept constant; S2 “Low variant” scenarios, low-fertility assumption; S3 “Medium variant” scenarios, medium-fertility assumption, and S4 “High variant” scenarios, high-fertility assumption (United Nations, 2017). The data sources of N loading predictions in China are described as following. The historical data on cropland and agricultural land areas were derived from the China Agriculture Yearbook (1980–1999) and World Data Atlas (Knoema, 2017) (Fig. 2a). Historical population data for China was from the China Statistical Yearbook (Fig. 2b). Raw data used to calculate fertilizer N application rate and manure N application rate 11

were derived from China Agriculture Yearbook, World Data Atlas and FAOSTAT statistical databases (FAO, 2012) (Fig. 3). The future populations were from World Population Projections 2017 revision (Fig. 4).

2.3 Portfolio optimization analysis The portfolio optimization model is an investing model where the investor attempts to take minimal level of market risk to capture maximum-level returns for a given portfolio of investments. In this study, portfolio optimization analysis is used to help planners make strategic management allocation that manage risk of uncertainty more effectively to increase the value of future Watershed N loading index (WNLI). The results of Portfolio analysis give exact portfolio weightings for the three stages. Output from the portfolio analysis can either be concrete or abstract management efforts. For instance, Ando and Mallory (2012) used capital investments as their major management effort, whereas Crowe and Parker (2008) used strategy as their major management effort to protect forests of white spruce. However, one can also consider the combination of capital investment, policy making, technology improvement and alike as an abstract form of management effort. In the present study, the abstract form is used since at this stage it is still too early to separate each specific management effort. To build the portfolio optimization model, a factor was created with the name Watershed N loading Index (WNLI). The equation to calculate the WNLI is shown below: 12

𝑊𝑁𝐿𝐼 = (𝑊𝑁𝐿𝑚𝑎𝑥 ― 𝑊𝑁𝐿𝑖)/(𝑊𝑁𝐿𝑚𝑎𝑥 ― 𝑊𝑁𝐿𝑚𝑖𝑛)

(3)

where i, min and max represent the values, minimum values and maximum values of the N loading of China (WNL) in the dataset under future scenarios. This factor is an original concept that has not been used elsewhere. WNLI is a normalized value, ranging between 0 and 1. It represents the watershed quality in this paper. When WNLI is higher, it means the catchment environment quality is better. WNLI will be modeled for the next 30 years and divided into three stages (first stage: year 2021-2030; second stage: year 2031-2040; third stage: year 2041-2050), including four plausible outcomes of population projection scenarios (future scenarios) which has been described in the section 2.2 (World Population Projections 2017 revision). In the portfolio optimization model, corresponding proxies for the variables are: (1) the expected return of each stage equals the value of WNLI under the future scenarios; (2) the expected variance of each stage equals its variance in return over all considered future scenarios; and (3) the expected covariance of each stage’s WNLI with other stages equals its covariance in return across all future scenarios (Hua et al., 2015). The objective function of the portfolio optimization model is to minimize the expected variance and covariance of the portfolio of projected WNLI for each stage over the future scenarios evaluated. In the analyses, the return is defined as WNLI and the expected WNLI in each stage is defined as below E[𝑊𝑁𝐿𝐼𝑖] =

∑ 𝑃 × 𝑊𝑁𝐿𝐼 𝑗

𝑖𝑗

for all future scenarios 𝑗

𝐽

13

(4)

where the expected WNLI in stage i is the sum of the probabilities (P) of each future scenario times the realized WNLI in stage i for future scenario j; the sum of the probability (Pj) of each future scenario to occur is 1. Ando et al. (2012) explain this approach in more details. In this study, it is difficult to say how likely each of the future scenarios is to occur.

Two sample probability distributions are considered to illustrate

the probabilities (P) of four future scenarios: “historic trend” has weights that assume the “no change” scenario has the highest probability to occur; the other called “uniform,” assumes each future scenario is equally likely to occur.

3. Results 3.1 Distribution of the nitrogen loading in China N loading is denoted as annual weight per hectare (kg/ha/year), and is calculated for 2005 and 2010, respectively, due to data availability. The distribution of N loading at the watershed scale is seen for year 2005 and 2010 (Fig. 5). Comparing these two years, it can be seen there are not obvious differences of N loading trends in terms of spatial distribution, ignoring the change of the N loading value. For both years, it can be seen that the watersheds with the highest N occur in the North China Plain, the Middle-lower Yangtze Plain and the region of southeastern China, as expected. These regions are featured with cropland area or high population density, which are considered main sources for N loading. Moreover, the distribution of N loading has an increasing trend from the northwest to southeast. Human and natural factors such as population distribution, agricultural area, and vegetation characteristics are the main 14

reasons for this trend.

3.2 Results and validation of the statistic model Data from 2010 were used for model calibration for parameters a-e in Equation (2). The regression coefficients for the parameters in the statistical model for N loadings in 2010 are listed in Table 2 (R2 = 0.848, P < 0.005). Data from 2005 were used for validation based on parameters obtained in Table 2. The correlation between the statistical model predictions and the Mass Balance Model in each of the watersheds is high (R = 0.942) and follows the 1:1 line (Fig. 6). This means this statistical model can be used to predict the past and future N loading of China.

3.3 Simulation of nitrogen loading in China in the past and future The N loading in China from 1949 to 2010, as calculated by the statistical model, Equation (2), is shown in Fig. 7. In this process, it was assumed that natural factors are fixed and only the anthropogenic factors, namely population growth and agricultural land use change, were used as independent variables. N loading in China had large anomalies between year 1950 and 1961, which is consistent with historical events. The two major events that caused this to happen were: First, China began to heavily use fertilizers starting from 1955, and second there was a severe famine in China during 1959 and 1961. After 1962, N loading is seen with steady growth. Over the past 60 years, N loading to landmass in China increased from 22.2 Tg (109 kg) in 1950 to 80.8 Tg in 2010, an increase of almost four times. The annual growth rate of N loading is 15

approximately 1 Tg N/year. The total N loading to landmass was calculated to 3914 Tg from 1949 to 2010. In general, the inter-annual variation of N loading is positive. The fluctuation of N loading decreased with the time. In order to effectively manage the N loading in the future, the variation of N loading from 2020 to 2050 under 4 different management scenarios (Fig. 8) was predicted. In the four scenarios, the difference of N loading variant in the early period (2020-2024) is not very obvious. However, the variation of the change is very significant in the later period. Under S1, S2, S3, which are no change, low variant, and medium variant scenarios respectively, the N loading first increased and then decreased. While N loading has continued to rise in high variant scenarios (S4). These results are consistent with the demographic trends in the four scenarios (shown in Fig. 4).

3.4 Portfolio optimization analysis for reducing nitrogen loading in China WNLI was calculated for every decade starting from 2020 and for the four population growth scenarios (Table 3). The results of the portfolio optimization analysis were obtained by using the assumed probabilities of four future scenarios, as in Table 3, and the calculated WNLI values across the three stages as input. In Fig.9, the horizontal axis is standard deviation of WNLI (risk of uncertainty), and the vertical axis is the expected return of WNLI (the value of WNLI in the managed stages). These two curves are called the efficient frontiers. The efficient frontier is the portfolio set of three temporal stages with the best expected return achievable for a given level of risk. The efficient frontier describes the relationship between the return 16

that can be expected from a portfolio and the riskiness (volatility) of the portfolio. The efficient frontier gives the best return that can be expected for a given level of risk or the lowest level of risk needed to achieve a given expected rate of return. In this paper, two probability distributions of population growth scenarios are described, they are “uniform” probability and “Historic Trend” probability. The dotted curve and solid curve represent the efficient frontier for “uniform” probability and “Historic Trend” probability, respectively. The points A to E (F to J) (risk and expected have been marked in the figure) are the sample points to explain the efficient portfolio for three temporal stages. Table 4 shows the detailed information of sample points from Fig. 9. The table includes the Uncertainty (standard deviation of WNLI), E[R] (expected return of WNLI) and Portfolio weights (the proportion of an investment portfolio to reach the expected). The Portfolio weights get from the efficient frontiers. This table gives the detailed portfolio weightings for each managed stage to maximize the WNLI in the managed stages under a certain uncertainty. Under the uniform condition, from point A to point E, the expected value of WNLI increased by 14% accompanied by the risk increased by 49%. Similarly, under the historic trend condition, the expected value of WNLI increased by 12% with the corresponding risk increased by 87%. For both probability distributions of population growth scenarios, from low point to high point on efficient frontiers, the maximum of WNLI increase is 14% which is seriously unbalanced between the expected and uncertainty. Therefore, it is advised to choose the portfolio weightings for each management stage with the minimized risk. The best 17

portfolio weightings of management efforts among three temporal stages are 28.55% (2021-2030), 71.45% (2031-2040) and 0 (2041-2050). Furthermore, the probability distribution of future scenarios affects the weight and position of an efficient portfolio frontier in risk/expected benefit analysis. Among the same uncertainty (standard deviation of WNLI), “historical trend” the probability distribution scenario holds much higher expected return of WNLI compared to the “uniform” scenario. The probabilities of future scenarios also influence the choices of portfolio weights for each stage. For instance, considering points J and E in the two probability distributions, one should put all effort into the third stage for the “historic trend” probability distribution (point J); whereas in the “uniform” probability distribution, all effort should be put into the first stage (point E). This way, they give the maximum expected value of WNLI over a maximum risk.

4. Discussion The goal of this paper is to calculate the N loading in China since its establishment in 1949. In order to achieve this goal, MBM was applied to calculate the distribution of N loading. This result shows good agreement with the result of Green (2004). For both results, it can be seen that the watersheds with the highest N occur in the North China Plain, and the region of southeastern China. However, spatially explicit estimates with MBM require detailed maps of many data over 62 years, which are not available for most layers. Accordingly, using watersheds as the study unit, this study constructed a multiple linear regression equation to estimate N loading for China at the multi-decadal 18

time scale. It was calculated that N loading to landmass in China increased from 22.2 Tg in 1950 to 80.8 Tg in 2010 (Cui et al. (2013) estimated a total 79 Tg N of China in 2010). The surplus N was distributed in terrestrial and aquatic ecosystems which led to severe impacts on soil and water quality. The following information can be obtained according to the regression coefficients in Table 2. A 10% increase in population leads to a 3.27% increase in N loading (holding all other variables constant). Similarly, a 10% increase in crop area × fertilizer N application rate or agricultural area × manure N application rate leads to a 0.39% or 0.48% increase in N loading. It is clear that N loading is more sensitive to population. If N loading in watersheds need to be reduced rapidly in a short amount of time, population-related human dietary choices will need more attention by publicizing a low-nitrogen diet. However, considering the long-term N management, it is not enough to only manage population and human dietary choices. Decision makers can also restrict the rate of N manure and N fertilizer applications. N can be greatly reduced if nutrients are applied at rates that match their uptake by crops, and if fertilizers are applied when crops are growing rapidly (Zhu and Chen, 2002, Kanter et al., 2015). Overall, the development and spread of N-efficient technologies should be a high priority for planners. It is known that the management of N loading is complex. N is widely and unevenly distributed in terrestrial and aquatic ecosystems rather than directly into waste treatment systems. Prediction results show China still faces excessive emission of nitrogen in the future. These make the treatment of surplus N more difficult. Through 19

the portfolio optimization analysis, management plans can be implemented to achieve better cost-benefit ratios. The results of portfolio optimization analyses give the percentage of capital investments, policies, and technological improvements for each temporal stage.

Decision makers can choose the best portfolio of management efforts

on the efficient frontier over an uncertainty about the likelihood of population change. In this study, it will be a better solution to allocate the management effort in each management stage when the risk of uncertainty is smallest. The best allocations of management efforts among three temporal stages are 28.55% (2021-2030), 71.45% (2031-2040) and 0 (2041-2050). If management efforts are constrained not to exceed some fixed resource from 2020 to 2050 (fixed money, reasonable policy, etc.), 28.55% in the first stage means the percentage of management efforts in the first stage within this total fixed resource is 28.55%, The percentage of management efforts in the second stage within this total fixed resource is 71.45%. The third temporal stage comes to zero. However, this does not suggest inaction in the third stage, it means we should manage excessive N as early as possible. New investments of infrastructure may not be needed, but maintenance work introduced in the first and second stage should continue. The management efforts include capital investments, new policies, improving technology, etc.. The capital investments in the first and second temporal stages are used to construct wastewater treatment plants, improve operational efficiency of sewage treatment, manage contaminated rivers and soils, etc.

These new management efforts

can help control pollutant emission, decrease N based fertilizer application, and use 20

crop straw more efficiently. In the third temporal stage, it is unnecessary to reinvest in the construction of new wastewater treatment plants or develop new policies. The wastewater treatment plant and management policies introduced in the first and second temporal stages, if they can be implemented according to plan, can already meet the N management in the third temporal stage. So, the percentage of the third temporal stage comes to zero. Furthermore, management of the future N loading, has the following characteristics: First, population change poses significant uncertainty regarding the future N quality for the management activity being targeted. Second, action to adapt to population change should be taken long before that uncertainty is resolved. Third, to take reasonable action to make the best N quality in managed stages at a given risk of uncertainty. In view of these characteristics, MPT can help planners make strategic management allocation that manage risk of uncertainty more effectively to increase the value of WNLI. Portfolio optimization analysis is used in this study to cope with future uncertainties which were heavily impacted by population. There are also some other studies (Van Drecht et al., 2003; More et al., 2013) with objectives to analyze the uncertainties of the result of the N loading, but they focused primarily on the historical estimates or model parameters. Studies that involved future prospects and management opinions are relatively rare. Our study made the first attempt to give future scenario analysis regarding N loading at national scale with uncertainty taken into consideration.

21

5. Conclusions In the present study, the N loading at 211 watersheds in China was calculated using a mass balance method. The result showed that the watersheds with the highest N occurred in the North China Plain, the Middle-lower Yangtze Plain and regions along the coast of southeastern China. The connection between the anthropogenic activities (population and agricultural land uses) and N loading at watershed scale was constructed by a statistical model. The model parameters were first calibrated using data from 2010 and then validated using data from 2005. The validation showed good agreement between the mass balance model and statistic model when using the population and the agricultural land use as the independent variables. Based on the statistical model, the past (1949–2010) and future (2020-2050) N loading in China at national scale were analyzed. Over the past 60 years, N loading to landmass in China increased from 22.2 Tg in 1950 to 80.8 Tg in 2010, an increase of almost four times with an annual growth of 1 Tg N/year. The total N to landmass was calculated to be 3914 Tg over the past 60 years. Moreover, N loading can reach up to 88 Tg (S4 “High variant” scenarios) by 2050. Subsequently, a portfolio optimization model was used to distribute management efforts according to three temporal stages (2021-2030, 20312040, 2041-2050). Our results found that the probability distribution of future population projection scenarios affected the weight and position of the efficient portfolio frontiers. For the same uncertainty level, the “historical trend” probability distribution scenario exhibited much higher expected return of WNLI compared to the 22

“uniform” scenario. Additionally, our results showed that the best strategy to allocate the management efforts for the three temporal stages in the future is 28.55% (first stage), 71.45% (second stage) and 0 (third stage) for the both probability distribution scenarios in order to achieve the best expected return of WNLI.

Acknowledgements This work was supported by the National Key R&D program of China (project number 2016YFC0402806). Additional support was provided by the National Natural Science Foundation of China (grant no. 41890852 and 41701042).

References Abdelnour, A., Mckane, R. B., Stieglitz, M., Pan, F., & Cheng, Y. (2013). Effects of harvest on carbon and nitrogen dynamics in a pacific northwest forest catchment. Water Resources Research, 49(3), 1292-1313. Anderson, I. C., Tobias, C. R., Neikirk, B. B., & Wetzel, R. L. (1997). Development of a process-based nitrogen mass balance model for a virginia (usa) spartina alterniflora salt marsh: implications for net din flux. Marine Ecology Progress, 159(8), 13-27. Ando, A. W., & Mallory, M. L. (2012). Optimal portfolio design to reduce climaterelated conservation uncertainty in the Prairie Pothole Region. Proceedings of the National Academy of Sciences, 109(17), 6484-6489. 23

Bernhardt, E. S. , Band, L. E. , Walsh, C. J. , & Berke, P. E. . (2008). Understanding, managing, and minimizing urban impacts on surface water nitrogen loading. Annals of the New York Academy of Sciences, 1134(1), 61-96. Beusen, A. H. W. , Dekkers, A. L. M. , Bouwman, A. F. , Ludwig, W. , & Harrison, J. . (2005). Estimation of global river transport of sediments and associated particulate c, n, and p. Global Biogeochemical Cycles, 19. Beusen, A. H. W. , Van, B. L. P. H. , Bouwman, A. F. , Mogollón J. M., & Middelburg, J. J. . (2015). Coupling global models for hydrology and nutrient loading to simulate nitrogen and phosphorus retention in surface water – description of image–gnm and analysis of performance. Geoscientific Model Development, 8(12), 4045-4067. Bo, L., Dong-Guo, S., Min, S., & Yu-Long, L. (2013). Simulating the Effect of Reducing the Non-point Source Pollution by;Buffer Zone with SWMM. International Conference on Measuring Technology & Mechatronics Automation. Carpenter, S. R., Caraco, N. F., Correll, D. L., Howarth, R. W., Sharpley, A. N., & Smith, V. H. (1998). Non-point pollution of surface waters with phosphorus and nitrogen. Ecological Applications, 8(3), 559-568. CAY, 2011. China Agriculture Yearbook. China Statistical Publishing House, Beijing. CESB, (2010). China Environmental Status Bulletin. Ministry of Environmental protection of the People’s Republic of China, Beijing. Cleveland, C. C., Townsend, A. R., Schimel, D. S., Fisher, H., Howarth, R. W., & 24

Hedin, L. O., et al. (1999). Global patterns of terrestrial biological nitrogen (n2) fixation in natural ecosystems. Global Biogeochemical Cycles, 13(2), 623-646. Crowe, K. A. , & Parker, W. H. . (2008). Using portfolio theory to guide reforestation and restoration under climate change scenarios. Climatic Change, 89(3-4), 355370. CSY, (2016). China Statistical Yearbook. China Statistical Publishing House, Beijing. Cui S, Shi Y, Groffman PM, Schlesinger WH, Zhu YG (2013). Centennial-scale analysis of the creation and fate of reactive nitrogen in China (1910-2010). Proc Natl Acad Sci USA, 110(6), 2052–2057. Dentener, F. J., & Crutzen, P. J. (1994). A three-dimensional model of the global ammonia cycle. Journal of Atmospheric Chemistry, 19(4), 331-369. Diamond, M. L. (1995). Application of a mass balance model to assess in-place arsenic pollution. Environmental Science & Technology, 29(1), 29. (EDC) Earth Resources Observation Systems (EROS) Data Center (2000). Global Land Cover Characterization. EROS Data Center Distributed Active Archive Center. http://edcdaac.usgs.gov/. Fekete B.M., Vorosmarty C.J. and Grabs W. (2002). High-resolution fields of global runoff combining observed river discharge and simulated water balances. Global Biogeochem. Cy. 16(3): 1042. Figge F (2004) Bio-folio: Applying portfolio theory to biodiversity. Biodivers Conserv. 13:827–849. 25

Food and Agriculture Organization of the United Nations (FAO) (2012). FAOSTAT: FAO Statistical Databases. http://www.fao.org/faostat/en/#data. Galloway, J. N. , Dentener, F. J. , Capone, D. G. , Boyer, E. W. , Howarth, R. W. , & Seitzinger, S. P. , et al. (2004). Nitrogen cycles: past, present, and future. Biogeochemistry, 70(2), 153-226. Galloway, J. N. (2008). An earth-system perspective of the global nitrogen cycle. Nature, 451(7176), 293. Green, P. A., Vörösmarty, C. J., Meybeck, M., Galloway, J. N., Peterson, B. J., & Boyer, E. W. (2004). Pre-industrial and contemporary fluxes of nitrogen through rivers: a global assessment based on typology. Biogeochemistry, 68(1), 71-105. Gu, B. , Ju, X. , Chang, J. , Ge, Y. , & Vitousek, P. M. . (2015). Integrated reactive nitrogen budgets and future trends in china. Proceedings of the National Academy of Sciences, 112(28), 8792-8797. He, B., Kanae, S., Oki, T., Hirabayashi, Y., Yamashiki, Y., & Takara, K. (2011). Assessment of global nitrogen pollution in rivers using an integrated biogeochemical modeling framework. Water research, 45(8), 2573-2586. Hua, S., Liang, J., Zeng, G., Xu, M., Zhang, C., Yuan, Y., ... & Huang, L. (2015). How to manage future groundwater resource of China under climate change and urbanization: An optimal stage investment design from modern portfolio theory. Water Research, 85, 31-37. Kahrl, F., Yunju, L., Roland-Holst, D., Jianchu, X., & Zilberman, D. (2010). Toward 26

sustainable use of nitrogen fertilizers in China. ARE Update, 14(2), 5-7. Kanter, D. R. , Zhang, X. , & Mauzerall, D. L. . (2015). Reducing nitrogen pollution while decreasing farmers’ costs and increasing fertilizer industry profits. Journal of Environment Quality, 44(2), 325. Klein Goldewijk, K., A. Beusen, J.Doelman and E. Stehfest, New anthropogenic land use estimates for the Holocene; HYDE 3.2, in prep. ftp://ftp.pbl.nl/hyde/hyde3.2/ Knoema (2017). World Data Atlas: China. https://knoema.com/atlas/China. Koellner T, Schmitz CJ (2006) Biodiversity, ecosystem function, and investment risk. Bioscience, 56:977–985. Li, X., Wellen, C., Liu, G., Wang, Y., & Wang, Z. L. (2015). Estimation of nutrient sources and transport using Spatially Referenced Regressions on Watershed Attributes: a case study in Songhuajiang River Basin, China. Environmental science and pollution research international, 22(9), 6989. Liang, J. , Zhong, M. , Zeng, G. , Chen, G. , Hua, S. , & Li, X. , et al. (2017). Risk management for optimal land use planning integrating ecosystem services values: a case study in changsha, middle china. Science of The Total Environment, 579, 1675-1682. McDonald, R. I., Weber, K. F., Padowski, J., Boucher, T., & Shemie, D. (2016). Estimating watershed degradation over the last century and its impact on watertreatment costs for the world's large cities. Proceedings of the National Academy of Sciences of the United States of America, 113(32), 9117. 27

Mockler, E. M., Deakin, J., Archbold, M., Gill, L., Daly, D., & Bruen, M. (2017). Sources of nitrogen and phosphorus emissions to irish rivers and coastal waters: estimates from a nutrient load apportionment framework. Science of the Total Environment, 601-602, 326. More, A. L. , Beusen, A. H. W. , Bouwman, A. F. , & Willems, W. J. . (2013). Exploring global nitrogen and phosphorus flows in urban wastes during the twentieth century. Global Biogeochemical Cycles, 27(3), 836–846. Mouri, G., Shinoda, S., & Oki, T. (2010). Estimation of total nitrogen transport and retention during flow in a catchment using a mass balance model incorporating the effects of land cover distribution and human activity information. Water Science & Technology A Journal of the International Association on Water Pollution Research, 62(8), 1837-47. National Research Council (NRC) (1985). Nutrient Requirements of Domestic Animals Series. Committee on Animal Nutrition. National Academy Press, Washington, DC. Organization for Economic Co-Operation and Development (OECD) (2015). OECD Environmental Data. http://stats.oecd.org/ Pathak, D., Whitehead, P. G., Futter, M. N., & Sinha, R. (2018). Water quality assessment and catchment-scale nutrient flux modeling in the ramganga river basin in north india: an application of inca model. Science of the Total Environment, s 631–632, 201-215. 28

Peterson, B. J., Wollheim, W. M., Mulholland, P. J., Webster, J. R., Meyer, J. L., & Tank, J. L., et al. (2001). Control of nitrogen export from watersheds by headwater streams. Science, 292(5514), 86-90. Santos, R. M. B., Fernandes, L. F. S., Cortes, R. M. V., Varandas, S. G. P., Jesus, J. J. B., & Pacheco, F. A. L. (2017). Integrative assessment of river damming impacts on aquatic fauna in a portuguese reservoir. Science of the Total Environment, 601602, 1108. Sharpe, W. F. (1970). Portfolio theory and capital markets. McGraw-Hill College. Smith, V. H. (1998). Cultural Eutrophication of Inland, Estuarine, and Coastal Waters. Successes, Limitations, and Frontiers in Ecosystem Science. Springer New York. Smil, V. (1999). Nitrogen in crop production: an account of global flows. Global Biogeochemical Cycles, 13(2), 647-662. United Nations, Department of Economic and Social Affairs, Population Division (2017). World Urbanization Prospects: The 2017 Revision, Online Edition. Van Drecht, G & Bouwman, Alexander & Knoop, J.M. & Beusen, A & R. Meinardi, C. (2003). Global modeling of the fate of nitrogen from point and nonpoint sources in soils, groundwater, and surface water. Global Biogeochemical Cycles, 17. 261. Vitousek, P. M., Aber, J. D., Howarth, R. W., Likens, G. E., Matson, P. A., & Schindler, D. W., et al. (1997). Technical report: human alteration of the global nitrogen cycle: sources and consequences. Ecological Applications, 7(3), 737-750. 29

Wu, Y., & Chen, J. (2009). Simulation of nitrogen and phosphorus loads in the Dongjiang River basin in South China using SWAT. Frontiers of Earth Science in China, 3(3), 273-278. Yang, X., Liu, Q., Fu, G., He, Y., Luo, X., & Zheng, Z. (2016). Spatiotemporal patterns and source attribution of nitrogen load in a river basin with complex pollution sources. Water Research, 94, 187-199. Yu, C., Huang, X., Chen, H., Godfray, H.C.J., Wright, J.S., Hall, J.W., Gong, P., Ni, S., Qiao, S., Huang, G., Xiao, Y., Zhang, J., Feng, Z., Ju, X., Ciais, P., Stenseth, N.C., Hessen, D.O., Sun, Z., Yu, L., Cai, W., Fu, H., Huang, X., Zhang, C., Liu, H. and Taylor, J. (2019) Managing nitrogen to restore water quality in China. Nature 567(7749), 516-520. Zhang, J., Gong, H., Li, X., & Ross, M. (2009). Assessment of Nitrogen Loading in Miyun Reservoir Beijing Using Bayesian Decision Network. International Conference on Bioinformatics and Biomedical Engineering (Vol.7, pp.1-4). IEEE. Zhang, X., Davidson, E.A., Mauzerall, D.L., Searchinger, T.D., Dumas, P. and Shen, Y. (2015) Managing nitrogen for sustainable development. Nature 528, 51 Zhou, Q., Gibson, C. E., & Foy, R. H. (2000). Long-term changes of nitrogen and phosphorus loadings to a large lake in north-west Ireland. Water Research, 34(3), 922-926. Zhu, Z. L. , & Chen, D. L. . (2002). Nitrogen fertilizer use in china – contributions to food production, impacts on the environment and best management strategies. 30

Nutrient Cycling in Agroecosystems, 63(2-3), 117-127.

Fig. 1. 211 watersheds in China for constructing the relationship between anthropogenic activities and N loading.

31

Fig. 2 Tendencies over 60 years in China. (a) time series of agricultural area and crop area. (b) population changes with time

32

Fig. 3 Variation of fertilizer N application rate and manure N application rate from 1949 to 2050 in China (these rates represent the N loading per hectare in China, they were calculated by the function of rate= total fertilizer N (manure N) of China/total inland of China). The rate from 1949 to 2016 were calculated by China Agriculture Yearbook and the future rate by the linear trend of history.

33

Fig. 4 Future population projection change in China: (1) “No change” variant, both fertility and mortality are kept constant; (2) “Low variant”, low-fertility assumption; (3) “Medium variant”, medium-fertility assumption (4) “High variant”, high-fertility assumption

Table 1 The elements and sources of MBM to estimate the spatially distributed N loading in 2010 and 2005. All fields were resampled to 5 minutes (longitude × latitude) spatial scale. Element (N loading)

Source/achieve

N deposition

Dentener and Crutzen (1994)

Undisturbed vegetation

land cover type data (http://www.landcover.org/data/lc/, data for 2005 and 2010)

Fertilizer consumption

Knoema databases (https://knoema.com/atlas/China)

Animal intake and excretion

China Agriculture Yearbook (CAY, 2011)

Crops production data

FAOATAT databases (FAO 2012)

Distribution of population

HYDE version 3.2 (Klein Goldewijk et al., 20xx (in prep))

Distribution of cropland

HYDE version 3.2 (Klein Goldewijk et al., 20xx (in prep))

Distribution of agriculture

HYDE version 3.2 (Klein Goldewijk et al., 20xx (in prep))

Runoff

Global Runoff Data Center river discharge database (http://www.grdc.sr.unh.edu)

Precipitation

WorldClim-Global Data (http://worldclim.org/)

Sewage treatment

OECD (2015) https://stats.oecd.org/ 34

groundwater storage, soil moisture

available at ftp://podaac-ftp.jpl.nasa.- gov/

Table 2 Regression coefficients describing the relationship between predictor variables and N load Variable*

Estimate

T value

P

a (Crop area × fertilizer N application rate)

0.041

3.223

0.001

b (Agriculture area × manure N application rate)

0.050

3.772

0.000

c (Population)

0.338

12.395

0.000

d (Watershed area)

0.535

20.173

0.000

e (Intercept)

3.211

19.360

0.000

*All variables are standardized by log10-transformed

Table 3. Basic parameters for optimal portfolio analyses No change

Low variant

Medium variant

High variant

Probabilities of future scenarios Historic Trend

0.40

0.30

0.20

0.10

Uniform

0.25

0.25

0.25

0.25

Average watershed nitrogen loading index (WNLI) a

a

First stage b

0.678

0.724

0.612

0.503

Second stage

0.656

0.722

0.471

0.226

Third stage

0.825

0.876

0.463

0.058

WNLI values are derived from N loading predictions of the future in different population increase

scenarios; variances and covariances are calculated by the authors. b

First stage: year 2021-2030; Second stage: year 2031-2040; Third stage: year 2041-2050.

35

Fig. 5 The distribution of watershed N loading in China for 2005 and 2010 by mass balance model. It was classified 5 grades by equal intervals.

36

Fig. 6 For all the watersheds in China, the correlation between the estimates of N loading in this study were predicted by the Mass Balance Model and Statistical Model for 2005 and 2010 (Tg =109 kg). The red line is the 1:1 line. The number of point (N) is 211.

Fig. 7 The estimate of N loading in China over more than 60 years (1949-2010). Including the variation of N loading with time and the inter-annual variation of N loading. The dash represents 37

the base line (i.e., zero inter-annual variation of N loading).

Fig. 8 Future N loading of China under four population projection variants (no change, low variant, medium variant and high variant) predicted by statistical models

38

Fig. 9 Results of portfolio optimization analysis for two probability distributions. Points A to E are the expected values of WNLI for the “uniform” probability. Points F to J are the expected value of WNLI for the “Historic Trend” probability.

Optimization of management strategies for reducing nitrogen loading in China Shanshan Hua a, b, Xin He c, Chunmiao Zheng b, a School



of Environment and Energy, Peking University, Shenzhen 518055, China

Corresponding author, Chunmiao Zheng, School of Environmental Science and Engineering, Southern University of Science and Technology, Shenzhen 518055, China. Tel.: +86 0755-88010020; E-mail address: [email protected] 39

b

State Environmental Protection Key Laboratory of Integrated Surface Water-

Groundwater Pollution Control, School of Environmental Science and Engineering, Southern University of Science and Technology, Shenzhen 518055, China c Department

of Water Resources, China Institute of Water Resources and Hydropower

Research, Beijing 100038, China

Table 4 Selected results of optimal portfolio analyses Portfolio weights a

Poin t on

First

Secon

Outcomes Third

Uncertaint E[R]

figure

stage 0.285

d stage 0.714

stage

y

0.000

A

0.550 0.1574

5

5

0

4

0.464

0.535

0.000

0.570

B

0.1633 1

9

0

1

0.642

0.357

0.000

0.589

C

0.1798 7

3

0

7

0.821

0.178

0.001

0.609

D

0.2044 0

0

0

4

1.000

0.000

0.000

0.629

E

0.2345 0

0

0 40

0

0.285

0.714

0.000

5

5

0

0.483

0.450

0.065

9

8

3

0.538

0.232

0.229

1

5

4

0.592

0.014

0.393

3

3

5

0.000

0.000

1.000

0

0

0

0.1574

0.614

F 4 0.1707

0.633

G 6 0.1896

0.652

H 8 0.2086

0.672

I 0 0.2938

0.691

J

a

2

The portfolio weights of investment in different managed stages under the

different expects and risks.

41

1. A statistical model is constructed to estimate the long-term N loading for China. 2. Portfolio optimization analysis is used to manage the N loading in next 30 years. 3. The past and future N loading in China at national scale is estimated. 4. N loading to landmass in China increased from 22.2 Tg in 1950 to 80.8 Tg in 2010.

42