Predictions of soil and nutrient losses using a modified SWAT model in a large hilly-gully watershed of the Chinese Loess Plateau

Wenhai Shi ,Mingbin Huang

a The State Key Laboratory of Soil Erosion and Dryland Farming on the Loess Plateau,Institute of Soil and Water Conservation,Northwest A&F University,Yangling,Shaanxi Province,712100,China

b School of Water and Environment,Chang’an University,Xian,710054,China

c CAS Center for Excellence in Quaternary Science and Global Change,Xian,Shaanxi,710061,China

d University of Chinese Academy of Sciences,Beijing,100049,China

Keywords:Sediment Runoff Total N Total P

A B S T R A C T Hydrological models are effective tools for assessing the effects of soil and nutrient losses on land degradation.SWAT(Soil and Water Assessment Tool)model is widely used to simulate soil and nutrient losses caused by various management regimes.However,its performance of predicting nutrient loss has not been assessed adequately on the Loess Plateau.This study proposed a modified SWAT model by incorporating the modified Soil Conservation Service curve number method,the storm-based Chinese soil loss equation and the nutrient loss model.The observed daily data of runoff and sediment over 16 years and the monthly soluble phosphorus(P)and nitrate losses over 9 years and 4 years,respectively at the outlet of the upper Beiluo river(UBR)basin were used to assess the model performances.Global sensitivity and uncertainty analyses of parameters to runoff,sediment,soluble P and nitrate in the modified SWAT were conducted.The findings during calibration and validation showed that the modified SWAT was highly accurate in terms of model efficiency(calibration:0.83,0.83,0.48,and 0.49;validation:0.58,0.57,0.53,and 0.65)for runoff,sediment,soluble P loss and nitrate loss,respectively.High model efficiency indicated that the modified SWAT could accurately predict soil and nutrient losses at the river basin scale for the Loess Plateau.Moreover,the temporal variations from month to year and the spatial variations at the sub-basin scale for soil loss and the total N and P losses were analysed using the data simulated by the modified SWAT.The results indicated that the critical loss period occurred in July and August,and the Grain for Green project significantly affected the hydrological behaviour and reduced the soil and nutrient losses in the UBR basin.

1.Introduction

The Loess Plateau region,which has semi-humid,semi-arid,and arid climate with heavy rainfall intensity,highly dissected topography,sparse vegetation coverage and loose loessial soil,is one of areas in China and the world that soil erosion is most severe,and almost 90% of the sediments generated by this erosion flows into the Yellow River(Fu et al.,2004;Huang et al.,2003;Huang &Zhang,2004).In addition to runoff and eroded surface soil,a large amount of soil nitrogen and phosphorus,especially nitrate and soluble phosphorus(P),are transported into the Yellow River and finally into the marine ecosystem(Zhao et al.,2013),which leads to soil degradation,non-point source pollution(NPS)and water body eutrophication(Fu et al.,2011).

Due to the temporal and spatial variation in nutrient loss,it is hard to identify and assess nutrient loss in response to hydrological processes at the watershed scale(Zhang&Huang,2011).Modelling provides a useful analysis tool to assess nutrient loss and develop best management practices,which ultimately reduces soil degradation and NPS loading into water bodies(Shen et al.,2012).The existing nutrient loss models can generally be divided into two categories:empirical and process-based models(Bi&Wang,2012).The empirical models,such as the export coefficient method(Johnes,1996),have advantages such as ease of application and simple parameterization;however,these models often lack a theoretical basis and have low accuracy(Zhang et al.,2007).In contrast,the process-based models are efficient tools to analyse the runoff,erosion processes and impacts of nutrient loss on water quality(Liao &Tim,1997).Because soil and water losses result in soil nutrient loss,runoff and sediment are the most prominent processes that also require precise studies to solve water quality problems.Therefore,a number of basin-scale models have been developed with the water,sediment,and nutrients components,but these models differ in complexity and data entry requirements(Borah&Bera,2004).The annualized agricultural non-point source(AnnAGNPS)model(Baginska et al.,2003),the Hydrological Simulation Program-Fortran(HSPF)(Bicknell et al.,1996),the European Hydrological System model(MIKE SHE)(Refsgaard&Storm,1995),and the Soil and Water Assessment Tool(SWAT)(Abbaspour et al.,2015;Arnold et al,1998,2012)are among these commonly used models.In the Loess Plateau region,nutrient loss has been simulated using these existing models.For example,Wu et al.(2016)used the export coefficient model to predict the total nitrogen(TN)and total phosphorus(TP)load in the Beiluo river basin,while Jia(2005)and Li et al.(2008)used the AnnAGNPS model to simulate the respective runoff,sediment and nutrient loss from the Zhuanyaogou and Heihe river basin.However,most of these studies have only predicted the nutrient loss at an annual scale.

Some researchers have evaluated and reviewed different watershed models(Borah&Bera,2004;Gassman et al.,2007)and have concluded that the SWAT is a widely used model for continuous runoff,soil loss,and nutrient loss simulations(Behera &Panda,2006;Jeong et al.,2010;Pile et al.,2016;Clement &Djebou,2018;Ayivi & Jha,2018;Briak et al.,2019;Guse et al.,2019;Pandey et al.,2019;Malago’et al.,2019;Tan et al.,2019;Abbasi et al.,2019).The SWAT model has been tested and has estimated nutrient loss at sites from around the world(Abbaspour et al.,2007;Douglas-Mankin et al.,2010;Ouyang et al.,2009)as well as in China(Yang et al.,2011).However,most of the previous studies that have applied the SWAT model in the Loess Plateau region mainly focused on the hydrological processes(Li et al.2009,2017;Yan et al.,2017;Zuo et al.,2016);few previous studies have comprehensively evaluated the accuracy of SWAT in terms of predicting nutrient loss at the watershed scale(Qin et al.,2009).This may reflect the weakness of SWAT in terms of its application to the Loess Plateau,which are beyond the spatial area where the model was originally developed(Yu et al.,2013);thus,specific improvements are needed(Gassman et al.,2014).

