SWAT modeling of hydrology, sediment and nutrients from the Grand River, Ontario

GRW

Grand River Watershed

SWAT

Soil and Water Assessment Tool

SUFI2

Sequential Uncertainty FItting

SWAT-CUP

Calibration and Uncertainty Procedures

TP

Total phosphorus

TN

Total nitrogen

Temperate freshwater lakes are most commonly phosphorus limited (Schindler 1974, 1977), and during the 1960s, excessive phosphorus was responsible for harmful algal blooms (HABs) in Lake Erie (United States Environmental Protection Agency (USEPA) 2011). Although HABs continued through the 1970s, a 60% reduction of phosphorus by 1985 has resulted in a corresponding 89% decrease in Aphanizomenon (USEPA 2011). However, Microcystis blooms returned in the 1990s and have continued to date, especially in the western basin (USEPA 2011). The largest algal bloom in the lake's recorded history occurred in 2011, although the 2015 bloom was the most severe on record (based on biomass over peak 30-day period) (National Oceanic and Atmospheric Administration (NOAA) 2015). Agricultural practices coupled with changing climatic conditions were identified as the main drivers of the eutrophication of Lake Erie (Michalak et al. 2013).

The Grand River watershed (GRW) covers 25% of the Canadian land area within the Lake Erie watershed (and 10% of the total Lake Erie basin) while contributing around 5% of total phosphorus (TP) loading to the lake (Scavia et al. 2014). Over the years, water quality in the Grand River (based on dissolved oxygen, TP and suspended solids) has been consistently rated as poor when measured at the Dunnville Dam monitoring site near the river's mouth (Cooke 2006).

The GRW has been greatly altered over the past 150 years, with large areas converted to agriculture, growth in population density, and modifications of the stream network. Due to their dendritic patterns and high width:depth, upland headwater streams are important sites for biogeochemical transformation of nutrients, making them critical in controlling nutrient export downstream (Alexander et al. 2007). Changes in land use have resulted in changes to the natural meandering features of many of these smaller streams, altering hydrology and ecological functions with ramifications for sediment and nutrient transport, water quality, and biological characteristics, as well as creating negative economic and social consequences (Alexander et al. 2012).

Modeling is an important tool, widely used to enhance our understanding of physical, chemical or biological processes occurring at the local or watershed scale. The Soil and Water Assessment Tool (SWAT) is a freeware that is a continuous-time, semi-distributed, physically based (i.e. deterministic) watershed model operating on a daily time step (Arnold et al. 1998). SWAT was developed by the United States Department of Agriculture – Agricultural Research Service (USDA-ARS), incorporating a number of previous models into its structure. Due to SWAT being physically based, it provides the unique opportunity to simulate the hydrology and water quality of ungauged streams and to quantify the relative impacts of alternative input data on hydrology and water quality in watersheds. Today, SWAT is widely used to assess the environmental impacts of land use on water quantity and quality in small agricultural fields to continental-size watersheds with varying soil types, topography and land uses (e.g. Abbaspour et al. 2007, 2015; Yang et al. 2011; Booty et al. 2014; Huang & Xiang 2015). Furthermore, SWAT, being a continuous time model, is capable of simulating periods ranging from 1 to 100 years, providing output on daily, monthly or annual time scales.

Im et al. (2003) compared the applicability and performance of SWAT and HSPF (Hydrologic Simulation Program Fortran) to predict stream flow, sediment, and nutrient loading from the Polecat Creek watershed in Virginia (∼120 km²). HSPF was slightly better than SWAT at predicting stream flow and water quality parameters, however the researchers indicated HSPF is less user-friendly than SWAT, due to complexity of the hydrologic cycle, sediment and nutrients transport model components, and their numerous parameters that must be optimized. Zhang et al. (2016) compared SWAT and the Distributed Large Basin Runoff Model (DLBRM) to predict daily runoff in the mountainous Heihe River watershed in northwest China (∼128,000 km²). Both models performed well, but DLBRM performed better than SWAT due mainly to data requirements and constraints, differing interpolation structures, and spatial scales of the models. Yet, both models performed poorly in predicting low flows due to inability to capture impacts of snow-melting processes and the presence of frozen soils in the watershed. Golmohammadi et al. (2014) compared the applicability of three models: MIKE SHE, the physically based, distributed European Hydrological System Model; APEX, Agricultural Policy/Environmental Extender, a small-scale physically based model developed by USDA and Texas A&M University; and SWAT, to predict hydrology in the Canagagigue watershed (52.6 km²), a sub-catchment within the GRW. Their results showed that MIKE SHE was slightly better than SWAT at simulating and predicting daily and monthly stream flow at the outlet of the watershed, while APEX performed less well, perhaps due to its intended field-scale of operation. Land use may be an important determinant of model performance. For example, Bosch et al. (2011) demonstrated that SWAT is very effective when applied to agricultural and forested lands as compared to urban areas. Although SWAT has not consistently outperformed all other models, these studies have consistently shown that SWAT is robust and performs well in a variety of scenarios, environments, and scales (temporal and spatial).

