Assessment of Impacts of Climate Change on Tile Discharge and Nitrogen Yield Using the DRAINMOD Model

Abstract

The detrimental impacts of agricultural subsurface tile flows and their associated pollutants on water quality is a major environmental issue in the Great Lakes region and many other places globally. A strong understanding of water quality indicators along with the contribution of tile-drained agriculture to water contamination is necessary to assess and reduce a significant source of non-point source pollution. In this study, DRAINMOD, a field-scale hydrology and water quality model, was applied to assess the impact of future climatic change on depth to water table, tile flow and associated nitrate loss from an 8.66 ha agricultural field near Londesborough, in Southwestern Ontario, Canada. The closest available climate data from a weather station approximately 10 km from the field site was used by the Ontario Ministry of Natural Resources and Forestry (MNRF) to generate future predictions of daily precipitation and maximum and minimum air temperatures required to create the weather files for DRAINMOD. Of the 28 models applied by MNRF, three models (CGCM3T47-Run5, GFDLCM2.0, and MIROC3.2hires) were selected based on the frequency of the models recommended for use in Ontario with SRA1B emission scenario. Results suggested that simulated tile flows and evapotranspiration (ET) in the 2071–2100 period are expected to increase by 7% and 14% compared to 1960–1990 period. Results also suggest that under future climates, significant increases in nitrate losses (about 50%) will occur along with the elevated tile flows. This work suggests that climate change will have a significant effect on field hydrology and water quality in tile-drained agricultural regions.

1. Introduction

Field agriculture is a significant land use in many areas of the world and has been documented as a contributor of non-point source pollution to surface waters. In a recent study in Uruguay, [1] researchers found a strong correlation between total phosphorus (P) and agriculture land use in a river. In central Asia, remote sensing was used to analyze effects of land use on surface water contamination [2]. Although the relationships were not straightforward, they did find that cropland was a significant contributor to surface water organic pollution. In another study in Europe, [3], analyses of nutrients and other components in river water found that the proportion of arable land in the watershed did affect water pollution levels. To further investigate agricultural land use as a non-point source of pollution, the hydrologic pathway of a pollutant from field to surface water body must be identified and quantified. This study focusses on quantifying tile flow and accompanying nitrate load to surface waters under future climate conditions of the Great Lakes region of North America.

Subsurface (tile) drainage is considered an effective and important practice in field agricultural water management to facilitate farm operations, especially in cold climatic regions. In Ontario, Canada, more than 50% of its arable land is under artificial drainage systems [4]. Tile drainage increases productivity of land by increasing crop yields, but also impacts the environment by degrading water quality. Indeed, it is widely recognized that nitrogen (N) and P leaching results in elevated nutrient concentrations in surface waters representing both environmental and economic (fertilizer loss) impacts. High levels of nitrate may even remain in tile drainage several years after nitrogen fertilizer reduction ([5,6]). Hence, to reduce contamination of surface waters, the impacts of changing cropping systems that may, for instance, occur under future climate regimes are essential. Although climate change is anticipated to play an important role in subsurface drainage, it is not possible to characterize the potential impacts of climate change on field hydrology and/or water quality in agricultural watersheds using field observations. However, modeling of these potential impacts does have potential to at least give acceptable trends in future tile flow and associated nitrate loss. For instance, [7] used the Root Zone Water Quality Model (RZWQM) to estimate the impact of agricultural management systems adaptations on gaseous and drainage nitrogen (N) losses in Iowa, USA. They found that the optimal N rate to minimize loss and maximize maize production was 120 kg N ha⁻¹. Using DRAINMOD, researchers [8] were able to suggest when to block tile outflow to keep the groundwater table at an acceptable depth in central-western Poland. In Illinois, [9], assessed impacts of changing atmospheric conditions using RZWQM along with Support Decision for Agrotechnology Transfer (DSSAT). They found that increasing rainwater nitrate concentration had a moderate impact on nitrate loss in drainage tiles.

