Evaluating the influence of hydrological condition on the phosphorus loads in an agricultural river basin using the SWAT model

Nonpoint source (NPS) pollution is an important source of water environment contaminants (Ongley et al. 2010). In the US, although agricultural best management practices (BMPs) have been widely implemented, approximately 50% of water pollution problems are caused by NPS pollution, mainly including sediments and nutrients from agricultural areas (Teshager et al. 2017). In China, water contamination caused by NPS pollution has become increasingly serious (Li et al. 2016a). Currently, one focus of water pollution control is phosphorus loads in water environments, as phosphorus is a limiting nutrient for eutrophication (Rocha et al. 2020). Agricultural NPS pollution in the flood season is an important factor influencing overall NPS pollution (e.g. Chen et al. 2005; Ding et al. 2020). Therefore, it is necessary to formulate effective agricultural NPS pollution prevention and control measures under different hydrological conditions.

Due to the spatial and temporal variations in NPS pollution, it is difficult to conduct regular long-term and large-spatial-scale monitoring (Easton et al. 2008). Moreover, the findings of small-scale field experiments are difficult to extend to the watershed scale (Chen & Lu 2014). With the applications of remote sensing, geography information systems and global positioning systems (3S technology) in the field of hydrology, many hydrological models have been developed and used to simulate hydrological processes and water contamination (e.g. Jha et al. 2007; Lai et al. 2011). The Soil and Water Assessment Tool (SWAT), a distributed hydrological model, has been developed to predict the water yield and the circulation of runoff, sediment and NPS pollution under different land use, soil and management conditions on a watershed scale (Liu et al. 2016a). The SWAT model has been widely used in the study of various watersheds worldwide (e.g. Song & Zhang 2012; Wang et al. 2016; Yang et al. 2016). Studies on NPS pollution conducted using the SWAT model have mainly focused on the sources of NPS pollution (e.g. Grizzetti et al. 2005; Qu et al. 2015; Zhang et al. 2017) and on understanding the mechanisms of spatial and temporal variations in pollution loads in river systems (e.g. Zhang et al. 2013; Duncan 2014).

Models are simplified numerical descriptions of real world, and the simulation uncertainty inherent in a model comes from the model structure, parameters, natural randomness, input data and calibration data (Gupta et al. 2005). For NPS pollution at the basin scale, the influences of the model structure and parameters on the simulation uncertainty have been widely explored (e.g. Arabi et al. 2007; Shen et al. 2010). Regarding the uncertainties arising from the data used for modeling, previous studies have focused on the influences of the spatial resolutions of the digital elevation models (DEMs) (Xu et al. 2016), forcing precipitation data (Schurz et al. 2019), and streamflow and water quality data used for model calibration (Chen et al. 2018). The pollution load in a river system comes from different sources. For point sources, the amount of pollutants entering a river system is relatively easy to measure. For NPS pollution, estimating the pollution loads from crop production and the livestock and poultry industries entering a river is difficult, and these loads are influenced by many factors. Some studies have shown that livestock and poultry breeding is an important cause of NPS pollution (e.g. Wu et al. 2012; Wang et al. 2018a). The livestock and poultry amount data used in SWAT modeling usually include the statistical data of different administrative regions, for which the boundaries are usually different from those of the subbasins discretized for SWAT simulation. However, the influence of this mismatching on model simulation results is still unclear.

NPS pollution processes are affected by many factors. Previous studies have mainly focused on the impacts of climate change and land use change on NPS pollution loads (e.g. Yang et al. 2017; Wang et al. 2018b; Liu et al. 2020). NPS pollution processes are also influenced by rainfall-runoff processes in river basins. However, studies exploring the characteristics of the spatial distributions of pollution loads under different hydrological and seasonal conditions are rare. The currently employed calculations of water environmental capacity and estimations of pollution loads use design flows to simplify the hydrological conditions, without taking into account the effects of seasonal changes in the hydrological conditions on the migration and transformation of water pollutants in the basin. To achieve reasonable and dynamic management of the quality of a water environment in a given basin, it is necessary to estimate the load of pollutants into the river by considering different hydrological conditions (Ide et al. 2019; Xuan et al. 2019).

