GIS Spatial-Temporal Modeling of Water Systems in Greater Toronto Area, Canada

INTRODUCTION

The Oak Ridges Moraine in the GTA is a unique area for natural resources, environment and social-economic developments. The area is covered by up to 300 m glacier deposits which form a unique glacial landscape extending approximately 200 km from Rice Lake in the east to the Niagara Escarpment at Caledon Hills in the west. The moraine area varies in width from about 1 km to about 24 km. The ORM has been recognized as a major aquifer complex within Ontario that provides one of the most heavily used ground water sources in Canada. More than 6 million populations live in the GTA. 15 % of the population in urban areas, small towns and rural areas north of Toronto, use groundwater for drinking water. Groundwater has proven to be a reliable and economical resource in the GTA. Understanding the spatial-temporal distribution of the groundwater resource on a regional scale and the interactions among precipitation, ground water recharge/discharge and surface water runoff is essential for urban planning, social-economical development and water resources management in the area.

A number of government programs have been conducting in the area which have generated a vast amount of data available for hydrological research. Among them, the Geological Survey of Canada together with the Ontario Geological Survey has carried out a National Mapping (NATMAP) project (Sharpe et al., 1996). This project has created a comprehensive geodatabase for the area including a surfacial geology, the DEM, water well database and some remote sensing images. Environment and Meteorological Service of Canada has been carrying out a program monitoring the river flow and sediments data in rivers cross Canada. It has published a HYDAT database which contains the daily, monthly and instantaneous streamflow data, mean monthly water level data for over 2 900 active stations and some 5 100 discontinueous sites across Canada (Environment Canada, 1996). About 100 gauging stations are located in the study area. The streamflow data can be linked to the topographic and hydrological properties of watershed by means of GIS technology. Environment Canada has been collecting precipitation data from the ORM since 1985 which makes daily precipitation information available from about 40 weather stations located in the study area. These data have been used in the current study to characterize the spatial-temporal distributions of surface water and groundwater in relation to precipitation and topographical variances.

Watershed was derived from the DEM of 30 m resolution by means of ArcInfo (ESRI Inc., 1997) with verification with stream lines from Ontario Base Map (OBM) Sheets. A database can be created for each watershed with river gauging stations as seed point. The database contains the geometrical information of watershed (area, perimeter, slope), land use types (fraction of areas of each type of land types), stream patterns (order and length of stream segments, density of stream per watershed) etc.. Based on the information of watersheds, a number of analyses have been conducted to associate the streamflow to other properties of watershed.

(1) A fractal/multifractal model was developed to the stream patterns to define an index measuring the relative randomness of the evolution of the landscape or stream networks. It was found useful for characterizing surface stream patterns in relation to topographic properties, geological and hydrological properties (Cheng et al., 2001).

(2) An inverse mode was developed for estimating the hydraulic parameters by integrating water levels from water wells and soil types and watershed properties (Han and Cheng, 2000).

(3) Weights of evidence model was applied to model the flowing wells with respect to topographical and geological and hydrological properties. It shows that the flowing wells, springs and seepages are main sources for supplying groundwater into streamflow.

(4) An integrated spatial-frequency model was developed to identify groundwater related streamflow component (Lu, 2001). It was found that the stream flow recorded at the gauging stations in the HYDAT database can be separated into seasonal variable component which is due to rain fall and baseflow which has long wavelength variability which can be further decomposed into watershed size-dependent and independent component. The latter is most likely caused by groundwater discharge.

(5) An inverse model has been recently proposed by Ko and Cheng (2004) for associating precipitation to stream runoff. A number of statistical indexes were applied to characterize the behaviour of response of river flow to the precipitation event in relation to watershed and land use types occupying watersheds.

This paper introduces some of these models used for the above objectives; some of these models may be applicable to the similar problem in other areas of the world.

Water Systems and Their Interactions in ORM

As a landform made of sand, gravel and silt de-posited by receding glaciers, the ORM is a ridge forming the regional surface water division between water flowing south to Lake Ontario and water flowing north to Lake Scugog, Lake Simcoe, Rice Lake and Georgian Bay. As a huge water recharge, discharge and storage area, it is essential for maintaining base flow in the stream systems and water levels in kettle lakes in the GTA. In recent years, there is increasing demand for the moraine's surface water and groundwater resources for residential, commercial, industrial and recreational uses (Storm, 1997).