For the direct runoff calculation module,the Soil Conservation Service curve number(SCS-CN)method was used in the SWAT model(SCS,1972).However,many studies from around the world have reported several drawbacks in this method,especially regarding the applicability of the SCS-CN method,which needs further improvement(Jiao et al.,2015;Mishra et al.,2007;Singh et al.,2015).One of the main drawbacks is that the constant initial abstraction(Ia)is used as an intrinsic parameter,which may result in inconsistency of the soil moisture accounting(SMA)(Michel et al.,2005).Another major limitation is that rainfall intensity and duration are ignored in the SCS-CN method(King et al.,1999).To eliminate these drawbacks,Sahu et al.(2012)derived an enhanced SCS-CN model by incorporating storm duration,and Singh et al.(2015)improved the method by embedding a SMAbased function for antecedent soil moisture.More recently,Shi et al.(2017)incorporated the antecedent moisture content prediction and rainfall intensity to estimate runoff in the Loess Plateau,which led to an improved model.

In the sediment estimation module,the Modified Universal Soil Loss Equation(MUSLE),developed by Williams and Berndt(1977),was adopted in the SWAT model.However,Liu et al.(2002)reported that the systematic practices of soil conservation used in China were distinctive to those areas where the original MUSLE model was developed.Subsequently,according to the typical management practices for soil conservation employed in the Loess Plateau,Liu et al.(2002)developed the Chinese soil loss equation(CSLE)to predict gross annual soil erosion in China(Li et al.,2012;Wang et al.,2015;Zhu et al.,2012).However,in order to predict storm-based soil erosion,Shi,Huang,and Barbour(2018)introduced a new factor based on the CSLE,which accounted for both rainfall and runoff characteristics;moreover,this modified method has been validated using a series of measurements collected from plots with different types of soil,vegetation and land use in the Loess Plateau.

In the sediment-bound nutrient transport modules,sedimentbound nutrient loss mainly depends on the original soil nutrient concentration,the soil loss,and the enrichment ratio(ER);furthermore,a logarithmic equation between ER and sediment nutrient concentrations was developed by Menzel(1980)and adopted in SWAT.However,Shi,Huang,and Wu(2018)found that the soil erosion,rainfall characteristics,steepness and land-use have significant influence on ER in loess soils;thus,a sedimentbound nutrient loss model was developed to account for the combined effects of these factors.For the runoff-associated nutrient transport modules,the amount of solutes migrated into runoff is predicted using the runoff volume,the solute concentration in the top soil of 10 mm,and an empirical coefficient(Neitsch et al.,2011).However,Wang et al.(1999)reported that the nutrient transport module of original SWAT was not applicable to the Loess Plateau;thus,the authors developed the equivalent convection model,in which the solutes transport is assumed to be dependent on the mass transfer coefficient,km.Shi,Huang,and Wu(2018)further modified the model by replacing kmwith the sediment and runoff rates to predict the losses of soluble P and nitrate.

This study focused on whether the integration of these modules developed in local areas into SWAT could improve the model’s efficiency to predict runoff,soil loss and nutrient loss and identify the critical area and period of nutrient loss at the river basin scale in the Loess Plateau.Thus,the objectives of this study are to:(1)develop a modified version of SWAT by embedding the modified SCS-CN,CSLE and nutrient loss modules into the SWAT model;(2)compare the performance of the modified and original SWAT models in terms of predicting runoff,soil loss and runoff-associated nutrient loss on the Loess Plateau;and(3)identify the spatial and temporal variation of the soil and nutrient losses.

2.Materials and method

2.1.Brief description of original SWAT

SWAT is a semi-distributed,continuous watershed simulator,which was developed to assess the impacts of management practices and climate on flows,soil loss,and nutrient loss in watersheds.SWAT(version 2012)was adopted within ArcGIS 10.3 in this study(Dile et al.,2016).The major components of SWAT include hydrology,climate,erosion,land cover and plant,nutrients and pesticides,management practices and channel process.An interface in the ArcGIS tool to define hydrological characteristics of a watershed and store relevant spatial and tabular data,have been provided by SWAT,which can divide watershed into sub-basins for high-level spatial detail simulation;moreover,based on unique combinations of land cover,soil and management,these sub-basins are then further discretized into hydrological response units(HRUs).

The SWAT model comprises two phases:a land phase solved at HRU level,and a stream phase solved at reach(sub-basin)level.The land phase comprises the computation of HRU daily water balance,sediment yields and nutrient loadings.The HRU daily water balance considers precipitation,irrigation,evapotranspiration,surface runoff,lateral flow,and percolation to shallow and/or deep aquifers.Daily outputs of all HRUs of a sub-basin are routed through the stream network(stream phase).The stream phase comprises the routing runoff,sediments and nutrients in the cascading sequence of reaches composing the stream network.The runoff,soil loss,and nutrient loss of each sub-basin are routed to the reach and enter the outlet of the watershed along the stream network(Neitsch et al.,2011).

2.2.The improvements of SWAT

2.2.1.Modified surface runoff module

In the original SWAT,the SCS-CN method(SCS,1972)was used to simulate surface runoff,which can be expressed as follows:

where Rdayis the daily rainfall depth(mm);Q is the daily runoff depth(mm);Iais an initial abstraction(mm),which is commonly approximated as 0.2 S;and S is the potential maximum retention(mm),which is defined by Eq.(2):

where CN is a non-dimensional value that varies from 0 to 100.

In the original SCS-CN,the incorporation of antecedent moisture conditions(AMC)using three discrete levels caused a discontinuous relationship between CN and soil water content,which often results in unreasonable abrupt spikes in the predicted runoff.Moreover,the constant initial abstraction(Ia)is used as an intrinsic parameter,which may result in inconsistency of the soil moisture accounting(SMA)(Michel et al.,2005).Another major limitation is that rainfall intensity and duration are ignored in the SCS-CN method(King et al.,1999).

To eliminate these limitations,Shi et al.(2017)modified the SCSCN based on the three hypotheses:i)the runoff coefficient is equal to the degree of saturation(Singh et al.,2015);ii)the threshold soil moisture Sais equal to the sum of initial soil moisture V0and initial abstraction Ia(Michel et al.,2005);and iii)the cumulative infiltration(F)is equal to the sum of the static infiltration(Fc,occurring largely due to gravity)and dynamic infiltration(Fd,occurring mainly due to capillarity).Therefore,the modified SCS-CN model of Shi et al.(2017)was used for the runoff estimation in the modified SWAT,which is expressed as:

where Sais the threshold soil moisture(mm);Sbis the sum of S and Sa(mm);Fcdenote the static infiltration which can be calculated by the product of rainfall duration and the minimum infiltration rate(fc,mm h-1);β is an empirical parameter;and V0is soil moisture store level of the 0-15 cm soil layer(mm)at the beginning of simulation time-step of a continuous model,which was obtained by the calculation of the original and modified SWAT models.

2.2.2.Modified sediment estimation module

In the original SWAT,sediment yield is estimated with the MUSLE method which can be depicted in the following general form(Williams,1995):

where sed is the soil loss(metric tons);Q is the volume of surface runoff(mm);qpeakis the peak flow rate(m3s-1);areahruis the area of an HRU(ha);and KUSLE,CUSLE,PUSLEand LSUSLEare the soil erodibility(t h MJ-1mm-1),crop management,conservation practice and topography factors,respectively.Finally,CFRG is the coarse fragment factor.In the MUSLE method,the systematic practices of soil conservation were distinctive to the areas of China(Liu et al.,2002).Subsequently,according to the typical management practices for soil conservation employed in the Loess Plateau,Liu et al.(2002)developed the Chinese soil loss equation(CSLE)to predict gross annual soil erosion in China.For a storm event,the soil loss is determined from the product of the runoff amount and the bulk sediment concentration.Kinnell(1997)reported that sediment concentration for a storm event is dependent on the kinetic energy per unit quantity of rain Ek(given by ERdivided by rainfall amount P,MJ ha-1mm-1)and I30.Because both the direct runoff Q and EkI30exhibit a good correlation with storm-based soil loss,Shi,Huang,and Barbour(2018)introduced the term QEkI30as a rainfall-runoff erosivity factor based on the CSLE,which accounted for both rainfall and runoff characteristics to predict storm-based soil erosion.The modified CSLE model of Shi,Huang,and Barbour(2018)was adopted in the modified SWAT,which can be expressed as follows:

where Ekis the kinetic energy per unit quantity of rain(MJ ha-1mm-1);I30is the maximum 30-min rainfall intensity(mm h-1);K,L,Se,B,E,and Teare the soil erodibility factor(t h MJ-1mm-1),slope length,steepness,biological control(cover),supporting practice and tillage factor(dimensionless),respectively;and a and b are the empirical parameters.

2.2.3.Modified nutrient transport modules

In the original SWAT,the organic N and P and the nitrate and soluble P that are transported in the surface runoff are calculated as(Neitsch et al.,2011;Williams & Hann,1978):

where ERiis defined as the enrichment ratio of organic nutrient i(N or P)(dimensionless);ORGiis the loss of organic nutrient i(kg ha-1);Coniis the concentration of original soil nutrient i(g Mg-1);concsed,surqis the concentration of sediment in the runoff;Psurfand NO3,surfare the soluble phosphorus and nitrate migrated into the runoff(kg ha-1),respectively;Psolutionand concNO3,mobileare the concentration of phosphorus and nitrate,respectively,in the mobile water of the top soil of 10 mm(kg mm-1);βNO3and kdare the nitrate percolation and phosphorus soil partitioning coefficient(m3mg-1),respectively;ρbis the dry soil bulk density(mg m-3);and H0is the depth of the surface soil layer(10 mm).

In the original sediment-bound nutrient transport modules,a logarithmic equation between ER and sediment nutrient concentrations was adopted in SWAT.However,Shi,Huang,and Wu(2018)found that the soil erosion,rainfall characteristics,steepness and land-use have significant influence on ER in loess soils;thus,a sediment-bound nutrient loss model was developed to account for the combined effects of these factors.

For the original runoff-associated nutrient transport modules,the number of solutes migrated into runoff is predicted using the runoff volume,the solute concentration in the top soil of 10 mm,and an empirical coefficient(Neitsch et al.,2011).However,Wang et al.(1999)reported that the nutrient transport module of original SWAT was not applicable to the Loess Plateau;thus,the authors developed the equivalent convection model,in which the solutes transport is assumed to be dependent on the mass transfer coefficient,km.Shi,Huang,and Wu(2018)further modified the model by replacing kmwith the sediment and runoff rates to predict the losses of soluble P and nitrate.

