Thomas W Charnock, John Elgy, and Peter D Hedges,

Application of GIS Linked Environment Models over a Large Area


Abstract

Linking environmental models to a GIS makes it practical to use the models over a large region whilst maintaining the scale of detail that the model was designed for. This paper discusses the use of GIS and linked numerical models to evaluate the poten tial effects of groundwater abstraction upon crop production. Existing models of the groundwater and soil moisture system have been selected and linked with the GIS, crop models are currently being evaluated.

The most appropriate models for the soil moisture component of the system are 1-D soil moisture models that can simulate the response of the soil moisture profile to the falling watertable. To build up a picture of the soil moisture deficit as it devel ops through the growing season, the model is run repetitively across the region of interest.

Such a modelling arrangement leads to large data requirements and issues of calibration, particularly when the soil moisture model is linked to models describing the other environmental subsystems. During this project experience has shown that one must be pragmatic about addressing these issues.

INTRODUCTION

Environmental process level models have a scale range over which it is appropriate to apply them and there are pitfalls in using a model at a scale or resolution greater or smaller than that for which it was designed. GIS has the data handling capabili ty that makes it practical to perform modelling operations over a large area whilst still maintaining the scale of detail for which the model was designed. Furthermore GIS provides a medium for linking different environment models together.

This paper concerns the use of GIS and numerical models for Environmental Impact Assessment. The purpose of the methodology is to answer one simple question; what is the best estimate of the spatial variation of crop yield reduction due to water extrac tion from the Shropshire Groundwater Scheme? Clearly implicit is that cost of answering the question is to be less than the expected value of the crop loss itself. This assumption drives the methodology. Thus we cannot derive unnecessarily complex models from scratch nor sample ground truth at very fine resolutions. To this end a modelling system is being developed by linking several numerical models via a GIS. Whilst the strategy of overall system development and the problems encountered are discussed in Charnock et al (1996), this paper looks principally at the soil moisture and groundwater components of the system, in particular the data requirements.

Of prime concern is the movement of moisture within the soil profile in response to watertable drawdown. The models which best represent this are 1-D vertical models of the unsaturated zone, which adequately represent soil moisture conditions over a fe w tens of metres. As the area of potential effects is over a hundred square kilometres, it is therefore necessary to run the 1-D model at points across the landscape to simulate the spatial development of soil moisture deficit through a growing season.

Although this approach avoids the pitfall of using a model at a scale for which it was not designed other issues are raised. These include data requirements and model calibration, as well as choosing an appropriate regular or irregular spacing at which to run the soil model.

BACKGROUND

The Shropshire Ground Water Scheme is being developed in order to augment flow in the river Severn during periods of drought. The area of interest is the Tern river basin, Figure 1.

a

b A field of potatoes c potato plants

d Tern River with typical pasture field in the riparian zone
Figure 1,(a) location map, (b) field of potatoes next to Heath House Bore Hole, (c) potatoes next to Heath House Borehole with rain-gauge and soil moisture monitoring tube, (d) the River Tern with typical pastural riparion zone.

Here the watertable is close to the surface and, though there is some clay drift, for the most part the soil is in hydraulic continuity with the aquifer. Under these conditions it has been shown that pumping can potentially reduce the soil moisture avail able for crops (Hedges and Walley, 1983).

soil moisture theory,(a) idealised tension conductivity curve, (b) critical height, (c) idealised soil moisture profile at field capacity
Figure 2,(a) the relationship of conductivity to tension, (b) the idealised tension profile at field capacity, (c) idealised moisture content profile.

When a soil with excess water drains it approaches a state at which tension is proportional to height above the watertable. However, a critical tension is reached at which the unsaturated conductivity drops rapidly and further drainage is negligible, Figure 2a. The soil moisture profile achieved at this tension is the field capacity profile. Hedges and Walley introduced the concept of critical height, defined as the height above the water table below which the position of the water table controls the soil moisture content (Figure 2b and c). If the watertable drops then any plants whose roots previously penetrated this zone will experience a loss of moisture. As illustrated in Figure 3.

loss of soil moisture available to crops due to groundwater drawdown
Figure 3, the loss of soil moisture due to groundwater drawdown.

Crops are vulnerable wherever the soil is in hydraulic continuity with the aquifer and the depth to the groundwater table is smaller than the critical height plus the crop rooting depth.

D<ho+dr (1)
D is the depth to watertable.
ho is the critical height.
dr is the crop rooting depth.