The Oak Ridges Moraine has a humid continental climate, warm summers and mild winters. It has a long growing season, about 180-200 days per year. Severe temperature was minimized by the influence of the Great Lakes. The north and south slopes of Oak Ridges Moraine form two different climate regions. North slope has the mean daily temperature 6.7 ℃, whereas in the south slope, the mean temperature is 7.8 ℃. It receives an average annual precipitation of 710-720 mm (from 1931 to 1960) (Sharpe et al., 2002; Brown et al., 1980; Phillips and McCulloch, 1972).

Sharpe et al. (1997) have established a geological model consisting of six principal stratigraphic elements for the area, which has served as the geological benchmark for characterizing spatial variability of river flow, ground water aquifers and ground/surface water system interactions. Among these six lithology elements, the ORM is an extensive stratified glaciofluvial-glaciolacustrine deposit 150 km long, 5 to 15 km wide and with thickness up to 150 m. It forms a prominent ridge of sand and gravel running from near Rice Lake to the Niagara Escarpment (Fig. 1). The lower contact of the ORM is an irregular channeled surface of Newmarket Till. The channels may be confined within, or have eroded through, the Newmarket Till into the lower drift below. It is generally recognized that the ORM is the main source of recharge in the region. Studies have indicated that recharge may also take place through till units adjacent to the moraine. The channels developed in the Newmarket Till contain mainly sandy sediments related to the ORM complex and some channels contain thick gravels. These channels may be hydrologically significant as high yield aquifers (Sharpe et al., 1997). Therefore, moraine deposits generally serve as recharging layers and tills as resistant layers for storing groundwater. The contacts of moraine and till deposits particularly at those local areas along the contacts intersecting with the paleochannels developed underneath the moraine deposits are the potential areas where ground water may flow into lakes, ponds and rivers to maintain the base flow in 65 rivers. The quantity and quality of storm runoff, which affects stream chemistry, also play a critical role in functioning the ecosystem in the Oak Ridges Moraine area. To test the influence of ORM on the ground and surface water interactions, several models have been implemented as will be discussed in the following sections. One model serves to characterize the spatial distribution of groundwater aquifers, especially for the location of artesian aquifers through the study of locations of flowing wells, springs and seepages. The second model is to analyze river flow data in order to associate low river flow to potential ground water supplies. The third model is to analyze the association between precipitation and runoff in the area. Understanding the associations and balances of the water systems (precipitation, surface water and groundwater) and their interactions is the main objective of this study.

Figure  1.  Surfacial geology of the ORM (Sharpe et al., 1997).

DownLoad: Full-Size Img PowerPoint

Ground and surface water distributions and their interactions

It is generally recognized that groundwater discharge is one of the main sources to the low flow of streams in the Oak Ridges Moraine area. Groundwater discharge to streams may be through springs, flowing wells, and seepages where the water level can reach the land surface. A number of large spring sites have been well documented, but economically less significant springs and most areas with seepage remain unknown. Artesian wells are wells with water levels above the land surface due to hydrostatic pressure. From about 57 000 water wells from the newly revised Ontario Ministry of Environment and Energy (MOEE) dataset compiled by the Geological Survey of Canada (Russell et al., 1998), 353 were selected with water levels above the land surface. The spatial distributions of these wells are shown in Fig. 2, superimposed on the DEM with 30 m resolution (Kenny, 1997). It can be intuitively seen from the locations of wells in Fig. 2 that most of the flowing wells occur in areas with negative relieves and in the vicinity of steep slope zones that may create hydrostatic pressure due to significant gravity difference. Therefore, spatial correlation between the locations of flowing wells and distances from the high slope zones was tested (Cheng, 2001). These flowing wells can be chosen as samples representing areas with watertable above the land surface. The locations of the flowing wells show an irregular and clustered distribution. One of the tasks to be conducted in the current research is to statistically test the significancy of spatial association between the locations of these flowing wells and other geological, topographical and hydrological features. For this purpose, various features were extracted and defined with the aid of GIS. For example, the distance from the ORM, the slope and the drift thickness were constructed using the ArcView GIS. To calculate the spatial correction between point features and line or polygon features, the weights of evidence method (Bonham-Carter, 1994) was applied. The method can be used to find optimum cut-off values to construct binary evidence. It has been demonstrated that the locations of flowing wells in the ORM area are correlated with the distances from the outline of ORM, thick drift layers and steep slope zones within the vicinity of 500-5 000 m to ORM, 500-4 000 m to thick drift layers and 1 500-2 500 m to the steep slope zones, respectively (Cheng, 2001). The combination of distances from the ORM and from steep slope zones, yielding higher correlation and cells with this combination may have posterior probabilities of possessing at least one flowing well three times higher than the prior probability of randomly selected cells from the study area. The posterior probability was further integrated to include other layers of data to study the interactions of ground and surface water systems. More detailed study by Zhang (2001) in his M. Sc. thesis was devoted to check into the spatial association of flowing wells with several other evidential layers including layers created from water well data such as elevation of sand and gravel layers recorded in water wells. Due to the interdependency of the multiple evidential layers, Zhang (2001) used a multiple logistic regression model with 9 binary evidential layers to predict the potential areas where flowing wells may occur (Fig. 3). The high potential areas with potentiometric surface of the upper artesian aquifer above the land surface are mainly located in the Laurentian Channel area, south of the boundary of the Oak Ridge Moraine, and creek valleys. The study confirmed that regional geological, topographical and hydrogeological features are important influence factors on the spatial distribution of groundwater in this area.