The objective of this research was to build and test the applicability of a SWAT hydrological model for the GRW under current climate, land management, stream network, and physiographic conditions. The GRW is an important agricultural region in southern Ontario. The GRW is predominantly rural, but has several rapidly growing urban centers. The GRW is a major basin draining into Lake Erie, contributing to nutrient loading and likely to associated ecological problems experienced in Erie. This study describes the calibration and validation of the model for quantifying surface water hydrology in the GRW, and uses the model to predict sediment, phosphorus, and nitrogen loading from the GRW to Lake Erie. This model provides a valuable tool for predicting sediment and nutrient export from the watershed in response to interventions in the watershed or to scenarios of land use and climate change.

The GRW covers an area of 6,800 km², and is the largest watershed in southwestern Ontario (Figure 1). Spanning a length of 290 km and having an elevation differential of ∼362 m from source to mouth, the Grand River originates at the Dundalk Highlands and flows through Port Maitland into Lake Erie. Land use is dominated by agriculture (76%), with about 6,000 farms (Cooke 2006). Forested areas and urban centers (most centrally located) cover 17% and 5% of the total watershed area, respectively (Cooke 2006). In intense agricultural sub-watersheds in southern Ontario, corn is the main crop type with 25% coverage, while soybeans, other grains, and forage and fodder cover 16%, 12%, and 11%, respectively (Ontario Ministry of the Environment (OMOE) 2012). In the absence of crop coverage data specific to the GRW, these percentages were assumed to apply to agricultural areas of the GRW.

The duration of the SWAT run was set from 1996 to 2010, divided into three five-year phases: a warm-up period from 1996–2000, a calibration period from 2001–2005, and a validation period from 2006–2010. The initial testing period was used to perform a water mass balance in order to check the accuracy and stability of model outputs. The model was run on monthly and annual time scales.

A Digital Elevation Model (DEM) for Ontario was obtained via the Scholars GeoPortal network. This DEM (version 2.0.0) is a 3-D raster dataset capturing terrain elevations, with cell resolutions of 10 m in southern Ontario. The Grand River Conservation Authority's (GRCA) virtual stream network was burned onto the raster during the watershed delineation process to produce more accurate sub-basin delineations in subsequent steps. The virtual stream layer is a ground-truthed, fully connected network that represents the inferred flow through watercourses and water bodies.

A land cover shapefile was obtained from the GRCA. This land cover classification was based on Landsat 7 TM Imagery in 1999, updated in 2005. Pixel sizes are 25 m in the land cover grid. To date, SWAT has a library of 97 plant types and 8 urban land uses in its database.

Geospatial data were obtained from the Soil Landscapes of Canada's online geospatial database, maintained by the Canadian Soil Information Service (CanSIS) (http://sis.agr.gc.ca/cansis/). A soil database used in SWAT has to identify and link the physiochemical properties of the top 10 layers of the soil profile to the soil groups in the watershed. The names and characteristics of the different soil layers were downloaded from the CanSIS website in the Soil Name Table and the Soil Layer Table of SWAT.