Hedges (1989) used this simple rule within a GIS in order to delineate vulnerable areas, and used a map of predicted drawdown to identify those areas that would actually have been affected by that particular event. However, there are drawbacks to this approach; the groundwater table can be expected to drop naturally in a dry year; the soil moisture profile will be naturally depleted by crop root abstraction and, the soil moisture profile will have a complex response to watertable drawdown depending on the physical properties of the soil profile.

To improve upon the basic approach of Hedges, numerical models are being applied to simulate the groundwater movement and the response of the soil moisture and crop systems more realistically.

SYSTEM DESIGN

The overall system design and the issues raised by this approach to development are discussed in Charnock et al (1996), and will only be briefly reviewed here. Because of constraints of time and money the strategy of system development has been one of linking existing components together with a network of communicating programs. Minimal changes are made to existing code. It is doubtful whether any other approach would be cost effective.

GRASS is the GIS used (see Shapiro et al, 1993). This is a raster based system with some vector capability. The groundwater system is modelled by Modflow (McDonald and Harbaugh, 1988) which is a three dimensional finite difference formulation. T he soil moisture system is modelled with SWMS_2D (Simunek et al, 1994), which simulates unsaturated flow and solute movement in 1 or 2 dimensions. For this project only 1 dimensional moisture flow is required. A suitable crop model is still being s ought.

link between SWMS_2D and GRASS
Figure 4, the link between the GRASS raster database and SWMS_2D input files.

Figure 4 illustrates the link between the SWMS_2D and the raster database of GRASS. Rasters are used to represent the spatial variation of input parameters across the landscape. The linking program samples these raster at a specified point, for which it g enerates a finite element grid of the soil column and writes the input files for SWMS_2D. A picture of the spatial variation of soil moisture deficit development is built up by calling the model successively at different points across the landscape; the a ttribute of interest is extracted from the SWMS_2D output and written into the appropriate cell of an output raster. In order to facilitate this, the link is embedded within a script which controls the calling of the model link and the build up of the out put rasters

Once any model has been linked into the GIS, channels of communication are effectively opened between the model and all the other operations and functions of the GIS, including other models. The models can thus be used in a similar manner to any other function of the GIS.

The links to the models have been designed to be implemented quickly and easily, and so they contain none of the logic of the problem to which they are going to be applied. However, once linked into the GIS they can be embedded within an appropriate ma cro language (in GRASS this is a UNIX shell script), and it is here that the logic of the problem, in terms of the subsystems to be modelled and the assumptions made about the interaction between them, are specified.

the modelling system
Figure 5, the modelling system viewed as component layers and as a hierarchy of control.

Figure 5 illustrates the logical arrangement of the models, which can be viewed as a layered approach. The groundwater model defines the water table as it varies across the landscape, SWMS_2D is called at an appropriate regular or irregular grid spaci ng, and the crop model will, when it is included, act either like the SWMS_2D model at points across the landscape or be applied to larger spatial units. Figure 5 also illustrates the hierarchy of control. In the simplest situation the models are run sequ entially but more complex arrangements can simulate interaction between the environmental subsystems by repetitive calls to the models (Charnock et al, 1996).

DATA REQUIREMENTS

As with system development, there are constraints of cost and time on the collection of data. The approach adopted has been to use secondary and surrogate sources of data wherever possible. These include geological and soil maps, standard meteorological d ata, remotely sensed data, data generated during the planning, development and operation of the Shropshire Groundwater Scheme and other documented studies.

Data requirements for the whole system include; soil properties, watertable position, potential evaporation and transpiration, precipitation, crop rooting depth and uptake parameters, aquifer transmissivity, storage and vertical conductivity, drift dis tribution and river leakage parameters.

The large scale variation of meteorological processes is an ongoing subject of research. However the system under development is intended for application to scenarios which are dry with very little precipitation, so it is reasonable to employ data from local meteorological stations utilising simple procedures such as voronoi polygons or simple interpolation.

For this project the spatial distribution of crops is required, and this will change year by year. Remotely sensed data is the most convenient way to provide timely crop information. The literature provides typical values for crop parameters, such as r ooting depth, and thus a raster of crop distribution can be classified into a distribution of rooting depths.

Aquifer properties, such as transmissivity and storage coefficient, are derived from the many pumping tests of abstraction boreholes undertaken during the history of the scheme. Watertable levels have proved to be problematical. There is an extensive n etwork of monitoring boreholes and tube wells but few have been drilled specifically for observation purposes. Thus the network is non-random, and the application of interpolation routines leads to an over estimation of watertable heights in areas where t he network is sparse.

To illustrate the problems of deriving spatially varying model input parameters from secondary and surrogate sources of data, the example of soil data will be discussed in detail.

Soil data