Figure  2.  Artesian wells in the Oak Ridges Moraine area extracted from MOEE dataset (Cheng, 2001). Background map is the shaded relief DEM with 30-meter resolution (Kenny, 1997).

DownLoad: Full-Size Img PowerPoint

Figure  3.  Posterior probability map calculated by logistic regression with 9 binary evidence layers (Zhang, 2001).

DownLoad: Full-Size Img PowerPoint

River flow analysis

Maintaining river flow is important not only for supply to human beings but also for maintenance of environmental sustainable development in the GTA. Water supply for river flow can be in either groundwater discharge or surface runoff. Groundwater discharge provides a long-term sustainable component maintaining the low flow (baseflow) that is relatively stable and independent of seasonal precipitation. The surface runoff and overland flow usually correspond to precipitation and are influenced by topographical features such as soil types, vegetation cover, and stream network and drainage systems. Studies based on stream discharge records obtained from about 70 gauging stations in HYDAT CD-ROM (1996) for more than 30 years, show that baseflow, the part of stream discharge from groundwater seeping into the stream, is lower in western watersheds and higher eastern watersheds. Differences in baseflow per unit of watershed area suggest spatial differences in groundwater flow across topographic divides (Hinton, 1996). Significant spatial differences in baseflow are also noted within individual basins (Sharpe et al., 1998). In order to conduct thorough spatial analysis of river flow in the ORM area, more recent river flow records from the HYDAT database (Environment Canada, 1996) were extracted and 98 gauging stations were selected from the database.

To decompose the river flow into separate components reflecting baseflow and surface flow, an integrated statistics and spatial analysis was developed. The general model proposed by Cheng (2001) is illustrated as a flow chart in Fig. 4. The data processing consists of two main steps: the first step involves frequency analysis and the second step spatial analysis. In the first step, river flow is separated into low flow and residual flow on the basis of flow frequency. The low flow is the basic component reflecting "baseflow" variability. In the second step, the low flow was further separated into watershed-dependent and independent low flow. The watershed-independent low flow is most likely caused by ground water discharge potentially those paloechannels developed in Newmarket Till beneath the ORM.Lu (2001) hasimplemented the model in her M.Sc. thesis. The residual flow is the component closely corresponding to seasonal precipitation changes. Lu (2001) calculated the low flow using two methods, one of which is to calculate the average of 30-day minimum flow values over many years. The low flow calculated in this manner may still contain flow component due to surface runoff, which can be confirmed by plotting the low flow against drainage basins and their attributes such as size of drainage basins, perimeter of drainage basins and total length of stream segments in drainage basins. If we assume that groundwater discharge to the rivers is limited to the small areas where aquifers are exposed at or near surface, then these types of groundwater discharge will be relatively independent of the size of drainage basins. Therefore, in the second step of the process, the low flow was plotted against drainage area size related attributes. A multiple regression was then applied to the low flow values and drainage area size related attributes. Significant correlation has been demonstrated to exist between the low flow and the drainage-related attributes. The residual values of the regression, independent of all the drainage-related attributes, can be treated as the low flow component, mainly caused by groundwater supply (Fig. 5). The result of Lu's work (2001) has demonstrated that a significant component of low flow in the area may be due to groundwater discharge. More detailed implementation and discussion of spatial variability and association with geological, hydrological and topographical features can be found in Lu (2001).

Figure  4.  Flow chart showing the two-step processing for decomposition of river flow. The first step involves frequency analysis and the second step deals with drainage-based spatial analysis. The residuals obtained from a regression analysis can be treated as flow component related to ground water discharge.