The Luanhe River basin in the semiarid region of North China is a typical agricultural basin experiencing NPS pollution; in this region precipitation exhibits high inter- and intra-annual variations (Wang et al. 2015). Previous studies applying the SWAT model have mainly focused on inspecting the influences of land use and climate changes on hydrological processes (e.g. Li et al. 2016b; Yang et al. 2019) and drought assessments (e.g. Yan et al. 2013; Wang et al. 2020). Knowledge about how variations in hydrological conditions affect NPS pollution at the basin scale is very limited. The objective of this study is to explore the influences of hydrological conditions on the amounts and spatial distributions of phosphorus loads using the SWAT distributed hydrological model. Moreover, the influences of the input data on NPS pollution simulations and, more specifically, the influences of the spatial scale of the applied livestock and poultry production data on the phosphorus load estimations, are examined through model simulations using two datasets comprising different spatial administrative levels. The findings of this study are expected to provide decision-support knowledge regarding phosphorus pollution control at the basin scale under different hydrological conditions and improve the understanding of the uncertainties of phosphorus transport process simulations that occur due to the choice of input data related to NPS input.

The Luanhe River rises at the Zhangjiakou City of Hebei Province, China, flows through 27 counties of Hebei Province, Inner Mongolia and Liaoning Province and ends up in Bohai Bay in Laoting County after spanning about 880 km; eight main tributaries flow into the river. The area of the Luanhe River basin is approximately 44,750 km², and 98.2% of the total area is mountainous (Jiang et al. 2015). A total of 75.1% of the basin area is located in the Chengde City of the Hebei Province. The Luanhe River basin is in the northern region of the North China Plain in the Haihe River basin, ranging from 39–43°N and 115–120°E (Figure 1). The Luanhe River basin can generally be divided into three geomorphological units, i.e. plateaus, mountains, and plains, in terms of the different geology, geomorphology, and formation types. The plateau topography type is mainly located in the northern part of the basin with elevations ranging between 1,400 and 1,600 m. The mountainous regions are located south of the plateau and north of the plain regions, with slopes ranging between 20 and 40°. The plains are distributed in the southern region of the basin, with the longitudinal slope along the river channel ranging from 1/300 to 1/1,000 (Bi et al. 2018).

The climate of the Luanhe River basin is of the temperate continental monsoon climate type. From the upper to lower basin, the average annual temperature varies from −0.3 to 11 °C. Precipitation has strong interannual and intra-annual variability in the basin. The basin receives an average precipitation total of 560 mm. The average annual runoff of the river is 4,694 million m³, 70–80% of which occurs in the period from June to September (Shan et al. 2021). The major vegetation types in the Luanhe River basin include cultural vegetation, bush wood, plateau and broad-leaved forests, and a small area of coniferous forests and meadows (Li et al. 2015). To meet the needs of the NPS pollution simulation in the SWAT model, the soil types in the study area were generalized into seven types (Figure 2). The land use types mainly include four categories: farmlands, forests, grasslands and urban areas (Figure 3). In this study, we selected the Wulongji Hydrological station (hereinafter referred to as the WLJ station) as the outlet of the SWAT simulation. The total phosphorus data from a nearby regular water quality monitoring station operated by the Chengde Municipal Ecology and Environment Bureau, which is about 1 km downstream from the WLJ station, is also used for model calibration.

The SWAT model was developed by the Agricultural Research Center of the United States Department of Agriculture (USDA). It uses a daily time scale and is based on a GIS-based distributed watershed hydrological model with physical mechanisms. Application of the elevation data (DEM), remote sensing (RS) and geographic information systems (GIS) enable the model to perform continuous simulations and predictions of large and medium watersheds over long time series (e.g. Olivera et al. 2006; Wu et al. 2019). The SWAT model can simulate many different physical processes related to the water cycle in a basin, including water movement, sediment transport, plant growth, and nutrient migration management. The model is mainly composed of eight components: hydrology, meteorology, sediment, crop growth, soil temperature, agricultural management, nutrients and pesticides (Arnold & Fohrer 2005).

For SWAT modeling, the watershed is first divided into a certain number of subbasins according to topography. Each subbasin is further divided into hydrological response units (HRUs) according to land use type, soil type and slope. Using the built-in function of the SWAT model and the corresponding algorithm, the Luanhe River basin was divided into 93 subbasins and 335 HRUs based on the inputs of the selected land use, topography, soil, and slope maps.