Soil data is needed principally for the soil moisture component. SWMS_2D uses a modification of van Genuchtens (1980) equation to describe the hydraulic properties of the soil. The moisture-suction and conductivity-suction curves are defined using inpu t parameters which include; saturated hydraulic conductivity, residual and saturated water content, as well as empirical curve fitting parameters. Evans (1990) used regression analysis to derive relationships that enable many of these parameters to be cal culated from soil properties that are easily measurable or obtainable from secondary sources such as maps.

extent of available soil data
Figure 6, the availability of soil property data.

Figure 6, shows the extent of available soil data. The whole region is covered by a small scale (1:250,000) soil map, and part of the region is also covered by a larger scale (1:50000) map. The soil series is the basic mapping unit in the UK, with clas sification based on the soil profile from the surface to a depth of about 1.50m (see Soil Survey of England and Wales, 1984). The characteristic profile of a soil series is based on the particle size distribution, mineralogy, calcium carbonate content, or ganic content and the pedogenic (soil forming) processes of the soil horizons. Thus, using relationships such as those devised by Evans (1990) general purpose soil maps can be translated to maps of the spatial variation of model input parameters.

However, maps, particularly small scale maps, can be expected to have considerable inaccuracy. The types and causes of error and the importance of scale in maps generally have been well discussed elsewhere, see for example Burrough (1986). Soil maps al so have their own distinct set of limitations. The concept of sharp soil boundaries is an understandable fiction; boundaries are indistinct and within a soil series there is variation in the soil profile.

There are other soil data in addition to the two available maps. During the history of the Shropshire Groundwater Scheme several sites have been extensively studied, and there was also a very intensive survey over a small area in the Tern river valley, (see Figure 6). The point data takes a variety of forms: a detailed description of the soil profile obtained from a trial pit or augered sample; a description of soil hydraulic properties determined by direct measurement in the laboratory or field; or a record of the soil moisture regime over time. For the intensively studied area, the augered soil profile, obtained on an approximate 120m grid enabled the profiles to be classified into soil series, but with more detailed subdivisions being defined.

Crop yield varies non-linearly with respect to soil moisture parameters. Within a single soil series important parameters, such as saturated hydraulic conductivity, vary significantly. If we assume representative figures for parameters within an homoge neous unit, for example the arithmetic mean, then misleading overall yields could be produced. We must represent the naturally occurring variability within the model system. Burrough (1986) described geostatistical methods for predicting soil parameters a nd made the point that;

" if interpolation techniques are to be of any value then they should be used if possible in conjunction with, and not instead of, conventional landscape mapping methods."
Rogowski and Wolf (1994) used a GIS to facilitate the process of combining spatially interpolated values with delineated soil map units, and produced a map;
"that preserves the map unit boundaries but incorporates spatial variability of attributes within units".
Kiefer (1993) uses a stochastic representation of soil heterogeneity. Other authors have investigated the stochastic properties of soil heterogeneity, ( e.g. Binley et al 1989). Rather than review their work here, one very simple stoc hastic rule that has been used to generate the soil properties layers is included as an illustration. It is based on the following assumptions:

  1. There is no trend within a soil series.
  2. There is a weak spatial auto correlation over the grid size that has the form:
    p(d)=k^d (2)
  3. Individual parameters are correlated with each other in a simple linear fashion.
  4. Our measured points are augered samples not absolute values over the whole grid. We use these to estimate the spatial statistics and then generate representative values for the whole cell containing the pit.
  5. Different realisations of the same stochastic set will give different values of yield loss and this will enable reliability figures to be applied to the loss estimates.

The standardised value for the most critical or least variable parameter can be calculated for each grid square:

v(i,j)=a[v(i-1,j) + v(i,j-1)] + bv(i-1,j-1) + cn (3)
v(i,j) is the standardised value for the property for cell row i column j.
a,b and c are parameters for each soil series and soil property to be established by method of moments using the available large scale data.
n is a random(0,1) number.

The standardised value is multiplied by the standard deviation and added to the mean of each soil series. Subsequent parameters are calculated using a similar method or generated from the single core distribution, whereby cross-correlation is accounte d for, but the very weak spatial auto-correlation is not specifically maintained

a soil series delineations

b raster generated using stochastic properties of soil series
Figure 7,(a)an example of delineated soil series, (the real map cannot be shown because of copyright),(b)a realisation of the stochastic variation of saturated moisture content of the top horizon, this will be used as one of the raster layers in Figure 4.

This is of course a very simple formulation but it has the advantage that it is easy to implement using a map algebra facility, and serves well as a first step in handling soil series heterogeneity. To go any further may be unwarranted bearing in mind the accuracy of the model and other input parameters. More sophisticated approaches include those discussed by Gessier et al(1996, this publication)