In Canada mainly watershed-scale studies on impacts of climate change on water resources have been completed by [10,11,12,13,14]. In their study on a watershed in Quebec, [10], simulated future streamflows finding a slight decrease in annual runoff. In another study, [11], used a weather generator to recommend local water resource management adaptations for a watershed near the present study’s site in Ontario. Using streamflow data and the LARS-WG weather generator, the Soil and Water Assessment Tool (SWAT) model was used by [12] to show that streamflow in a southern Ontario watershed could increase by 12% compared to base period 1961–1990. Also, in southern Ontario, [13] used model simulations to estimate the increase in winter streamflow in several large watersheds. A future assessment of groundwater nitrate concentrations in a sub-watershed of southern Ontario was conducted by [14]. They found that an agricultural BMP was very effective in reducing nitrate in groundwater under a crop rotation system.

In the Great Lakes region, where this study was conducted, there have been several recent studies on impacts of climate change on tile flow and nutrient losses at the watershed scale. The WEPP-WQ model was used by [15] to estimate future N and P losses in two small watersheds finding increases in losses of both nutrients. In another study [16], used DRAINMOD to simulate future tile flow in western Lake Erie basin. They found an average decrease of about 9% in subsurface drainage and recommended controlled drainage to retain more water in the soil profile as a BMP for crop production. In the same basin [17], using the SWAT model found an increase in subsurface drain flow but a decrease in P load from drains. On the other hand, relatively few field-scale studies, which use field data to first calibrate and then validate a tile-flow model to assess climate change impacts have been attempted in the region. Field-scale studies have the advantage over watershed-scale studies of a better-defined source area with associated soil and cropping properties for assessing impacts of future climate on tile flow and accompanying nutrient losses. Using DRAINMOD [18], examined impacts of different tillage practices on future tile flow and nitrate loss in eastern Ontario. They estimated greater future nitrate loss under no-tillage than conventional tillage. At the same field site as the work presented here, The SWAT model was applied to future climate data similar to the approach used in this study [19]. In addition, their results indicated a shift in seasonal water balance and an associated increase in annual P losses by year 2100. Main differences between this study and [19] are use of a different model and modeling approach and analysis of N loss in this study as opposed to P loss in the [19] study.

The DRAINMOD model was developed by [20] as a process-based distributed field-scale model. Many studies have assessed and applied DRAINMOD including, to name a few [21,22,23,24]. The DRAINMOD model has been extensively revised and updated over the years and is still undergoing improvements [25,26]. The main reasons why it was chosen for this study are three-fold: it includes impacts of winter-season processes on soil hydrology along with soil nitrogen dynamics and it has been previously tested and proven satisfactory in Ontario and other cold regions (e.g., [18,27,28,29,30,31,32,33,34,35,36,37,38]). For information on the structure and details on the applied processes in DRAINMOD see publications by [22,39].

Previous work at the same field site as used in this study includes analyses of P transport by [19,40,41,42,43,44]. As well, [45] examined the contribution of preferential flow to tile drainage.

Based on the preceding literature review, this study was undertaken to address the paucity of field-scale studies of how tile flow and nitrate loss may change under future climate regimes in the Great Lakes region. As climate changes there is the potential for shifts in agricultural practices, which may lead to increased nitrate leaching [46]. Hence more studies such as herein are needed to establish a baseline nitrate loading using current field crop practices. In a recent review by [47], it was noted that it is important to improve field-scale hydrology and water quality models under cold climates found in Canada and elsewhere in high-latitude regions. Studies that apply field data to calibrate and validate models, as done in this study, will lead to improvements in model accuracy in the future.