The DEM of the study region was obtained from 30 × 30 m SRTM DEM data set V1 provided by the Cold and Arid Regions Sciences Data Center at Lanzhou, China. For saving computational time, the DEM was upgraded to 60 m resolution. The land use map of the Luanhe River basin in 2015 was derived from MODIS terrestrial level 3 standard data products. The soil data source for the study area was the 1:1 million soil dataset provided by the Institute of Soil Science at the Chinese Academy of Sciences. Seven soil types are considered for simulation, including brown soil, cinnamon soil, gray forest soil, castannozem soil, aeolian sand soil, meadow soil and fluvo-aquic soil. All the spatial data were converted to the same coordinate system and were used to construct the SWAT model.

The meteorological data used in this study were obtained from three weather stations during 1956–2016 and were provided by the National Meteorological Science Data Center; these data included the daily average precipitation, daily average temperature, and daily average maximum and minimum temperatures. The 2009–2013 daily streamflow data of WLJ station from the Hydrological Yearbook were used to calibrate and validate the simulated streamflow in the SWAT model. The monthly water quality data of total phosphorus were also used by model calibration.

The point source data included quantity of pollutants sourced from the discharge volume of urban sewage treatment plants, emissions from industrial and mining enterprises and rural domestic wastewater discharge. A point source input port was set for each subbasin, the acquired point source pollution source data were matched with every subbasin according to its spatial location, and the data were transformed into the format that the SWAT model required.

The NPS pollution that occurs due to crop production and livestock and poultry production is considered in the SWAT model. Both these NPS pollution types are inputted to each subbasin of the SWAT model. For crop production, the county-level planting area and fertilizer amount information were obtained from a statistical yearbook. Then, the fertilizer amount per unit area was computed at the county-level. The fertilizer input for each subbasin was assigned as the computed value of the county with the largest area in the subbasin.

For livestock and poultry production, two sets of livestock and poultry production amount data for Chengde City were used to run the model independently to examine the influence of the spatial scale of the livestock and poultry data on the model simulations. The total phosphorus export coefficients of livestock and poultry simulated in this study are set according to Wang et al. (2006) and are shown in Table 1. The first dataset is based on livestock and poultry production amount data obtained at the county-level from a statistical yearbook. These data were converted into fertilization amount per unit area and were input into the model as top dressing in the agricultural management measures of each subbasin. The second dataset is based on livestock and poultry production amount data at the town-level provided; this dataset contains many more details on the spatial heterogeneity of livestock and poultry production than the first dataset. The spatial distribution of the subbasins used in SWAT modelling and the county boundaries and town boundaries are shown in Figure 4.

Based on a literature review, seven and five parameters were selected for model calibration based on the streamflow and total phosphorus load observations, respectively. The Nash–Sutcliffe efficiency (NSE) and coefficient of determination (R²) were used to quantify the performances of streamflow simulations. The NSE and relative error (RE) were used to quantify the performances of the total phosphorus load simulations. Based on data availability, for the streamflow simulations, the 2009–2011 period was selected as the calibration period and the 2012–2013 period is used as the validation period. For total phosphorus, the model calibration was conducted over the 2011–2012 period, and observation data from the 2013 period were used for model validation. The model calibration was carried out using SWAT CUP software.

To explore the influences of hydrological conditions on total phosphorus loads, typical wet, normal and dry years need to be defined. Considering that long-term hydrological data are unavailable in the studied basin, based on the long time series of observed precipitation data spanning from 1956 to 2015, these typical hydrological years were determined by statistical analyses. The Pearson type III curve (P-III curve) was used to derive the empirical cumulative frequency of annual precipitation in the basin. Based on the fitted frequency curve obtained (Figure 5), 2013 was selected as the typical wet year (cumulative frequency: 13.6%), 2012 was selected as the typical normal water year (cumulative frequency: 50.1%), and 2009 was selected as the typical dry water year (cumulative frequency: 90.0%). The SWAT model simulations output for these three representative years were used to inspect the influences of hydrological conditions on the spatial and temporal distribution of total phosphorus loads in the basin. The intra-annual variations in precipitation in the basin are shown in Figure 6. Within a given year, there are distinct flood and dry periods, respectively. For the flood period of June–August, precipitation accounted for 66.8% of the annual precipitation.

The SWAT model was calibrated and validated for the average monthly streamflow and total phosphorus load by using the data recorded at the basin. The agreement between observed and simulated streamflow was measured using the NSE and R² values for both the calibration (2009–2011) and validation (2012–2013) periods and 2008 was used as a warm-up period (Table 2). In the calibration period, the statistics of the simulation output of the model were NSE = 0.69 and R² = 0.77; in the validation period, the statistics of the simulation results of the model were NS = 0.70 and R² = 0.73. The calibrated model can simulate the streamflow well, and the simulation results are considered reliable (Figure 7).