Fig.1.Locations(a),DEM(b),soil types(c)and land use(d)of the upper Beiluo river(UBR)basin.

In the modified SWAT,the nutrient loss method proposed by Shi,Huang,and Wu(2018)was adopted and modified as follows,in which Eq.(8)was also used for the organic N and P estimations,as in the original SWAT,and the runoff duration was replaced by the time of concentration for runoff tconc:where Siand Ciare the slope and land-use factor for organic nutrient i,respectively;cs0,pand cs0,ndenote the initial concentration of solute P and nitrate within H0(kg ha-1),respectively;θsis the saturated soil moisture(cm3cm-3);ρbis the dry soil bulk density(g cm-3);k is the soil adsorption rate(cm3g-1);tconcis the time for concentration(h);t is the runoff duration(h);aiand biare the empirical parameters for organic nutrient i;λP,λN,bPand βP,βN,bNare the empirical coefficients for solute P and nitrate;and 100 is the unit conversion factor.

2.3.Application

2.3.1.Study area

This study area is located in the upper reaches of the Beiluo river(UBR)basin.The Beiluo river basin is a first-grade tributary of the Weihe river basin,which is a secondary tributary of the Yellow River.The studied UBR basin(latitude:36°44′N-37°19′N;longitude:107°33′E-108°32′E;elevation:1270-1889 masl),which is monitored by the Wuqi hydrological station with an area of 3436 km2,and accounts for 12.7% of the entire Beiluo river basin(Fig.1).Wuqi and Dingbian counties located in the basin,where the former covers 77%,while the latter covers 23%of the total area(Ran et al.,2006).

The UBR basin is a typical hilly-gully area with highly dissected topography and gully densities of 2-3 km km-2.The climate of this region is semi-arid with a mean annual temperature and precipitation of 7.5°C and 448 mm,respectively,and 78% of this precipitation falls during the flood season i.e.,between May and October(Chen et al.,2016).Loess soil is the dominate soil type,which comprises 84.6% of the study basin.This loessial soil is the aeolian deposit that is loosely cemented by calcium carbonate with less than 20% clay,more than 50% silt,and 5% sand.

2.3.2.Data collection

Daily meteorological data such as solar radiation,temperature,relative humidity and wind speed were collected from Wuqi,Dingbian and Jingbian climate stations.The three stations are in and around the study basin,and the data of which from 1999 to 2016 were provided by the China Meteorological Administration(CMA).Daily precipitation records during the corresponding periods at 14 rain-gauge stations were used in this study(Fig.1).

Digital elevation model(DEM)data,soil map(1:1000,000)and land-use map in 2006 in the river basin with the spatial resolution of 30 m were collected from the National Earth System Science Data Sharing Infrastructure(http://loess.geodata.cn/)(Fig.1).According to the soil types and land-use derived from the soil and land-use maps,23 sub-basins were classified,and 861 HRUs were discretized in this river basin.The soil property data including soil texture,available water content,bulk density,hydraulic conductivity and organic carbon content for each soil type were provided by the Chinese Soil Database’s website:http://vdb3.soil.csdb.cn/,and the nitrogen and phosphorus content data were obtained from Zhao(2017),who investigated the spatial distribution of soil nutrients in this basin by field sampling.The land-use was divided into five main categories:farmland,grassland,woodland,urban and water.The major crops in the watershed are millet,soybean and potato.Millet,soybean and potato were planted in mid-May and followed by a fertilizer application of 120 kg N ha-1and 60 kg P ha-1before harvest occurred from late September to early October.

2.3.3.Data analysis

In this study,the performances of the models were evaluated using the indices of the percent bias(PBIAS),Nash-Sutcliffe efficiency(NSE)(Nash& Sutcliffe,1970),P-factor and R-factor.

The P-factor is the percentage of observation bracketed by the 95% prediction uncertainty(95PPU),while the R-factor is the average thickness of the 95PPU band divided by the standard deviation of the observed data.A satisfactory model performance and efficiency was that the P-factor approaching 100% while the Rfactor approaches one(Singh et al.,2013).

The indices of NSE and PBIAS can be calculated as follows:

where Ojand Pjare the observed and simulated values,respectively;Omis the average of all observations;and n is the total number of paired values.PBIAS is used to determine how well the model simulates the average magnitudes for the output response of interest,with small magnitude of PBIAS values are preferred.Moriasi et al.(2015)suggested that the satisfactory performance of the model was that NSE values for flow,sediment,nutrient loss should be greater than 0.50,0.45,0.35,respectively.

2.3.4.Parameter estimations

Previous studies have shown that the SUFI-2 program(Sequential Uncertainty Fitting,ver.2)(Abbaspour,2007)is very efficient for calibration(Abbaspour et al,2007,2015;Arnold et al.,2012).Therefore,SUFI-2 was used for the parameter estimations of the original and modified SWAT models in this study.The simulation period for the calibration of the daily runoff and soil loss was from 2001 to 2010;moreover,the validation period was 2011-2016,and the first 2 years(1999-2000)were used as a warm-up period.The monthly soluble P loss data from 2006 to 2010 and 2013-2016 were used for model calibration and validation,respectively.For the monthly nitrate estimation,the simulation period for calibration was 2013-2014,and the validation period was 2015-2016.The calibration of model was implemented in the order of 1)stream flow,2)sediment,and 3)nutrients,because model outputs of nutrient yields were dependent on sediment yields,and sediment yields are driven by the water balance.The sediment(or nutrient loss)calibration was implemented by including only parameters directly related to sediment(or nutrient)while flow(or flow and sediment)parameter values were fixed to their calibrated value.