The concept of a soil series is of course an attempt to handle soil heterogeneity but in a general purpose way. The question that remains to be answered is whether the use of soil series leads to sufficient specification of spatial variability of soil hydraulic properties to be of use for the current application. If the answer is no, then the use of some method, as described above, is called for.

Similarly, do the other sources of surrogate and secondary data describe the spatial variation adequately. The most obvious feature of the generated rasters is that they will have a spatially varying quality, depending on whether the soil series is tak en from the large or the small scale maps, which has to be taken into account when decisions are based on model outputs. The approach followed here is to go as far as is reasonable to capture spatial variability. Reasonable is defined in terms of the cos t in time and money of improving the quality of data with respect to the improvement in the overall modelling system and the estimates being made. Where accuracy is limited by other parts of the modelling system such as; the accuracy of other data within the model component, the accuracy of the model; the accuracy of the model components or the way the components are linked, then improving the parameters estimates for one particular subsystem are not justified.

CALIBRATION

An important stage in any numerical modelling is that of calibration (see Sorooshian and Gupta, 1995 for a good review). However the modelling system outlined above is not necessarily amenable to this process. Calibration can be applied to subsystems o f models, to the individual models or to the whole modelling system of GIS linked models.

Realistically, reaching an optimal set of inputs is probably impossible. However, obtaining a satisfactory set of parameters for a particular purpose and under a constrained set of conditions (in this case drying conditions with a falling watertable) m ay be possible, once it is accepted that the model will need to be recalibrated if used for some other purpose. The model system components can be calibrated individually, or the modelling system can be calibrated as a whole - again the soil system is des cribed by way of illustration

Calibrating the soil moisture component

Calibrating the soil moisture component is the most problematical aspect of the modelling system. Each point on the landscape is being modelled independently, i.e. lateral flow is not considered. However the input parameters are generated as spatially varying rasters. This a more realistic way to treat soil input parameters and, because the model is to be run at thousands of points across the landscape, specifying the input parameters individually is impractical.

The Shropshire Groundwater Scheme has a few sites at which the soil moisture regime has been monitored for a number of years. These form a good record against which to calibrate the model system. Here it is practical to vary the parameters individuall y. However once you have a set of parameters at a point the question is how to apply them over a large region. Once again we have the choice of treating the point data as representative of an homogeneous unit or attempting to take into account some of th e spatial variability of soils.

Hopmans and Guttierez-rave (1988), discussed the problem of calibrating the root uptake function of the soil moisture model SWATRE over a large area with a heterogeneous soil properties (NB, the root uptake function in SWMS_2D is based on the SWATRE fu nction). They found that using a monte-carlo analysis with a trial and error method to optimise parameters gave better transpiration results than with single site calibration. However, they also found a considerable error between model output and measure d transpiration in particularly dry years.

Fitting the model to monitored sites or point data that provide a very detailed description of the soil moisture profile through time but are spatially very sparse is one approach. An alternative is to use some measure with good spatial coverage but le ss detail in the profile. The obvious candidate for such a measure is remote sensing, particularly microwave images. The use of microwave remote sensing for soil moisture monitoring is a current research topic, see for example Rogows ki and Engman (1996, this publication). However this technique only interrogates the top layers which are primarily influenced by atmospheric processes. For this project it is the effect of the falling watertable that is important, and hence the use o f this approach should be treated with caution. Other possible candidates for an easily obtainable spatial indicator of moisture status include remotely sensed crop properties such as stress or temperature.

Calibrating the whole system

The object of the project is to estimate the environmental impact of groundwater drawdown on crop production. Thus crop yield is an appropriate measure against which the whole system can be calibrated. Though it is possible to gather spatially varying yield data, historically no data has been collected in this way within the Shropshire Groundwater Scheme. Some farm records are available, which give crop yields lumped either to the farm unit or in some cases to the field unit. Given the error and uncert ainty in the inputs, utilising this data for calibration is the most robust approach, as using a larger spatial unit acts to cancel out, but not remove, some of the smaller scale spatial error. There is also the potential for combining this data with remo tely sensed data to give a better estimate of the spatial variability of yield, through a simple overlay operation in the GIS.

Again the approach to calibration must be pragmatic and cost effective. Thus the method adopted is to calibrate both at the individual component level and at the modelling system level using all the available data either during the calibration or in th e verification, but not calibrating the system beyond what is required. Consequently a trial and error approach has been adopted.

SUMMARY

Linking large scale models to GIS with its data handling capabilities and its spatial manipulation functions makes it practical to apply the models over a large region.