Weather data (daily maximum and minimum temperature (°C), precipitation (mm d⁻¹)) were obtained from GRCA monitoring stations (Conestogo, Luther, Shand, Woolwich), and from Environment Canada datasets from stations located within or in proximity to the GRW (Fergus Shand Dam, Guelph Turfgrass and Guelph Turfgrass CS (combined), Waterloo Wellington A and Waterloo Wellington 2 (combined)) (http://climate.weather.gc.ca/historical_data/search_historic_data_e.html). The Global Weather Data tool was used to provide all missing data (http://globalweather.tamu.edu/). Given sparsely measured daily solar radiation, humidity and wind speed datasets, the Hargreaves' method (Hargreaves et al. 1985) of simulating potential evapotranspiration (PET) was selected for the model run. This method of inferring PET rates provided reliable estimates when compared to measured values in well-watered agricultural lands in the Midwestern USA (Jensen et al. 1990) and across 49 geographically diverse sites in the USA (Itenfisu et al. 2003). Given the close proximity of the Grand River of southern Ontario to the USA, and the applicability of the Hargreaves' method over geographically diverse sites throughout the USA, the Hargreaves' method is a reasonable approach to simulate PET in the GRW.

Daily mean flow data at Drayton and York stations were collected from the GRCA. Additional data were downloaded from Environment Canada's DataExplorer, a self-contained application with a browser and search engine for Environment Canada's HYDAT database. With the exception of York, which began collecting data in 2002, all other station data spanned the years 2001–2010.

Data on the extent and location of tile drains in Ontario were obtained from the Ontario Ministry of Agriculture, Food, and Rural Affairs, and the Ministry of Natural Resources (updated 2009). In the Grand River, a total area of 141,037 ha (or 28% of agricultural lands) are tile drained. SWAT required information on three parameters related to tile drainage in those areas where drains are located: depth to drain, time to drain soil to field capacity, and tile drain lag time. Based on common practice in similar watersheds in southern Ontario, the depth to surface drain (mm) was assumed to be 900 mm, the time taken to drain soil to field capacity was assumed to be 24 hours, and tile drain lag time was assumed to be 3 hours (Yang et al. 2013).

The GRW SWAT model was calibrated and validated using the semi-automated software, SWAT-CUP4 using the Sequential Uncertainty FItting (SUFI2) algorithm. Prior to calibration and validation, relevant parameters were examined to determine their sensitivity with respect to SWAT output. Given these parameters and suitable ranges (based on literature values (e.g. Abbaspour et al. 2007, 2015; Yang et al. 2013; Booty et al. 2014; Leon et al. 2014), direct knowledge of site measurements, or from one-on-one sensitivity analysis testing), random samples of parameter values from a multidimensional distribution were drawn by Latin Hypercube Sampling (LHS). Thereafter, the objective function was calculated for each simulation (having its own parameter set). Compared to random sampling such as Monte Carlo, LHS is iterative and can reduce sampling times and provide a 10-fold greater computing efficiency (Vachaud & Chen 2002).

The qualifiers v, r and a, are part of the parameter set-up process whereby (v__) refers to the substitution of a parameter by a value from the given range, (r__) refers to a relative change in the parameter where the current value is multiplied by 1 plus a factor in the given range, and (a__) refers to a given value being added to the parameter value. For watershed spatial parameters r or a was used while for non-spatial parameters, v was used.

SWAT delineated a total of 699 sub-basins within the GRW with a mean and median size of 961 ha and 779 ha, respectively and the first and third interquartile sizes were 505 ha and 1,293 ha, respectively. Hydrologic response units (HRUs) are a function of slope, land use and soil type and represent the functional unit of the watershed. A total of 7,699 HRUs were delineated, ranging in size from 0.01 to 1,182 ha. The mean and median size of the HRUs were 30 ha and 6 ha, respectively; first and third interquartile sizes were 1 ha and 29 ha, respectively.

In the SWAT model, evapotranspiration was 58% of total precipitation, and represented the largest water loss from the watershed. Approximately 23% of precipitation was lost as surface runoff, entering stream flow during discrete precipitation events or snow melt, and the remaining 19% percolated to the shallow aquifer. Of the water entering the shallow aquifer, most entered stream flow under base flow conditions (16% of total precipitation, thus 52% of total annual stream flow was surface runoff). The remainder was either lost to evaporation from the shallow aquifer (2% of total precipitation) or recharged the deep aquifer (1% of total precipitation).

Currently, ∼1 million persons live in the GRW, forecast to increase by 60% by the year 2056 (GRCA 2013). Municipal water supply, dewatering and commercial uses account for 71% of consumed water (109 Mm³y⁻¹), and although only 1% of the watershed is irrigated, this use represents 7% of consumption in the GRW (10.8 Mm³y⁻¹). Water stress in the GRW seems inevitable. A 1% deep aquifer recharge (as predicted by SWAT) is worrisome given the current high rate of water extraction, and anticipated growth in population. Water balance in the GRW will also be vulnerable to climate change, with higher predicted atmospheric moisture content and rainfall intensities, and a projection of more frequent storm events, with concomitant flooding, erosion, and unstable hydrographs. An increase in evapotranspiration is also expected with warmer temperatures resulting in higher water stress in crops. The GRW SWAT model provides a tool to predict how scenarios of land use change, increasing human population densities, increased water consumption, and climate change may affect water balance in the GRW.

Twenty-two parameters were initially chosen for flow calibration in headwaters. The least sensitive were excluded after two successive iterations of 500 simulations each, with 14 parameters retained. A similar process was followed for calibration to stations in the main channel (Galt, Brantford, and York), beginning with 22 parameters, and culminating in 10 parameters retained. Two additional iterations consisting of 500 simulations were run to parameterize the watershed under its specific land use, slope and soil conditions. After the calibration process, the parameters, along with their minimum and maximum values, were re-written by SWAT-CUP preparing the model for validation.

The different sub-set of model parameters proved most sensitive with respect to predicting stream flow in the headwaters (Conestogo and Eramosa/Guelph) versus the main channel and lower Grand (Galt, Brantford, and York) (Table 1). The parameters used for calibration in Galt, Brantford, and York were applied to prediction of the ungauged station at Dunnville, near the mouth of the Grand River, assuming similar sensitivity to this set of parameters.

All modeling has associated uncertainty arising from uncertainty in the conceptual model itself, the input data, and the parameters used. Uncertainty is propagated SWAT output and is represented by 95% probability distributions (95PPU). Other commonly used functions relating to model performance include: p-factor, r-factor, R², PBIAS, bR², and the Nash-Sutcliffe efficiency (NSE) value. Criteria for assessing model performance based upon these functions are reviewed elsewhere. For brevity, the better a model performs in predicting a parameter, the nearer the p-factor is to 1 (Abbaspour et al. 2007), and the smaller is the r-factor (Abbaspour et al. 2007). As performance increases, R² and, more importantly, bR² approach 1 (Krause et al. 2005), as does NSE (Nash & Sutcliffe 1970). Moriasi et al. (2007) suggested that the model simulation is satisfactory with NSE > 0.5. While a PBIAS of 0 is ideal, Moriasi et al. (2007) suggested that if the value of PBIAS ranges from 15 to 25 (absolute values), the output from SWAT is evaluated as ‘satisfactory’, a range from 10 to 15 is evaluated as ‘good’ and less than 10 is ‘very good’.

The GRW SWAT model simulated flow with a high degree of accuracy, as measured by the objective's functions (Table 2), suggesting the model is very robust and would be able to predict flow at ungauged streams at different locations in the watershed, such as at Dunnville, near the discharge into Lake Erie. The model's performance compared well against similar studies within the Lake Erie basin. Bosch et al. (2011) investigated the applicability of SWAT to six different US watersheds of varying characteristics within Lake Erie basin. In model validation, R² and PBIAS values of 0.75–0.95 and −11 to 11, respectively, were obtained for stream flow. Golmohammadi et al. (2014) compared different models in simulating stream flow in the Canagagigue Creek watershed, a sub-basin of the GRW, for the period of 1990–1998, obtaining R² and PBIAS values of 0.64 and −12.5, respectively.

SWAT predicted hydrographs matched well to measured hydrographs over the ten-year calibration and validation period (Figure 2), although the model appeared to under-predict summer flows in the main channel of the watershed. One factor contributing to this may be the effects of dams within the watershed. In addition to 28 dams owned and operated by the GRCA, over 200 dams are owned by municipalities and private landowners. The seven major dams controlled by the GRCA, including the Shand, Conestogo and Guelph dams, can reduce spring flood peaks by at least 50% in the middle and lower Grand. However, water stored in the reservoirs is released primarily in summer to maintain minimum flows downstream. This is required to supply municipal drinking water systems and to dilute effluent from wastewater treatment plants discharging into the Grand. The daily volume of water released from all reservoirs, a required input for SWAT, is not available. However, this discharge may contribute substantially to flow in certain portions of the Grand at certain times of year (GRCA 2005). Another factor compromising SWAT's predictions of hydrology, particularly in dry seasons, is water extraction. By 2000, there were over 800 active water extractions in the watershed by municipalities, farm irrigators, golf courses, aggregate producers and industry. SWAT can account for extraction given the location, day, and volume extracted either from the reach or aquifer. Again, these data are unavailable, apart from a few scattered data in GRCA reports.

A summary of monthly average discharge (over the ten-year calibration and validation period) showed an annual cycle of high spring and low summer flows corresponding to periods of high precipitation and snowmelts (Figure 3). Again, there was a close relationship between model predicted and measured monthly average discharges at the five stations for which data were available. The parameterized SWAT model was used to simulate flow at Dunnville, predicting a maximum average flow of 179 m³ s⁻¹ in March and a minimum flow of 20 m³ s⁻¹ in August. Transferring the calibrated flow parameters from Brantford and York, base flow has greatly improved with maximum base flows of around 50 m³ s⁻¹ (Figure 3).

Many of the parameters that were most sensitive in modeling flow were also sensitive in modeling total suspended solids (TSS) and phosphorus transport, including snowfall and snowmelt temperatures, base flow alpha factor, groundwater delay time, curve number, Manning's n-value for the main channel, effective hydraulic conductivity in the main channel, soil available water storage capacity, soil hydraulic conductivity, soil bulk density, maximum canopy storage (Table 1). While flow data were collected at a daily time scale, available water quality data for the GRW, such as sediment and nutrients, were collected by the Ministry of the Environment 3–5 times annually, with some years missing. Grab sampling at such a low frequency is wholly inadequate for calibration and validation purposes. LOAD ESTIMATOR (LOADEST) produces a regression model for estimating a constituent load, estimates can then be used as surrogates for observed values in SWAT calibration (Runkel et al. 2004). However, LOADEST outputs were statistically rejected due to too few sampling points. Given no other option, additional sediment and nutrient parameters that were not included in the flow sensitivity analysis were adjusted by means of manual calibration. Some of these parameters include sediment sensitive parameters, such as a linear parameter for sediment routing factor in main channels (SPCON), a channel re-entrained exponent parameter (SPEXP), a channel re-entrained linear parameter (PRF), channel erodibility factor (CH_K2), and channel cover factor (CH-COV). TP sensitive factors included phosphorus availability index (PSP), an adsorption coefficient (PHOSKD), and a soil solution P concentration (SOL-LAB-P). These parameters were not changed by more than 10% of SWAT's default values.

After calibration, simulated TSS (Figure 4(a)) at all sites were within the expected range based on long term observations by the Ministry of the Environment and illustrated by the GRCA (Wong 2011). The seasonal pattern in TSS transport was similar to the pattern for hydrology, with higher TSS in spring and lower TSS in summer, related to seasonality of physical and or biological processes in the watershed (e.g. frequency and intensity of rainfall, coupled with land use that affects erosion; seasonal changes in algae populations). Currently, there is no Provincial Water Quality Objective (PWQO) for TSS applied to streams. However, the Canadian Water Quality Guidelines (CWQG) states that the maximum increase of TSS should be no more than 25 mg L⁻¹ from base flow conditions for a short-term pulse (e.g. 1-day period), and only a maximum average increase of 5 mg L⁻¹ from base flow conditions for longer term (e.g. up to 30 days). In addition, if the base flow TSS levels are between 25 and 250 mg L⁻¹, TSS should not increase by more than 25 mg L⁻¹ during periods of high flow (Canadian Council of Ministers of the Environment (CCME) 2002). Although intended for practical purposes, these recommendations are very difficult to implement since the base flow TSS values need high frequency sampling, not three times a year as is currently practiced. Further, sampling only during months of high flow or only during months of base flow will skew any inference on sediment interactions within the watershed (Tate et al. 1999). Furthermore, base flow is highly influenced by monthly and annual precipitation rates, thus decades of samples are required before one can accurately infer the base flow TSS. With this acknowledged, if August and September simulated values are assumed to represent base flow TSS, then values for the remaining months are >25 mg L⁻¹ above base flow TSS at almost every site.

Based on a calibrated/validated hydrology and a semi-automated calibrated TSS model, simulated sediment yields from the watershed were extracted. These estimates were restricted to the main channel (i.e. at stations Galt, Brantford, York, and Dunnville), as the purpose was to model sediment export from the river to Lake Erie. In addition to the modified universal soil loss equation (MUSLE), SWAT predicts sediment yield by accounting for bed sediment erosion and deposition. Thus, sediment yield is the sum of TSS transported and bed-load sediment. Sediment yield followed a seasonal pattern with higher export in spring and lower levels in summer (Figure 4(b)), suggesting that the most important drivers of sediment yield are precipitation and land use.

Simulated TP levels throughout the GRW far exceed the provincial objective of 0.030 mg L⁻¹. The model predicts average annual TP values of 0.90, 0.34 and 0.50 mg L⁻¹ (as P) at Galt, York and Dunnville, respectively. TP values were similarly high in headwaters draining row crops on tile drained lands. These findings of elevated TP levels within the watershed agree with GRCA water quality reports (e.g. GRCA 2005, 2013; Wong 2011). Simulated TP values decreased as the river reached Dunnville. Loomer & Cooke (2011) indicated that the Grand River tends to recover as it flows toward Paris from Cambridge, perhaps because the Grand flows through a steep valley with many riffles, enhancing oxygenation. Loomer & Cooke (2011) also suggested a significant contribution of groundwater to stream flow in this area helps to dilute the nutrient levels in the river.

Most predicted TP is in organic/bound form and would not be readily available for biotic assimilation. TP yields are very high during the months of November–March across the watershed, with the majority in the form of organic phosphorus (Figure 4(c)). However, predicted soluble reactive phosphorus (SRP) increased in proportion during the low TP yield months. This may be due to increased residence time in the headwaters, which will facilitate greater nutrient cycling and decomposition of organic nutrients by benthic organisms. Predicted TP was correlated with TSS (r = 0.72) and both were generally highest during significant runoff events such as spring runoff. TP on an annual basis was higher in the highly agricultural sub-basins. However, despite some recovery in water quality in lower reaches, the Grand River's annual export of TP into Lake Erie remains high (200–300 tonnes y⁻¹). This estimate is in accord with the estimate by Shaker (2014) of 285 tonnes y⁻¹ for 2000 to 2009, based on 8–12 grab samples per year from two Dunnville stations. Dolan & Chopra (2012) estimated total P loading to Lake Erie, including 6183 tonnes y⁻¹ from non-point sources including tributaries. The SWAT model estimate suggests the Grand River contributes 3–5% of this load, and a substantial portion of the load to the Eastern Basin.

Predicted nitrate concentrations in the GRW showed less pronounced seasonality, and a shift to highest concentrations in autumn, particularly October (Figure 4(d)). However, predicted total nitrate export was still highest in spring as a function of higher total discharge. Nitrate concentrations were generally greatest in sub-basins with high densities of livestock. Overall, modeled nitrate concentrations generally exceeded the PWQO (2.93 mg L⁻¹) in the headwaters and middle reaches of the Grand River. In the spring and early summer, nitrate concentrations were highest in the middle reaches, where the river flows through the portion of the watershed with most intensive row crop agriculture. The predicted decrease in nitrate near the mouth was due to the influx of groundwater, and to a lesser degree biological processes including assimilation into biomass and denitrification.

This study is among the few studies that applied SWAT or a similar ecohydrologic model to investigate hydrology, and sediment and nutrient loading from agricultural lands in Southern Ontario. A sub-basin-scale hydrologic representation of the highly agricultural GRW was constructed for the years 1996–2010 using the conceptual based, continuous time, distributed parameter SWAT model. The model effectively predicted hydrology at gauged stations within the GRW, and because the model is physically based, calibrated values from gauged stations can be transferred to predict water flow and quality in streams that are ungauged or that have limited available data. Predictions of TP export from the GRW were comparable to previous estimates extrapolated from grab samples near the mouth of the Grand River over the same time period. This suggests the SWAT model adequately predicts TP export. However, the real value in SWAT is not predicting TP export over the recent past. Rather, it provides a tool for hypothesis-driven research. The GRW SWAT can be used to predict how hydrology and material transport (sediment and nutrients) may be affected by changes in the watershed, for example, increased urbanization, changes in crop type, and increases in water extraction. It can be used to explore implementation of best management practices and alternative land use scenarios, and how they may affect sustainable water management, and mitigate (or exacerbate) nutrient export to Lake Erie. Finally, it can be used to model climate change scenarios to assess effects of increasing rainfall intensity and frequency on stream flow and material export from the GRW.