The calibration procedures were conducted using the Nash-Sutcliffe Efficiency as the objective function and running model 1000 times for each iteration.After each iteration,the range of each parameter was updated(normally narrowed down)based on both the new parameters suggested by the SWAT-CUP tool and their reasonable physical ranges.This procedure was repeated three to five times until no further improvement of model performance could be achieved.

3.Results

3.1.Parameter sensitivity analysis

Table 1 lists both the original and modified SWAT models’sensitive parameters and the final values that were optimized in thecalibration process.Twelve,seven,five and four model parameters were used for runoff,sediment,soluble P and nitrate calibration and validation of original SWAT model,respectively.An additional two parameters(β,fc)(see Eq.(4)),two parameters(a and b)(see Eq.(6)),four parameters(λP,βP,bP,k)(see Eq.(13))and three parameters(λN,βN,bN)(see Eq.(14))were also used for runoff,sediment,soluble P and nitrate calibration and validation of the modified SWAT model,respectively.Moreover,all the above parameters were applied for sensitivity analysis of the modified SWAT.

Table 1 Parameters of the two tested models that were fitted and their final calibrated values.

Fig.2 shows the parameter sensitivity analysis of the modified SWAT model based on p-value and t-stat.The p-value determines the significance of the sensitivity,and the value closer to zero is more significant,while t-stat provides the extent of sensitivity.

As presented in Fig.2,the SCS runoff curve number(CN2),sediment routing factor(ADJ_PKR),channel re-entrained linear parameter(SPCON),phosphorus uptake distribution parameter(P_UPDIS),nitrate percolation coefficient(NPERCO),nitrogen uptake distribution parameter(N_UPDIS),and the three parameters of λN,βN,and bNin Eq.(14)are the eight most sensitive parameters having very high t-stat with p-values close to zero(p<0.05).The sensitive of parameter CN2indicated that the characteristics of the watershed,which represented by CN2are significantly affecting runoff,while the sensitive of channel parameters(SPCON,ADJ_PKR)showed that sediment transportation and deposition is also one of the most important factors that impact on sediment yield.The parameters of nitrate movement(NPERCO,λN,βN,bN)have a larger influence on nitrate loss according to the sensitivity analysis.Moreover,the plant nutrient uptake distribution parameter(P_UPDIS and N_UPDIS)also affect nutrient losses of soluble P and NO3.

3.2.Calibration and validation

3.2.1.Runoff

Fig.2.Global sensitivity analyses of parameters to(a)runoff,(b)sediment,(c)soluble P and(d)nitrate in the modified SWAT.

Table 2 presents the comparative overall performance of runoff,sediment,soluble P and nitrate losses following the evaluation criteria obtained from the original and modified SWAT models.Fig.3a and b shows the time series of observed versus simulated daily runoff for model calibration(2001-2010)and validation(2011-2016)in the UBR basin using the original and modified SWAT models.It can be seen from Fig.3a that the simulated runoff did not match the trend of the measured runoff very well,as the magnitudes of the simulated runoff were lower than those of their measured counterparts during the peak events.However,the graphical representation of Fig.3b indicates that the temporal variation in the runoff simulated by the modified SWAT model was more closely matched with the observed values than those simulated by the original SWAT model.Moreover,the modified SWAT shows less uncertainty than the original SWAT,with 75%and 85%of the measured data bracketed by the 95PPU and the R-factor values of 0.90 and 1.59 during the discharge calibration and validation,respectively.Table 2 also shows that the performance of the modified SWAT model yielded a higher NSE value of 0.83 and a lower PBIAS value of 4.1% than those of the original SWAT model(NSE=0.69,PBIAS=9.4%)for calibration;additionally,the modified SWAT model continued to outperform the original SWAT model during validation period,with the NSE values of 0.58 vs.0.11 and PBIAS values of 12.2% vs.27.1%,respectively.

3.2.2.Sediment

The observed and predicted sediment values from the two tested models are presented and graphically compared in Fig.3c and d.Fig.3c shows that the sediment values are under-predicted by the original SWAT model during the calibration and validation,and most peak events were underestimated.Nevertheless,the predictions of sediment during the whole periods obtained from the modified SWAT model were in good agreement with its measured counterparts(Fig.3c).High NSE values(0.83 and 0.60)further indicate that the modified SWAT model simulated sediment well within the acceptable level of accuracy.Similar to the discharge,the sediment simulation of the modified SWAT also shows fewer uncertainties with higher P-factor values(0.58 and 0.60)and closer R-factor values to one(0.88 and 1.38)as compared to the original SWAT(P-factor=0.15 and 0.08,R-factor=0.66 and 1.94,respectively)in the calibration and validation periods(Table 2).Thus,the modified SWAT model adequately captured the characteristics of the watershed interactions with rains.

3.2.3.Soluble P

The comparisons between the measured and predicted soluble P loss for the original and modified SWAT models are presented in Fig.4a and b.Fig.4a shows that the soluble P was under-predicted by the original SWAT model in most years except 2013.The modified SWAT model exhibited a much better match between the estimated and observed soluble P values than the original SWAT model,with high NSE values(0.48 and 0.53,respectively)for the calibration and validation periods,respectively(Table 2).Moreover,the uncertainty results of modified SWAT were also better than the original model with more measured data bracketed(the P-factor values were 0.48 vs.0.23 for calibration and 0.40 vs.0.24 for validation,respectively),and most of the data missing the 95PPU band were from the small soluble P points,while all the peaks were bracketed.

3.3.Nitrate