The objectives of this study were: (1) to calibrate and evaluate the latest version of DRAINMOD 6.1 using data from a field site in southern Ontario; and (2) to assess the impact of climate change on tile discharge and nitrogen yields at the same field site using the calibrated model. Using the CGCM3T47 with SRA1B emission scenario, we estimate the impact of future climate on tile discharge and its accompanying nitrate loading for the next century. This research may assist the farming community in developing adaptation strategies to minimize negative impacts of tile drainage systems on water resources.

2. Materials and Methods

2.1. Description of Study Area

The study site (LON) is in southern Ontario, Canada (UTM 472219 E and 4767583 N) near Londesborough. The study was conducted on an 8.66 ha agricultural field with both overland flow and tile drainage restricted to the study field. The site was under a reduced overland tillage (RT) and data collection spanned a range of years and therefore experienced a range of climatic conditions. Tile drainage systems in the field are systematically drained at 90 cm depth through 10-cm diameter laterals (13.5 m spacing) that connected to a larger main tile (20-cm diameter) that exits at the edge of the field. The contributing area to surface drainage within the field is 7.79 ha.

Long-term average annual precipitation measured at Blyth, Ontario, the closest Environment Canada weather station to LON site, is 1247 mm and long-term mean annual temperatures are 7.2 °C (https://climate.weather.gc.ca/historical_data/search_historic_data_e.html). Monthly temperatures vary seasonally across the year, with warm summers and cold winters. Daily maximum and minimum temperatures are variable throughout the year; however, daily mean air temperatures are generally below freezing between December and March. Due to the climate in Ontario, spring planted crops are seeded in early May, whereas fall seeded crops (winter wheat) are planted in September. Harvests are completed by early August (wheat), mid-September (soy) or November (corn).

The topography of the site is gently undulating with slopes ranging from 0.2 to 3.5% [48]. The field consists of soils from the Perth Clay Loam association [49]. These soils developed on clay loam glacial deposits with imperfect drainage. Soil samples collected at the site suggest that the average texture in the top 15 cm is silt (clay 0.9 ± 1.3%, silt 75.7 ± 2.1%, sand 23.4 ± 3.3%) [48].

Table 1. Some soil Properties at Londesborough field site.

Soil Properties
Depth (cm)	0–10	10–20	20–35	35–150
Soil hydraulic conductivity (cm/h)	0.26	0.07	0.2	0.26
Soil pH	7.7 ± 0.3
Organic matter (%)	4.1 ± 0.7

gives some basic soil properties for the LON site.

2.2. Model Formulation

The DRAINMOD model requires different inputs including air temperature and precipitation data, infiltration parameters, soil properties, crop information, and drainage system parameters. These data were used to calculate the relationships among the drained volume and depth to water table (WTD), and relationships between WTD and maximum steady upward flux. Infiltration was simulated using the Green–Ampt equation. Within the model, a soil moisture retention curve for the soil above the tile was used to calculate infiltration parameters as a function of WTD.

The model provides options to the user to use observed ET data or apply daily maximum and minimum temperatures to calculate ET using the Thornthwaite equation. In this study, Thornthwaite approach has been used. Details related to crop management and timing for the study site are given in

Table 2. Management practices considered in the study.

Crop	Date	Activity Code
Soybeans	Mid-October 2010 (assume 15 October 2010)	Harvest
Winter Wheat	Mid-October 2010 (Considered 16 October 2010)	Plant
	April 11, 2011	Red Clover Air Seeded
	Late April 2011 (Considered 25 April 2011)	Fert App
	Late July 2011 (Considered 25 July 2011)	Harvest
Grain Corn	5 October 2011	Fert App
	5 October 2011	Spray
	11 November 2011	Tillage
	18 April 2012	Spray
	25 April 2012	Tillage
	26 April 2012	Plant
	26 April 2012	Fert App
	23 May 2012	Fert App
	9 June 2012	Spray
	21 October 2012 (estimate as no rain that day)	Harvest
Soybeans	9 November 2012	Tillage
	3 May 2013	Plant
	24 May 2016	Spray
	25 September 2013	Harvest
Winter Wheat	27 September 2013	Plant
	19 April 2014	Plant CC
	26 April 2014	Fert App
	1 June 2014	Spray
	7–10 August 2014	Harvest
Grain Corn	9 October 2014	Fert App
	9 October 2014	Spray
	17 October 2014	Tillage
	19 April 2015	Spray
	27 April 2015	Tillage
	28 April 2015	Plant
	28 April 2015	Fert App
	May 15/15 (estimate)—check weather station	Spray
	25 May 2015	Fert App
	15 July 2015	Spray
	17 October 2015	Harvest
Soybeans	9 November 2015	Tillage
	3 May 2016	Plant
	24/05/2016 (estimate)—check weather station	Spray
		Harvest

.

Many of the input parameters for DRAINMOD are transferable between sites from nearby cold regions. Hence, parameters related to snow melt, ice content etc. were obtained from the literature, such as [29,30] and [18]. Crop-related constants including nitrogen uptake and transformation factors and organic matter dynamics were taken from [5,31,33,34,50,51].

The lower boundary at 1.2 m was assumed to be impermeable and at a constant soil temperature of 7 °C, which is approximately the long-term average air temperature [52].

The model requires initial conditions such as concentrations of NO₃-N and ${ \mathrm{NH}}_4^{+}\mathrm{N}$, given in

Table 3. Initial nitrogen concentrations in four soil depth ranges at the study site.

Initial NO3-N concentration in soil (mg L−1)
0–0.15 m	5.17
0.15–0.30 m	2.05
0.30–0.60 m	1.62
0.60–1.20 m	1.58
${\mathrm{Initial} \mathrm{NH}}_4^{+}$-N concentration in soil (mgl−1)
0–0.15 m	3.22
0.15–0.30 m	2.62
0.30–0.60 m	1.74
0.60–1.20 m	1.41
${\mathrm{NH}}_4^{+}$ sorption distribution coefficient (cm3g−1) ${\mathrm{K}}_{\mathrm{d} }$
0–0.15 m	2.17
0.15–0.30 m	2.32
0.30–0.60 m	2.85
0.60–1.20 m	3.61

.

2.3. Model Calibration

To calibrate and assess the DRAINMOD model, measurements of tile discharge collected in the field from 2011–2016 were divided into two parts: data from the first two years (fall of 2011–summer of 2013), were used to calibrate DRAINMOD, and the remaining two years of the data (fall of 2013–early winter of 2016), were used to validate the model. The chosen calibration parameters were manually adjusted by minimizing differences between observed and simulated tile flow at first visually and then statistically. The most sensitive parameters were selected for calibration based on previous work using DRAINMOD [18,27,28,32,34,35,36,50]. As in other studies, hydrology was first used to calibrate DRAINMOD using K_sat of soil profile and restrictive bottom layer and maximum surface storage which controls water runoff (

Table 4. Parameters Used in Calibration of DRAINMOD model.

Parameter			Value
Saturated hydraulic conductivity (4 layers, cm h−1)			Varies
Vertical hydraulic conductivity of restrictive layer (cm h−1)			0.0025
Maximum surface storage (mm)			1.6
Rooting depth (cm)	Month	Day
	1	1	5.0
	5	5	5.0
	6	15	24.0
	7	15	36.0
	7	30	45.0
	8	31	45.0
	10	15	5.0
	12	31	5.0
Mineralization rate (d−1)			NT:0.00002 CT:0.00003
Maximum Nitrification rate (${\mathsf{\mu }\mathrm{gNg}}^{-1}{\mathrm{d}}^{-1})$			NT:9.5 CT:7.7
Nitrification optimum temperature (°C)			20
Nitrification half saturation constant (${\mathsf{\mu }\mathrm{gNg}}^{-1}{\mathrm{d}}^{-1})$			90
Maximum denitrification rate (${\mathsf{\mu }\mathrm{gNg}}^{-1}{\mathrm{d}}^{-1})$			NT:4.0 CT:3.0
Nitrification empirical shape factor			0.5
Nitrification Optimum temperature (°C)			23.0
Nitrification half saturation constant (${\mathsf{\mu }\mathrm{gNg}}^{-1}{\mathrm{d}}^{-1})$			30.0

).

Simulations were initially evaluated using visual comparisons between the observed and predicted values. These were subsequently assessed using various statistical parameters. The statistical parameters used were the coefficient of determination (R²), the percent bias (PBIAS) [53] and the Nash-Sutcliffe efficiency (NSE; [54]). These parameters were calculated according to:

$$R^2={\left[\frac{{\displaystyle \sum _{i=1}^{i=n}(O_i-\overline{O})}(P_i-\overline{P})}{\sqrt{{\displaystyle \sum _{i=1}^{i=n}{(O_i-\overline{O})}^2}}\sqrt{{\displaystyle \sum _{i=1}^{i=n}{(P_i-\overline{P})}^2}}}\right]}^2 (0 \le R^2 \le 1)$$

(1)

$$PBIAS=\frac{100\hspace{0.17em}{\displaystyle \sum _{i=1}^{i=n}{O}_i-P_i}}{{\displaystyle \sum _{i=1}^{i=n}{O}_i}}$$

(2)

$$NSE=1-\left[\frac{{\displaystyle \sum _{i=1}^{i=n}(O_i}-P_i{)}^2}{{\displaystyle \sum _{i=1}^{i=n}(O_i}-\overline{O}{)}^2}\right]$$

(3)

$\mathrm{where},$ n is the total number of compared values, $O_i$ is the ith observed value, $\overline{O}$ is the mean of observed values, $P_i$ is the ith predicted value, $\overline{P}$ is the mean of predicted values.

2.4. Climate Change Simulations

Global Climate Models (GCM’s) or Regional Climate Models (RCM’s) are used to generate future climate data needed to assess the possible impacts of future climate changes on drainage and runoff in a watershed. These circulation models are based on physics and provide accurate predictions under different greenhouse gas emission criteria groups, defined by Intergovernmental Panel on Climate Change (IPCC) in their Special Report on Emission Scenarios (SRES) [55] The Ontario Ministry of Natural Resources and Forestry, Canada provided a facility through its web portal Aquamapper (http://climate/aquamapper.com/) to generate future climate data using up to 28 GCMs and RCMs under three gas emission scenarios of SRB1, SRA1B, and SRA2, and SRA2, used in future climate impact studies [56]. Emission scenario SRB1 describes a converging future world with constant population experiencing rapid structural changes towards an economy of service and information by introducing clean and resource-efficient technologies; whereas SRA1B describes a future scenario with rapid economic growth where population peaks at the middle of the century, and rapid introduction of new energy efficient technologies while development is being balanced across energy sources. A highly heterogeneous world with continuously increasing population with a fragmented and slower processing regional economic growth is described by emission scenario SRA2 [57,58].

For this study, we selected the emission scenario presented by SRA1B, which has been recommended for most of the 12 regions of Ontario considered in [23]. SRA1B is also considered to be one of the scenarios across Ontario by Ontario Climate Change Data Portal [59]. This emission scenario assumes the same level of socio-economic–technological growth throughout the region with a rapid economic growth, introduction of efficient new technologies, and has a balance among all energy sources. However, with variable grid sizes of hundreds of kilometers from model to model [60,61], the model predictions may lack the precision needed for smaller sites [62,63,64].

Since the high coarseness of the spatial (45 km) and temporal (daily) resolution of the GCMs to represent the physical processes of convection, land atmospheric interactions, and especially in predicting future rainfall intensity-duration-frequency characteristics, downscaling may be needed for obtaining a more accurate picture of the future climate scenario [61,64,65,66], which we have not attempted in this study.

For this study site near Londesborough, we selected Blyth (43°43′, 81°23′, Elevation = 350.5 m), the closest weather station presented in Aquamapper, to generate future predictions of daily precipitation, maximum and minimum air temperatures required to create the weather files for DRAINMOD. Out of the 28 models in Aquamapper, three models (CGCM3T47-Run5, GFDLCM2.0, and MIROC3.2 MedRes) were selected based on the frequency of the models recommended for use in Ontario with SRA1B emission scenario. The base data from 1971 to 2000 was used to generate the climate data for 2011–2040, 2041–2070, and 2071–2100.

DRAINMOD was run using the future climate data described above and weather data from 1960–1990 as a comparison. The DRAINMOD input parameters from the calibration exercise, also described above, were used in all simulations. By using the pre-existing soil properties and agricultural management system, the future simulations should yield differences between past and future tile flow and nitrate loading only due to climate change.

3. Results and Discussion

3.1. Model Calibration

Several trial-and-error runs were performed by adjusting K_sat of soil profile and restrictive bottom layer and maximum surface storage. The model results were the most sensitive, in order, as given in Table 4. Crop rooting depth was also varied to calibrate DRAINMOD (Table 4). Following the soil hydrology calibration, various nitrogen factors were determined using the same approach. Calibrated values are listed in Table 4, which are similar to those given in previous work (e.g., [21,33,35,51]).

3.2. Model Evaluation

3.2.1. Field-Scale Tile Discharge

Figure 1 and Figure 2 show time series graphs of observed and simulated daily tile discharge and precipitation during the calibration and validation periods. Figure 1, showing the calibration period, indicates that measured tile-flow events correspond well with precipitation. DRAINMOD tile-flow peaks almost consistently match the timing of major precipitation events; however, it appears to underestimate more often than overestimate some peak flow values. Please note that DRAINMOD correctly shows no tile flow when precipitation occurs during dry soil conditions in the summer season. A similar pattern is shown in Figure 2 during the validation period. DRAINMOD does a very good job of when tile flow occurs.

The model performance during the calibration and validation stages at daily and monthly time steps are shown in

Table 5. Monthly and daily calibration and validation statistics.

Statistical Index	Monthly		Daily
	Calibration	Validation	Calibration	Validation
R2	0.78	0.74	0.65	0.54
NSE	0.76	0.70	0.62	0.62
PBIAS (%)	1.2	1.4	1.1	1.4

.

Observed and simulated drain discharge were in close agreement at both daily and monthly time scales during the calibration period. The statistical values given in Table 5 for the monthly time interval are within the acceptable ranges for both NSE (≥0.65) and PBIAS (≤ ±15%), as suggested by [67]. The PBIAS values are very similar for monthly and daily time intervals.

3.2.2. Field-Scale Nitrogen Yields

After DRAINMOD was calibrated and validated to successfully simulate tile discharge, a series of nitrogen simulations were done to calibrate the nitrogen component of the model. Ammonium yields were excluded from the model evaluation because both field measurements and model simulations exhibited very small $N{H}_4^{+}$-N yields in tile drain effluent.

The calibration and validation results for daily and cumulative NO₃-N yields are shown in Figure 3. In general, observed and predicted NO₃-N yields were in reasonable agreement. During the calibration period, there were two major nitrate-loss events corresponding to significant tile-flow events. The June 2012 nitrate loss occurred during the first major precipitation event following fertilizer application in the spring. However, the major loss of nitrate in April 2013 was not associated with a particular fertilizer application but followed a winter season after corn was grown (Table 2). Figure 3 shows that there is a strong linear relationship between simulated and observed with the coefficient of determination only slightly lower during the validation period than calibration.

3.3. Climate Change Analysis

3.3.1. Comparing Past and Future Climates

Figure 4a,d compares annual measured precipitation and DRAINMOD-estimated evapotranspiration for 1960–1990 to three different climate-model estimates of precipitation and DRAINMOD-estimated evapotranspiration for three future time periods. As well, Figure 4b,c compare past-measured and future-estimated daily average maximum and minimum air temperatures using the three different climate models. Please note that both precipitation and evapotranspiration increase in the future as determined in other studies. There are small differences between the three climate-model estimates of precipitation, but model CGCM3T47 does predict somewhat lower future temperature and corresponding evapotranspiration. Therefore, it has been selected as the scenario of future climate change assessment of tile flow and nitrate loss at the site. It represents the least temperature and evapotranspiration increase giving the most predicted future surplus water and, hence, presumably the worst-case scenario for increases in tile flow and corresponding nitrate loss.

Figure 5 shows past and future precipitation and temperature estimates from CGCM3T47 on annual and seasonal bases. Annual precipitation estimates increase with time with a difference of over 20% from the near- to far-future time periods. Precipitation is divided roughly evenly between seasons; however, summers become slightly drier and other seasons slightly wetter with time. These changes could have significant impacts on reducing crop production in summer and a greater potential for nitrogen loss through denitrification in non-cropping seasons. In terms of average daily temperatures, both minimum and maximum increase in every season almost without exception as time progresses. The summer season appears to show the greatest increase in temperatures from beginning to end of the modeled time period.

3.3.2. Effects of Climate Change on Water Balances

The DRAINMOD model was used to predict the impact of future climates on field-scale water balance (Figure 6) and the associated nitrate yields from the LON field site in Ontario.

Evapotranspiration and tile flow first decreased and subsequently increased due to increased precipitation and temperatures under future climate conditions. The comparison of simulated results for the 2071–2100 period (late century) with 1960–1990 period (historical) show that tile flow and ET were found to increase by 7% and 14%.

3.3.3. Effects of Climate Change on Tile Discharge and Nitrogen Yields

Mean seasonal and annual tile discharge for the historical and future time periods is shown in Figure 7. As noted above, the DRAINMOD-estimated annual tile discharge increased slightly during the future period (vs. historical). The winter season appears to dominate this annual increase. The increase in tile flow during winter likely occurs because the estimated maximum daily temperature rises above 0 °C in the future (Figure 5). This would result in more snow melt and more precipitation occurring as rain. In support of this postulation, [13] predicted that streamflow will increase during winter season in the future in several southern Ontario watersheds. As well, [68] found that tile flow is expected to increase based on their study of four Lake Erie watersheds.

Figure 7 also shows past and future nitrate loss as estimated by DRAINMOD. Increasing loss in winter season as time progresses matches well with tile-flow estimates. However, the steady increase in nitrate loss in spring does not correspond exactly with tile flow. Average annual estimated tile flow during spring season is variable with time. Hence, the concentration of nitrate must increase in tile water during spring to account for this steady increase in load. Other studies have revealed a range of results. For instance, [69] found very little change in future annual nitrogen loss from watersheds of various sizes in northeastern Indiana. In a study using SWAT model, [70], found that annual P loading would decrease in the future from a watershed in Lake Erie basin due to increased evapotranspiration and decreased snowfall. On the other hand, [14,18,19,68,71,72] all found that nutrient losses may increase in future in Ontario and Quebec watersheds.

4. Conclusions

This study has once again shown that DRAINMOD is an effective field-scale model for simulating tile flow in Ontario, Canada. The model calibrated and validated well in comparison to tile-flow measurements collected from a single field site. The calibrated model was used to predict future tile flow and nitrogen loss from the same study site. Tile flow is estimated to increase especially in winter in the future. This result is not uncommon when compared to other studies. However, nitrogen loss appears to be more complicated than tile flow. Although nitrate loss appears to increase consistently with increasing tile flow in the future, increases in nitrate loss through tiles does not correlate well with tile flow. Perhaps an increase in winter and spring season temperatures, especially minimum temperature in spring, is leading to enhanced nitrification and hence excessive nitrate loss during spring.