The five parameters that are most sensitive to total phosphorus and their value ranges and final values are listed in Table 2; these parameters were used to examine the performance of the model simulations. The simulated total phosphorus loads were compared with the observed loads. The period from 2011 to 2012 was regarded as the model calibration period, and 2013 is taken as the model validation period. The NSE and RE were selected as the indices with which to evaluate the simulation results. The total phosphorus simulation results of the two models are shown in Table 3. For the model based on the township level NPS source data, in the calibration period, the simulation results of the model showed NSE = 0.57 and RE = 9.1%, and the simulation results of the model in the validation period showed NSE = 0.76 and RE = 12.4%. The model can simulate the total phosphorus load well and the simulation results are reliable (Figure 8). The statistical results of the model obtained using county-level data performance were NSE = 0.54 and RE = 11.11% and NSE = 0.40 and RE = 31.19% during the calibration and validation periods, respectively. The simulation results of the second model are poor in the validation period (Figure 8). The simulation results of the model simulations based on town-level data show better performance and were thus used to analyze the influences of different hydrological conditions on the phosphorus loads.

Using the precipitation data of 2013 (wet water year), 2012 (normal water year) and 2009 (dry water year), the calibrated SWAT model was driven to simulate outputs of the total phosphorus pollution load in each hydrological year. The total phosphorus pollution load output in 2013 (wet water year), 2012 (normal water year) and 2009 (dry water year) were 2,400.4, 1,511.7 and 919.7t, respectively. This means that the annual total phosphorus pollution load to the river system is highest in the typical wet water year, followed by in the normal water year and then in the dry water year. Considering that the point source and NPS pollution input data being used for the three years are the same, this difference could be considered a sign of the influences of hydrological conditions on the total phosphorus loads input to the river system.

To examine the intra-annual variations in the total phosphorus loads in the three typical hydrological years, the contributions of each seasonal load to the annual load are shown in Figure 9. The contributions of the loads in each month and season were quite different among different hydrological years. The total phosphorus pollution load in each hydrological year mainly resulted from the flood season in summer. However, the contribution of the summer load to the annual loads varied significantly among the three hydrological years and were 46.04, 35.75 and 29.52% in wet water year, normal water year, and dry water year, respectively. It was also found that the contributions from different seasons become more evenly distributed when the annual precipitation decreased. To gain a more thorough understanding of the relationship between hydrological conditions and the total phosphorus load, the monthly relationships for the three years are shown in Figure 10. Positive correlations were found between the precipitation amounts and total phosphorus loads in the three years. However, the slopes of the best-fitted linear relationship are quite different, indicating that the variations in this influence were high among the three hydrological years.

The spatial distribution of the total phosphorus pollution load shows great differences among the representative hydrological years (Figures 11–13). In the dry year, the range of the total phosphorus pollution load per unit area was 0.04–9.22 kg/km². The region with the highest load was located in the limited upstream area in the northern part of the basin. The output range of the total phosphorus pollution load in the normal hydrological water years was 0.06–18.36 kg/km² and high-load regions were located in the downstream subbasins. The output range of total phosphorus pollution load in the typical wet years is 0.04–30.90 kg/km²; in this year the phosphorus loads per unit area of the subbasin along the mainstream were higher than those of the other subbasins. From the comparison of the spatial distribution of phosphorus pollution loads among different hydrological years, it is shown that the key areas making the greatest contributions to phosphorus pollution differ among these different representative years.

The pollution load in a river system comes from a variety of sources. With more precise source data, the accuracy of the pollution load simulated by the SWAT model is expected to improve. In this study, we examined whether using livestock and poultry emission data of a finer spatial unit could improve the simulation accuracy; livestock and poultry emissions are an important input source to the phosphorus load in the studied river. By comparing the total phosphorus load simulation results of the two different modeling settings, it can be seen that the simulation results of the model obtained based on town level pollution source data were better than those obtained based on county level pollution source data. Especially in the validation period, the model simulation conducted based on town level pollution source data could achieve satisfactory results, while the model simulation results obtained based on county-level pollution source data performed worse. Therefore, data of pollution sources at the town-level can reflect the spatial distribution characteristics of pollution sources more accurately than those at the county level.