The time series of the measured and simulated nitrate values for the two tested models are compared graphically in Fig.4c and d.The simulated results in Fig.4c indicate that the nitrate was wellpredicted by the original SWAT model in the years 2013 and 2016,but there were great differences between the estimated and observed values in 2014 and 2015.However,the statistical test results revealed better agreement between the measured and the estimated nitrate for the modified SWAT model than in that of the original SWAT model;moreover,the NSE values were 0.49 vs.0.41 for calibration and 0.65 vs.-0.22 for validation,respectively(Table 2).Similar to the other variables,the nitrate simulation of the modified SWAT is also better than original SWAT with small uncertainties,while bracketing 58%of the data for calibration and 67%for validation(Fig.4c).

Table 2 Performance of the two tested models for the calibration and validation datasets of the upper Beiluo river(UBR)basin.

Fig.3.Time series of observed versus simulated daily discharge using(a)original SWAT and(b)modified SWAT,and observed versus simulated daily sediment using(c)original SWAT and(d)modified SWAT for model calibration(2001-2010)and validation(2011-2016)in the upper Beiluo river(UBR)basin.

3.4.Temporal and spatial patterns of sediment and nutrient loss

3.4.1.Temporal variations of sediment and nutrient loss

The temporal distribution of the average monthly sediment and nutrient loss of TN and TP simulated by the modified SWAT are presented in Fig.5a.From Fig.5a,the critical period of sediment and TN and TP nutrient losses occurred in July and August.The sediment in July and August accounted for 76.7% of the average annual sediment.The TN and TP nutrient losses in July and August were more than 76%and 87%of the total nutrient loss,respectively.This result indicated that more attention should be focused on July and August in terms of efforts to reduce sediment and nutrient loss.

Fig.5b shows the average annual sediment and TN and TP nutrient losses as well as the component proportion of the average annual TN and TP losses.The annual sediment from the UBR basin showed great variability,ranging from 0.18 to 34.55 million tons,with an average value of 7.92 million tons;moreover,this value sharply decreased by 77.7% compared with the average annual sediment of 1980-2000(i.e.,35.52 million tons)(Zhang et al.,2018).The sediment in 2001 was the highest of all annual soil losses,which might be due to a major flooding period with the maximum daily stream flow of 425 m3s-1.The annual yield of TN and TP also had a similar variability;specifically,TN ranged from 244.50 to 6593.26 t with an average value of 1789.75 t,and TP varied from 12.58 to 2766.20 t with an annual mean of 678.31 t.Therefore,the maximum annual TN and TP loads were 27 and 220 times the minimum ones,respectively.In addition,the nitrate loss,which accounted for 23%of the TN loss,was the major source of the inorganic N loss;furthermore,soluble P accounted for only 1% of the TP.However,the organic N and P accounted for 71%and 99%of the average annual TN and TP,respectively.These results indicated that the organic N and P in the UBR basin were the dominant mechanism of total nutrient loss.This result was supported by Ernstberger et al.(2004),who claimed that,although the regulation measures of inorganic nutrients are concerned,the TN and TP flowing into streams were dominated by the organic forms.

3.4.2.Spatial variations in soil loss and nutrient loss

By using the simulated results of the modified SWAT model,the critical areas of significant sediment and nutrient loss can be identified according to the average annual sediment and nutrient loss for the whole period within each sub-basin.Fig.6 shows the spatial distributions of surface runoff,sediment and nutrient losses at the sub-basin scales of the UBR basin.The rate of sediment ranged from 5.3 to 19.4 t ha-1,and the nutrient TN and TP loss rates varied in the range of 2.9-16.6 kg ha-1and 1.2-5.6 kg ha-1,respectively(Fig.6).Among the 23 sub-basins within the watershed,the numbers 1,2,4,7,8,9,12 and 16 were identified as severely eroded areas(i.e.,up to 11 t ha-1).

Fig.4.Time series of observed versus simulated monthly soluble P using(a)original SWAT and(b)modified SWAT,and observed versus simulated monthly nitrate using(c)original SWAT and(d)modified SWAT for model calibration and validation in the upper Beiluo river(UBR)basin.

From the spatial variations in sediment and nutrient loss,the most severe soil erosion and nutrient loss were distributed in the upper sub-basins of the UBR basin.The main reasons are:(1)These sub-basins are located in the upper reaches of the basin with steep slopes,while downstream areas are mostly flat,which results in less soil erosion and nutrient loss;and(2)these sub-basins that have farmland as the main land-use type are located in the northwestern part of the UBR basin which belongs the area of Dingbian County,while the Grain for Green project was launched in Wuqi County in 1999,which led to significant decreases in farmland(22.8%)and increases in woodland and grassland(22.6%)relative to the values in 1990(Miao et al.,2012;Yan et al.,2017).These results indicated the Grain for Green project significantly affected the hydrological behaviour and reduced the soil erosion and nutrient loss in the UBR basin.

4.Discussion

Although the original SWAT model exhibited good performance in simulating runoff and sediment of calibration,it still lacked accuracy to simulate runoff and sediment of validation.In addition,the original SWAT model could not estimate the soluble P and nitrate losses during both the overall periods within an acceptable level of accuracy.However,the results showed that the modified SWAT model was highly accurate in terms of its NSE values(calibration:0.83,0.83,0.48,and 0.49;validation:0.58,0.57,0.53,and 0.65)for runoff,sediment,nitrate loss and soluble P loss of calibration and validation.The good performance of the modified SWAT model indicates that it could accurately estimated runoff,sediment,soluble P loss and nitrate loss using the optimized(i.e.,calibration)parameters in the UBR basin.

4.1.Comparison of runoff simulation between the original and modified SWAT