DownLoad: Full-Size Img PowerPoint

Figure  5.  Drainage basins colored according to river flow values. (a) log-transformed flow values; (b). residuals obtained from regression. The residuals correspond to the decomposed low flow values representing ground water discharges. Yellow dots are flowing wells (from Lu (2001)).

DownLoad: Full-Size Img PowerPoint

Precipitation and runoff

The previous discussions have shown that groundwater and surface water system interactions are significant and that water supplies from surface and ground to river flow are essential for maintaining river systems in the area. In order to comprehend the balance of precipitation, groundwater and surface water systems, the connection between precipitation and runoff has to be further examined. Ko and Cheng (2004) studied the spatial relationship between precipitation events and storm runoff using several statistical analyses including cross correlation coefficient, auto correlation coefficient and fractal R/S analysis (Hurst, 1951). The cross correlation analysis is to determine the delay responds of storm runoff to precipitation events by measuring the optimum cross correlation coefficient between the two time series (storm runoff and precipitation). Auto-correlation and R/S analysis are applied to characterize the short-term and long-term temporal dependence of the stream flow records. Visualizing the spatial pattern of the calculated statistics (only correlation coefficients shown in Fig. 6) in the drainage basin context provides a better view of the spatial distribution of storm runoff generation in relation to the physical properties of different watersheds in the ORM. River flow and precipitation of multiple year records (1985-1999) have been used for the statistical study.

Figure  6.  Shifted times (in days) reach the highest correlation between river flow (with low flow component removed seen in Fig. 5) and precipitation calculated by cross-correlation coefficient analysis for selected river gauging stations. The colors represent the average results for each drainage basin with its gauging station as seed point.

DownLoad: Full-Size Img PowerPoint

The study has demonstrated that the respond time of the stream flow to the precipitation events measured by cross correlation coefficient statistics ranges from 1 to 7 days varying from drainage basin to drainage basin. The optimum auto-correlation range and the Hurst exponent (H) measured from the river flow records react to the storm runoff patterns from short-term and long-term point of views, respectively. To explain the distributions of these characteristic values, the physical properties of the watersheds were associated to the storm runoff processes, generally produced either by Hortonian overland flow (HOF), saturation overland flow (SOF) or fast subsurface flow (SSF).

Drainage characteristics quantified include soil and land types, topographic features and drainage basin geometry. It has been found that drainage basins containing large urban area tend to have short response time whereas the watersheds occupied by high proportion of the vegetated area show long response time between precipitation event and high peak river flow. In other word, the streams in such drainage systems would need longer time to respond to the precipitation event. In the latter case the stream flow shows more persistence behavior. In checking the association between the statistics (cross correlation, auto correlation and R/S index) and topographic property of drainage basin, it has been demonstrated that the general slope of drainage basins has significant influence on the behavior of the river system and as well as the responding time to precipitation events, for example, in the north side of ORM, the general slope of watersheds is relatively low in comparison with the drainage basins in the south side and, as shown in Fig. 6, the drainage basins in the north side generally take longer respond time to precipitation event and the flows are more persistent. Drainage basin geometry was found having significant influence on the river flow. For example, the ratio of total channel length and basin area, circularity and compactness were calculated to each of the drainage basins and they were plotted against the statistics (cross correlation, auto correlation and R/S index). It has been shown that long narrow basin would expect a lower and gentler peak respond with shorter time-lag to the precipitation event whereas a circular basin would generally produce a higher and sharper peak in respond to the precipitation event with longer time-lag. More results about analysis of precipitation and runoff in the study area can be found in Ko and Cheng (2004).

CONCLUSIONS

Various spatial analysis models proposed and applied to the datasets from the ORM of Greater Toronto area, Canada, have demonstrated that spatial and temporal aspects of the water systems are essential for quantification of the interactions among ground water, surface water systems and precipitation events on the basis of point observations and measurements made at river gauging stations, water wells and weather network stations. Integration of GIS and statistics can provide powerful tools for solving complex problems such as separation of baseflow and surface flow which might not be possible without taking into account the spatial associations of flow data and other geological, hydrological and topographic features. The quantitative results have demonstrated that the ORM, as a unique glacial deposit, plays an important role in maintaining and balancing the water distribution in the area. ORM provides significant aquifers for the GTA and the channels occur at the contact between ORM and Newmarket Till have potential to become a significant ground aquifer for the area. The environmental sensitivity of water systems of the area should be carefully considered when urban development is planned for the ORM.