For the project described in this paper experience has shown that the model system accuracy and hence its applicability to decision making are limited by the data availability and quality, together with the capability to calibrate the system rather tha n by the sophistication of the model components themselves. This is likely to be the case for any similar project unless a significant amount of time and money are invested in data collection.

ACKNOWLEDGMENTS

The authors acknowledge the support of the Engineering and Physical Sciences Research Council and The Royal Academy of Engineering, and the co-operation of the National Rivers Authority Severn Trent Region.

REFERENCES

Binley, A.M., Beven, K.J., and Elgy J., (1989) A physically-based model of heterogeneous hillslopes. II. Effective hydraulic conductivities. Water Resources Research 22: pp 1531-1536.

Burrough, P.A., (1986) Principles of Geographical Information Systems for land resources assessment. Oxford Science Publications: Monographs on Soil and Resources Survey no 12.

Charnock, T.W., Hedges, P.D., and Elgy, J. (1996) Linking multiple process level models with GIS. Application of Geographic Information Systems in Hydrology and Water Resources Management, Proceedings of HydroGIS'96, Ed. K. Kovar and H.P. Nachtnebel, IAHS Press, Wallingford, UK. pp 29-36.

Evans, D.A., (1990) Development of a conceptual model of the soil-moisture-plant sub-system of the hydrological cycle PhD Thesis, The University of Aston in Birmingham.

Gessler, P., McKenzie, N., and Hutchinson, M., (1996) Progress in Soil-landscape Modelling and Spatial Prediction of Soil Attributes for Environmental Model s. Proceedings of the 3rd International Conference/Workshop on Integrating Geographic Information Systems and Environmental Modelling, Sante Fe, New Mexico, January 21-25th.

Hedges, P.D., and Walley, W.J. (1983) A study of soil moisture losses due to the drawdown of a shallow water table. Scientific Procedures Applied to Planning Design and Management of Water Resources Systems Proceedings of the Hamburg Symposium A ugust 1983: IAHS publication no. 147: pp 489-499.

Hedges, P.D. (1989) The impact on agriculture of the drawdown of shallow watertables. PhD Thesis, The University of Aston in Birmingham.

Hopmans, J.W., and Guttierez-Rave, E. (1988) Calibration of a root water uptake model in spatially variable soils. Journal of Hydrology. 103: pp53-65.

Kiefer, E.M. (1993) A conceptual-stochastic model of unsaturated flow in heterogeneous soils. Journal of Hydrology. 143: pp3-18.

McDonald, M.G., and Harbaugh, A.L. (1988) A modular three-dimensional finite difference ground-water flow model. Techniques of Water Resources Investigations of the United States Geological Survey. Book 6, Chapter A1.

Rogowski, A.S., and Engman, E.T., (1996) Using a Decision Support System in a GIS Framework to Model Spatial Distribution of Soil Water. Proceedings of the 3 rd International Conference/Workshop on Integrating Geographic Information Systems and Environmental Modelling, Sante Fe, New Mexico, January 21-25th.

Simunek, J., Vogel, T., and van Genuchten, M.Th. (1994) The SWMS_2D code for simulating water flow and solute transport in two-dimensional variably saturated media, Version 1.2. US Salinity Laboratory, Agricultural Research Service, U.S.D.A., Ri verside California: Research Report No. 132.

Soil Survey of England and Wales (1984) Soils and their use in Midland and Western England. Soil Survey of England and Wales: Bulletin no. 12.

Sorooshian, S. and Gupta, V.K. (1995) Model Calibration. Computer models of watershed hydrology: Singh, V.P. (ed), chapter 2, pp23-63.

van Genuchten, M.Th., (1976) A closed-form equation for predicting the hydraulic conductivity of saturated soils. Soil Science Society of America Journal 44: pp892-898.


Thomas W. Charnock
Research Student, Aston University
Department of Civil Engineering
Aston Triangle
Birmingham, B4 7ET, UK
Telephone: UK 0121 359 3611
FAX: UK 0121 333 3389
EMAIL: charnotw@sun.aston.ac.uk

John Elgy
Lecturer, Aston University
Department of Civil Engineering
Aston Triangle
Birmingham, B4 7ET, UK
Telephone: UK 0121 359 3611
FAX: UK 0121 333 3389
EMAIL: J.ELGY@aston.ac.uk

Peter D. Hedges
Lecturer, Aston University
Department of Civil Engineering
Aston Triangle
Birmingham, B4 7ET, UK
Telephone: UK 0121 359 3611
FAX: UK 0121 333 3389
EMAIL: P.D.HEDGESaston.ac.uk