In the CN methodology of the original SWAT,S and Iaare updated at the beginning of each simulation time-step in response to rainfall amount,simulated surface runoff,evapotranspiration and soil moisture from the entire soil profile at the previous timestep.Such updating of S and Iaenables the use of CN method in a continuous SWAT model,without violating the integrity of its basic assumptions.Despite of the regular updating at each time-step,the lumped insertion of S and Iain calculating runoff(Eq.(1)and(2))does not allow the method to account for the variation of soil moisture during an event and over a continuous period.Moreover,it is more appropriate to update the S values based on soil water content in the topmost soil layer,which would more closely reflect surface saturation condition during heavy rainfall periods(Qiu et al.,2012).This lacking can be addressed by incorporating the direct contribution of antecedent soil moisture amount of top 15 cm soil layer into Eq.(1).The modified SWAT is incorporated by replacing Eq.(1)with Eqs.(3)and(4)in the sub-routine for runoff volume estimation.The sub-routine that calculates retention parameter at individual time-steps and updates CN based on the antecedent moisture conditions is kept unaltered.This can remove the inconsistency associated with the updating process of S based on daily root-zone soil moisture in a continuous simulation of SWAT model(Gassman et al.,2007).The SMA-based SWAT configuration consistently produces improved goodness-of-fit scores and less uncertain outputs with respect to streamflow during both calibration and validation.

Fig.5.The temporal distribution of(a)multi-year average monthly and(b)yearly sediment and nutrient loss total N and P,and the component proportion of annual average total N and P loss(pie chart)in the upper Beiluo river(UBR)basin.

Another reason for the discrepancy is that the simulated runoff of the original SWAT is not sensitive to rainfall duration.Given the same amount of rainfall,the original SWAT produces the same amount of runoff regardless of event duration or the distribution of rainfall intensity during the rainfall event.Many studies show that storm duration plays is a key role in the watershed’s hydrological processes of both runoff generation and prediction(Fang et al.,2008;Reaney et al.,2010).Shi and Wang(2020)has tested the effect of storm duration on runoff prediction of SCS-CN method by dividing storm events of the Loess Plateau into three regimes based on rainfall depth and duration and found that the effect of storm duration decreased in the order:rainfall regime 3(low intensity and long duration)＞rainfall regime 2(middle intensity and middle duration)＞rainfall regime 1(heavy intensity and short duration).The original SCS-CN method often under-predicted runoff for most of storm events in rainfall regime 1 and over-predicted runoff for most storm events in rainfall regime 3.The modified CN method incorporating storm duration improved runoff predictions for all rainfall regimes.In this study,the calibrated models are also compared over the period 2001-2016 at three different flow regimes using flow duration curve(FDC)analysis(Fig.7).These regimes consist of high,mid and low flows with 0.0-0.1,0.1-0.6 and 0.6-1.0 exceedance probability,respectively.The results showed the modified SWAT had better agreement with the observed data than the original SWAT model for all flow regimes.

4.2.Comparison of sediment and nutrient loss predictions between the original and modified SWAT

Sediment and nutrient loss presented great variability during the past 16 years,and the critical period appeared to occur in July and August(Fig.5).The results were consistent with the fact that heavy storm events mainly occurred during the rainy season and fluctuated from year to year in the UBR basin,and the periods of large sediment and nutrient loss were found to agree with those of peak flow(Zhang et al.,2018).Thus,accurate model simulations of sediment and nutrient loss are largely dependent on accurate simulations of stream flow and peak flows.The original SWAT model has been tested to underestimate high flows in this study,which led to underestimations of sediment and nutrient loss.Improving the accuracy of peak flow estimations is an ongoing challenge for SWAT developers(Arnold et al.,2012;Gassman et al.,2007;Moriasi et al.,2015).However,the modified SWAT exhibits a better match with the observed data for all flow regimes when embedding the modified SCS-CN method in the SWAT model,which could improve model simulations of sediment and nutrient loss during high flow events.

The quality and reliability of the sediment simulations were not only directly linked to the runoff simulations but were also very sensitive to the rainfall characteristic.For rainfall erosion,it consists of the processes of soil materials detachment and transport.It is found that raindrop impact plays an important role in both processes,which not only contributes to soil detachment but also causes runoff disturbance,thereby largely contributing to sediment transport(Kinnell,2010;Zheng,1998).However,the rainfall characteristics,such as rainfall energy and maximum 30-min intensity I30,are ignored in the MUSLE method(Eq.(3))in the original SWAT model.Because of these limitations,some peaks in sediment simulated by the original SWAT model did not match their measured counterparts.In contrast,the product of Ek and I30was used in the modified CSLE model to represent the influence of rainfall impact on soil particle detachment and transport.Moreover,the equation for Ekcalculation was developed using the measured rainfall data of the Loess Plateau,which can accurately characterise the rainfall energy of this region.These factors showed the modified SWAT model can reasonably predict sediment.

Fig.6.The spatial distribution of(a)the sub-basins,the annual average of(b)surface runoff,(c)sediment yield,nutrient loss of(d)NO3,(e)organic N,(f)total N,(g)soluble P,(h)organic P and(i)total P at sub-basin scale of the upper Beiluo river(UBR)basin.

