Understanding nitrogen transport in the unsaturated zone with fluctuations in groundwater depth


 Fluctuations in groundwater depth play an important role and are often overlooked when considering the transport of nitrogen in the unsaturated zone. To evaluate directly the variation of nitrogen transport due to fluctuations in groundwater depth, the prediction model of groundwater depth and nitrogen transport were combined and applied by least squares support vector machine and Hydrus-1D in the western irrigation area of Jilin in China. The calibration and testing results showed the prediction models were reliable. Considering different groundwater depth, the concentration of nitrogen was affected significantly with a groundwater depth of 3.42–1.71 m, while it was not affected with groundwater depth of 5.48–6.47 m. The total leaching loss of nitrogen gradually increased with the continuous decrease of groundwater depth. Furthermore, the limited groundwater depth of 1.7 m was found to reduce the risk of nitrogen pollution. This paper systematically analyzes the relationship between groundwater depth and nitrogen transport to form appropriate agriculture strategies.


GRAPHICAL ABSTRACT INTRODUCTION
Nitrogen fertilizer is widely used on farmland to increase crop productivity. However, the utilization rate of nitrogen fertilizer for crops (wheat, rice, and corn) is only 28-41% (Zhu & Chen ; China ). Excessive nitrogen on the soil surface infiltrates into groundwater through the unsaturated zone with precipitation or irrigation, and becomes one of the main sources of groundwater pollution (Zhang et al. a, b). LS-SVM uses different loss functions where inequality constraints are replaced by equality constraints and overcomes a huge computing burden (Wang & Hu ). So far, LS-SVM has been applied to estimate evapotranspiration and precipitation (Kundu et al. ), predict groundwater depth (Tang et al. ), and analytical chemistry (Balabin & Lomakina ). Therefore, LS-SVM was adopted to predict the groundwater depth. However, the above literature rarely discussed the transport of nitrogen quantitatively for different groundwater levels or proposed a reliable groundwater depth limit for groundwater-ecological risk assessments. In this study, based on the field-measured and collected data, LS-SVM and Hydrus-1D were applied to evaluate the groundwater depth and nitrogen transport, respectively. Additionally, the transport of nitrogen in the unsaturated zone at varying groundwater depths was simulated. Finally, the ecological threshold for groundwater depth was determined to provide a sustainable limit for contaminated sites.

Field site
The study area is located at Liangjiazi Town, Da'an City, west of Jilin (123 52 0 48″ E-45 18 0 45″ N). The annual average temperature is 5.02 C and the annual average precipitation is 422 mm. The annual average evaporation is 1,749 mm, which is 4.14 times the precipitation. The location is shown in Figure 1 and the area is approximately 600 m 2 . In the field site, the fertilize points are distributed in an 'S' type (see the right panel in Figure 1) to avoid farming errors.
The research area is a paddy field modified by saline soil and the planted rice is Jijing 88. Shallow and deep irrigation are used in the farmland during different periods of crop growth. From regreening to tillering, shallow irrigation was applied; from booting to heading and flowering, deep irrigation was applied; and then shallow irrigation was applied during the ripening period. The irrigation was stopped until the rice was harvested. According to the usual application mass of nitrogen fertilizer for paddy fields in the research area (150 kg/hm 2 -230 kg/hm 2 ), the trial plot was fertilized in two patterns: A (180 kg/hm 2 ) and B (220 kg/hm 2 ), which were used for calibration and testing, respectively. The proportion of base fertilizer, tillering fertilizer, booting fertilizer, and granule fertilizer was 3:4:2:1. Potassium and phosphate fertilizers of 90 kg/hm 2 and 70 kg/hm 2 were applied once in the base fertilizer. The field experiments were conducted for 150 days. The rice was raised on May 15, transplanted on June 19, harvested on October 11.
During the period of regreening, tillering, booting, ripening, and harvesting, the soil samples were collected at a depth of 0-20 cm, 20-50 cm, 50-80 cm, and 80-100 cm. The soil samples were taken to the laboratory and tested for concentration of nitrate nitrogen, ammonium nitrogen, and moisture content.