NPS pollution processes are driven by many factors; subsequently, the data used in a model will influence the outcome of the pollution load estimations. Chen et al. (2016) showed that errors in soil data greatly impacted NPS pollution predictions. In addition to the recommendations described in our study, a finer spatial resolution of input data related to agricultural management practices was also recommended by Ekstrand et al. (2010) and Bossa et al. (2012). Livestock and poultry production are important sources of phosphorus pollution in China (Liu et al. 2016b). Our findings highlight the necessity of using livestock and poultry production data collected at a fine administrative level in SWAT modeling to improve the model accuracy; in particular, the model is expected to identify the key regions of focus for pollution control.

Soil erosion processes are greatly influenced by the runoff generation process at hill slope scale (e.g. Zhang et al. 2008; Yang et al. 2012); this process is an important driving factor of particulate phosphorus movement in a basin. The amount and temporal distribution of precipitation affect the soil erosion process and subsequently have impacts on NPS pollution processes. Our results indicate that the annual total phosphorus in the studied river decreased from a wet year to a dry year, and the total phosphorus amount in the wet year was approximately 2.6 times greater than that in the dry year; this difference is almost within the range of the ratio of nutrient loads modeled under a wet climate to a dry climate by Li et al. (2018) in a basin with similar climate and geographic conditions to those of the Luanhe River basin. By using a lumped monthly scale model, Du et al. (2016) also found that the NPS pollution loads of wet years and flood seasons were much higher than those of dry years and non-flood seasons in a tributary of the Luanhe River basin. Our results also indicated that the seasonal patterns of phosphorus load are also quite different under different hydrological conditions. In the representative wet year, the load in the summer season was almost half of the total annual load. These results are consistent with those of Du et al. (2014), who found that nutrient loads in a wet year and flood season were much higher than those in a dry year and non-flood season. As the annual precipitation decreases, the differences in the phosphorus loads among the four seasons decreased. When developing pollution control strategies, measures for both point-source and NPS pollution need to be designed simultaneously (Altenburger et al. 2015). High variabilities in phosphorus loads from NPS pollution that occur due to differences in hydrological conditions indicate that it is necessary to adapt dynamic management strategies for point-source pollutions to the variations in NPS pollution.

In addition, our study also found that under different hydrological conditions the NPS pollution loads present obvious spatial differences; this result is consistent with a study by Guo et al. (2020). From the obtained variations in phosphorus loads at the subbasin scale, it was also shown that the load decreases when the hydrological conditions become drier. However, the magnitude of this decrease is quite different among different subbasins. This implies that the critical source area of a phosphorus load changes with a change in hydrological conditions; this is in line with the findings of Wei et al. (2016). This phenomenon indicates that in addition to hydrological conditions and NPS input (fertilizer, livestock and poultry excrement), many geophysical properties, such as topography, soil type and land cover type, also impact NPS pollution loads. For basins such as the studied Luanhe River basin in the dryland region, the runoff generation mechanisms may be quite different between dry and wet years (Sun et al. 2017), which will bring about changes to the phosphorus cycle in the soil system. Our findings demonstrate that when identifying the critical source area of phosphorus loads for pollution control, the influences of hydrological conditions need to be considered, especially in basins located in arid and semiarid regions in which the precipitation has high variability.

Excessive phosphorus discharge is an important source of NPS pollution and a key factor in the eutrophication of water bodies. Taking the Luan River basin as an example, this study explored the impact of hydrological conditions on the total phosphorus load based on the SWAT model and studied the impact of the spatial scale of livestock and poultry amount data on the model simulation performance. The estimations of the total phosphorus load based on town-level livestock and poultry amount data were better than the estimations using the county-scale data, especially in the validation period, revealing that livestock and poultry data obtained at finer administrative management levels can effectively improve the performance of a model simulation. The contribution of each season to the annual total phosphorus load is highest in summer, and this contribution of summer rises with increasing wetness. The spatial distributions of the total phosphorus load were significantly different among the three typical years, and the critical areas for pollution prevention and control also differed among different hydrological years. All these research results have shown that interannual and seasonal changes in hydrological conditions have very important impacts on the temporal and spatial distributions of the total phosphorus load at the basin scale. Therefore, when developing pollution control strategies for such basins, it is necessary to take the high variations in precipitation into account to achieve reasonable river basin management strategies.

This research was supported by the National Key R&D Program of China (Grant no. 2018YFC0406502), the GEF Mainstreaming Integrated Water and Environment Management (P145897), the Chinese National Special Science and Technology Program of Water Pollution Control and Treatment (Grant no. 2018ZX07110006), the 111 Project (B18006).

Data cannot be made publicly available; readers should contact the corresponding author for details.