For nutrient loss estimation,the original SWAT model accounted for only a single parameter,i.e.,nitrate percolation coefficient or phosphorus soil partitioning coefficient,to predict nitrate loss or soluble P loss in surface runoff.The model might be applicable to an area with little or no soil erosion,but it is not applicable for the severely erodible region of the Loess Plateau(Wang et al.,1999).To consider the effect of soil erosion on nutrient transport,Gao et al.(2004)introduced a raindrop-induced water transfer rate erin the nutrient loss model.However,studies have shown that more severely erosion leads a greater interaction between the soil particle and runoff;thus,more solutes of the soil are removed to overland flow(Sharpley,2002).Therefore,considering that solute loss has a greater relationship with sediment and runoff,the runoff-associated nutrient loss module in the modified SWAT model uses the sediment and runoff rates instead of the rainfall impact to account for variations in the solute loss(Schiettecatte et al.,2008).These factors ensure the accuracy of modified SWAT model in estimating nutrient losses associated with runoff.

4.3.Limitations and future research

The modified SWAT model was useful for predicting runoff,sediment and nutrient loss at the basin scale in the Loess Plateau;furthermore,the temporal and spatial variation of the sediment and nutrient loss estimated by the modified SWAT model could be used to analyse the dynamic changes that reflect the trends of the sediment and nutrient loss.However,the modified SWAT still has several limitations:

1)CN2from the handbook tables can accurately predict the runoff of traditional agricultural basins in the United States which was obtained and organized from the monitored rainfall-runoff events of 150 watersheds across a wide range of terrain,soil,land use and management conditions,but when applied to other countries and regions,it may lead to the unreliability of actual flood parameters.Lian et al.(2020)adjusted CN2values for 55 study sites in China with available observed rainfallrunoff event data and concluded that only one-third of the study sites could use the SCS handbook table to achieve a satisfactory(±5%)CN2value.Therefore,it is necessary to verify the local application of CN method to reduce the uncertainty of modelling results and promote the more widely use of this method in practice.

2)The modified SWAT model was developed based on the measurement of the Loess Plateau.The model is suitable for regions with similar arid,semi-arid and semi-humid climate if the measured data are available to calibrate the parameters.Since the characteristics and processes of runoff,erosion and nutrient loss in humid regions are different from those in this paper,the modified SWAT model needs to be further verified when applied to humid regions.

5.Conclusions

In this study,the modifications of the SWAT model were proposed by incorporating the modified SCS-CN,CSLE and nutrient loss methods.The reliability of modified SWAT model was tested against the original SWAT model in terms of predicting runoff,sediment,soluble P loss and nitrate loss in the UBR basin of the Loess Plateau.The following conclusions are drawn:

Fig.7.Flow duration curve comparison of streamflow at three different flow regimes(a)0.0-0.1:high flow;(b)0.1-0.6:moist and mid flow;(c)0.6-1.0:dry and low flow of the upper Beiluo river(UBR)basin.

1)The modified SWAT model improved runoff predictions and yielded a higher NSE value of 0.83 and 0.58 than those of the original SWAT model(NSE=0.69 and 0.11)for calibration and validation,respectively.

2)The sediment values are under-predicted by the original SWAT model during the calibration and validation,and most peak events were underestimated.Nevertheless,the predictions of sediment during the whole periods obtained from the modified SWAT model were in good agreement with its measured counterparts.High NSE values(0.83 and 0.60)further indicated that the modified SWAT model simulated sediment well within the acceptable level of accuracy.

3)The soluble P was under-predicted by the original SWAT model in most years except 2013.The modified SWAT model exhibited a much better match between the estimated and observed soluble P values than the original SWAT model,with high NSE values(0.48 and 0.53,respectively)for the calibration and validation periods,respectively.

4)The nitrate was well-predicted by the original SWAT model in the years 2013 and 2016,but there were great differences between the estimated and observed values in 2014 and 2015.However,the statistical test results revealed better agreement between the measured and estimated nitrate values for the modified SWAT model than the original SWAT model;moreover,the NSE values were 0.49 vs.0.41 for calibration and 0.65 vs.-0.22 for validation,respectively.

5)The critical periods of sediment and TN and TP nutrient loss occurred in July and August.This result indicated that more attention should be focused on July and August in terms of efforts to reduce sediment and nutrient loss.The organic N and P accounted for 71% and 99% of the average annual TN and TP,respectively,which indicated that the organic N and P in the UBR basin were the dominant mechanism of total nutrient loss.From the spatial variations in sediment and nutrient losses,the most severe soil erosion and nutrient losses were distributed in the upper sub-basins of the UBR basin,which indicated the Grain for Green project significantly affected the hydrological behaviour and reduced the soil erosion and nutrient loss in the UBR basin.

In general,the original SWAT underestimated runoff,sediment and nutrient loss for some high-flow events.The main reasons could be attributed to its dependence on many semi-empirical equations(such as SCS-CN,MUSLE and nutrient transport modules)and some fitted parameters based on the field data measured in some specific areas.When the original SWAT model was used to different climatic regions,these semi-empirical equations and fitted parameters would limit the accurateness of SWAT to track specific peak runoff,sediment and nutrient loss due to the differences in land management practices,geomorphological features and soil textures.In this study,some of semi-empirical equations and parameters were modified based on large field measurements conducted in the Loess Plateau.The modified SWAT model improved the predictions of surface runoff,sediment,and nutrient loss for the studied basin.This study expended the application of original SWAT model in the Loess Plateau and similar regions.The modified SWAT model can be used as a decision-making management tool to design appropriate management strategies for stakeholders to reduce runoff and soil loss in the studied region;the model can also be used for fertilizer management in agricultural fields to reduce NPS pollution.

Data availability statement

Data available on request due to privacy/ethical restrictions.

Conflict of interest statement

Authors declare no conflict of interest.

Acknowledgements

This research was supported by the Strategic Priority Research Program of Chinese Academy of Sciences,China(No.XDB40000000)and the National Natural Science Foundation of China,China(No.41571130082).