Groundwater depth modeling
Groundwater depth data in the study area are estimated by LS-SVM (Levenberg ; Marquardt ). In the modeling process, the groundwater depth is simulated by mapping the input space to the high-dimensional space nonlinearly. Linear regression is then performed and the regression function by a vector form is represented as f(x) ¼ wΦ(x) þ b, where w is the plane weight vector, b is the threshold value, Φ(x) is the nonlinear transfer function which maps the input vectors to the high-dimensional space. In the LS-SVM algorithm, the corresponding optimized problem is shown in Equation (1): where J(w, ξ) is the objective function; c is the equilibrium constant; and ξ i is the relaxation variable. The Lagrange function is shown as follows: where δ i is the Lagrangian operator, and the differentiation of w, b, ξ, δ with the Lagrange function is described as follows: The kernel function is defined as K( and the regression function is described in Equation (5): In the LS-SVM algorithm, the radical basis function (RBF) is taken as the kernel function, which is pivotal for determining the mapping function and feature space. First, data for 180 kg/hm 2 were used to calibrate the parameters in the transport model of nitrogen. Then, data for 220 kg/hm 2 were applied to test the calibrated transport model. Finally, consistency between modeling and observation values was assessed using the evaluation indexes in calibration and testing processes.

Flow transport
The process of flow in the unsaturated zone is described by the modified Richard equation, which varies with soil matrix potential, as shown in Equation (7) where θ is the volume moisture content; t is time; h is the pressure head; k is the unsaturated hydraulic conductivity; k A ij is the dimensionless expression of anisotropic tensor; x i is the space coordinate; and S is the sources and sinks of root uptake. The soil moisture characteristics parameters are obtained using the van Genuchten model. The expressions are as follows: where θ s and θ r are the saturated and residual water content; n, m, and l are the empirical fitting parameters, m ¼ 1 À 1/n; α is the reciprocal of intake pressure, α ¼ 1=h b ; and k s is the saturated hydraulic conductivity. with measured moisture data and fitted parameters were obtained. The related parameters including the saturated water content (θ s ), the residual water content (θ r ), the empirical fitting parameters (α, n, l) and the saturated hydraulic conductivity (k s ) are shown in Table 1.

Solute transport
During the leaching process of nitrogen in the unsaturated zone, excessive organic nitrogen undergoes biological transformation processes of mineralization, biological immobilization, nitrification, and denitrification, and produces inorganic nitrogen of ammonium nitrogen and nitrate nitrogen. In addition, chemical processes covering volatilization, adsorption-desorption, fixation and release of nitrogen occur in leaching (Tillotson et al. ). The transformation process is described by zero-order and firstorder dynamic equations in Hydrus-1D and is shown as follows: where c 1 , c 2 , and c 3 are the concentration of urea, ammonium nitrogen, and nitrate nitrogen, respectively; q i is the soil moisture flux in i direction; s is the mass concentration of ammonium nitrogen in soils (adsorption mass), s ¼ k d c 2 ; k d is the adsorption coefficient of ammonium nitrogen; c N is the organic matter content of soils; k 0 , k 1 , k 2 , k 3 , and k 4 are the mineralization rate of organic matter, hydrolysis rate of urea, nitrification rate (k 2W and k 2s are the nitrification rates in solid and liquid phases), denitrification rate, and fixation rate of biology, respectively; ρ is the bulk density of soils; S is the term of root water uptake; and D W ij is hydrodynamic dispersion tensor in all directions.
Considering the concentration of nitrate nitrogen is much greater than that of ammonium nitrogen in the unsa-  Table 2.

Root uptake
The driving force of root uptake is mainly derived from crop transpiration. The mass of root uptake is closely related to the distribution of crop root. In this research, the root uptake model (Feddes model) (Reicosky & Ritchie ) and the root growth model (Simunek et al. ), which are calculated in Hydrus-1D, are used and shown as follows: where a(x, h) refers to the function of soil matrix potential or suction head; b(x 1 ) refers to the distribution of the function of relative root density; L r refers to the root depth; and T P refers to the potential transpiration rate.  Table 3. The corresponding function of b(x 1 ) was provided by the database in Hydrus-1D as follows: where x 1m refers to the vertical maximum distance of the root distribution; x Ã 1 represents the horizontal coordinate of the vertical maximum root density; P x refers to the empirical parameter indicating the asymmetry of the root in the longitudinal direction. The root nitrogen uptake is the product of the water uptake mass and nitrogen concentration of the calculation point.
In addition, the parameter L r in the calculation process of root uptake was determined by the root growth module in Hydrus-1D, which considered the actual growth situation. The model assumes that the root depth is the function of the maximum root depth and the root growth coefficient (Simunek et al. ) as follows: where L r (t) represents the root depth which changes with time; L m refers to the maximum root depth; f r (t) is the root growth coefficient; and L 0 and r are the initial root depth and the root growth rate, respectively. The rice crop in this study had a certain depth of root in transplanting. The initial root depth was 15 cm and the maximum root growth depth was 50 cm.

Boundary and initial conditions
In the simulation process, the initial moisture content and the concentration of nitrate nitrogen and ammonium nitrogen were evenly distributed in the horizontal direction. In the field site, basin irrigation was used and the depth of irrigation was 10 cm. The moisture content of soil surface was saturated after the rice transplanting and the moisture content at other points was defined by linear interpolation.          Nitrate nitrogen also leaches into groundwater more easily and the concentration in shallower layers is lower than in deeper layers.
The values of RMSE are lower than 1.5 mg/kg and the values of R are greater than 0.8, as shown in Table 6.
These results indicate that the Hydrus-1D model is viable for exploring the migration and transformation of nitrogen.
In different layers, the structure of soil pores, connectivity, and physicochemical properties of salinity, pH, and moisture affect the transport and transformation of nitrogen.
These related parameters change with crop growth and could not be modified in Hydrus-1D. The fitting accuracies therefore vary with different layers, as shown in Table 6.  Considering the uncertainty of the data and parameters, the prediction interval of t-distribution is applied to improve the accuracies of the simulation results. The calculation result of 95% confidence interval is obtained as in Nie et al. () and is shown in Table 7.
Based on the dynamic range of groundwater depth shown in Table 7, Plan 1 (the upper limit of the predictive interval) and Plan 2 (the lower limit of the predictive   Upper limit (m) 6.39 6.49 6.69 6.71 6.75 6.81 Groundwater depth (m) 6.31 6.32 6.36 6.42 6.47 6.35  respectively. The concentration of nitrate nitrogen in Plan 1 is slightly higher than that in Plan 2. It can be seen that the concentration of nitrate nitrogen at 100 cm is slightly affected by the groundwater depth of 5.48-6.47 m.
The cumulative leaching loss of nitrogen was calculated from the modeling results in Table 8

ACKNOWLEDGEMENT
The research reported here was funded by the National Key R&D Program of China, grant number 2018YFC1800400.