The main problems in watershed management are deforestation and the expansion of agricultural land, which together cause soil erosion and nutrient depletion. Satellite Remote Sensing (SRS) and Geographic Information Systems (GIS) are useful tools for watershed management and monitoring. SRS data have been found to provide timely and reliable information for monitoring forest resources and agricultural crops, as long as they are supported by ground data. Obtaining the same information by field surveys would take a long time and considerable human resources. In addition, SRS data allows us to see changes in the spatial distribution of various land resources, which is otherwise difficult to keep updated on a national scale. However, the digital classification of SRS data on slopelands and mixed cropping is generally not very accurate. The Digital Elevation Model (DEM) has been found to improve digital classification in land use and cover maps. DEM is created from contour lines, using GIS techniques. Land use and land cover are major factors in watershed analysis, particularly in relation to soil erosion and the nutrient balance. Other factors such as rainfall intensity (erosivity), soil erodibility and slope gradient and length can also be identified and analyzed by GIS techniques.
Satellite Remote Sensing (SRS) is the gathering and recording of information by a sensor on board a satellite orbiting in space without any physical contact with the object or area being investigated. SRS is used to collect information about the earth's features, such as its geology, vegetation, soil, atmosphere, water, ice surfaces and land-use.
Geographic Information System (GIS) is a system which deals with information related to the spatial distribution of features on the earth's surface. The system is designed to efficiently capture, store, update, manipulate, analyze and display all forms of geographically referenced information.
The Digital Elevation Model (DEM) is defined as "any digital representation of the continuous variation of relief over space" (Burrough 1986). It can be created from stereopairs derived from SRS data or aerial photographs, or can be generated from digital terrain elevation data. As well as its many other uses, DEM derived products can be readily combined with SRS data. For example, Guneriussen and Johnsen (1996) used DEM to calibrate and geocode satellite data for monitoring snow. DEM is also very useful for discriminating land use and ground cover classes during the digital processing of SRS data.
The interfacing of GIS, DEM and SRS provides a new and exciting capability to analyze the dynamics of land-use change. Nellis et al. (1990) used GIS to classify and assess changes in vegetation in a landscape which had been altered by fire. This analysis showed that post-classification comparison within the context of GIS is a valuable procedure in landscape assessment, using sets of data from diverse sensors.
Coleman (1992) proposed a three-dimensional modeling of an image-based GIS to assist in land-use planning and management operations. The model was based on several interrelated databases, including image databases generated from Landsat data, digital line graph files, and a spatial database which handled the attribute data relating to these spatial coverages. Image-based GIS, viewed in a three-dimensional display, provided the user with a more realistic view of landscape conditions. The results showed that the image-based GIS is a first step towards developing a more comprehensive management system for agricultural land.
The unique feature of SRS compared to other tools is that it can be used to collect data for baseline inventory and future monitoring purposes. SRS can also be integrated with other tools for further analysis. The present SRS technologies can be categorized into optical (visible), near infrared (NIR) and microwave remote sensing. Optical sensing is more widely used than microwave sensing.
There are a number of SRS systems, either now operating or soon to come into operation. They include Landsat, MOSS, SPOT, NOAA, ERS, RADARSAT, JERS, MODIS and ADEOS. Each of the sensors they use has different but similar spectral and radiometric characteristics. The potential uses of these sensors depends on variations in the land feature under study.
SRS is superior in many respects, including its cost effectiveness in covering a large area. Other advantages are its geo-referencing, and the fact that its monitoring is close to real time. The vegetation in certain areas has a strong correlation with the productivity of the land, thus illustrating watershed conditions. As vegetation has a strong spectral response in the NIR band, vegetation indices are usually developed, to show the components of vegetation cover and their health.
In most tropical climates, a major limitation in the use of V/NIR data is the high probability of cloud cover. Cloud-free data are available only during the dry season. Such data contain little information on the spectral response of vegetation.
The Sakae Krang watershed is situated in the central plain of Thailand (Fig. 1), and covers 503 thousand hectares. The climate is tropical Savannah, with an average mean temperature of 28.8°C and an average annual rainfall of 1,225 mm. The average number of rainy days each year is 99, while the average relative humidity is 73.5%.
The most common landforms in the watershed are alluvial terraces and undulating fans, which together comprise about 45% of the total watershed area. High and low hills with eroded surfaces occupy about 41% of the area, while the rest is floodplains, peneplains and levees.
The eastern part of the watershed is a floodplain. The middle part is an alluvial terrace with undulating topography. Small residual limestone hills run from north to south, separating the lowlands from the uplands. The western part of the watershed is mountainous and forested.
Forest covers about one-third of the watershed, about 14% of which is classified as head watershed class-1 A (i.e. the forest is protected, and logging is prohibited). Due to the wide range of elevation (20 to 1,684 m above sea level), three types of forest are found in the area. Dry dipterocarp forest is found at lower altitudes, mixed deciduous in middle altitudes, and dry evergreen forest above 800 m. Most of the dry evergreen forest is protected, since it lies inside the Huai Kha Khaeng wildlife sanctuary (a World Heritage site).
Two-thirds of the area is being used for agriculture. Rainfed rice is the major crop in the lowlands, although upland rice is also cultivated in some areas. In the uplands, sugarcane is the main crop. Other upland crops such as corn, beans, oilseed and cassava are also grown, but only on a small scale. By and large, monocropping is practiced, but multiple cropping is found in about 4% of the watershed area. In some irrigated areas, a second crop of rice is grown during the dry season.
Apart from three soil series which are classified as having a "Udic" (humid), moisture regime, most of the 52 soil series in the flood plain and levee basin have an "Aquic" (wet) moisture regime. The soil texture of the Udic soils is mostly coarse loam, while the Aquic soils vary from fine clay to coarse loam. Soils in the alluvial terrace belong to the "Ustic" (wet and dry) moisture regime. The texture of these soils varies from clay and skeletal clay to loam and silt.
Remote sensing by satellite is an important tool for monitoring and managing watershed resources. In many cases, digital image classification, including data manipulation preprocessing and enhancement, is used to process the satellite data (Cracknell and Hayes 1993). In classification, spectral values are grouped together to form classes to which a theme is assigned to produce a thematic map, such as vegetation cover on a land-use map. One of the main problems when generating land cover maps from digital images is the confusion of spectral responses from different features. Sometimes two or more different features with similar spectral behavior are grouped into the same class, which leads to errors in the final map. The accuracy of the map depends on the spatial and spectral resolution and the seasonal variability in vegetation cover types and soil moisture conditions (Campbell 1987, Congalton 1988). Hence, improving the accuracy of classification has always been a major concern in digital image processing.
Attemps have been made to improve the accuracy of image classification. One approach is the use of multi-temporal imagery to distinguish classes (Conese and Maseli 1991, Barbosa et al. 1996). Another is the piecewise linear classifier with simple post-processing (Kai-Yi and Mausel 1993). The integration of GIS with ancillary information has also been tested, to improve image classification. For example, Gastellu-Etchegorry et al. (1993) and Ortiz et al. (1997) integrated GIS with information about soils, topography and bioclimate. Similarly, Palacio-Prieto and Luna-González (1996) employed GIS rules with ancillary data on terrain mapping units and elevation data.
When we were working with optical data such as Landsat-TM from the Sakae Krang watershed and processing digital images, we found we had a major problem in discriminating the different types of land use or ground cover. It was difficult to acquire cloud-free data during the rainy season (December - April), when all types of vegetation were growing. Cloud-free data could only be acquired after the rainy season was over, at which time the data contained very little information about crops.
Landsat-TM data acquired in February 1995 were employed in this study. Other materials used were topographic maps (1:50,000), a land use map dating from 1993 (1:250,000), a soil map (1:100,000), and field survey information. The topographic map was used for elevation data to generate the DEM, so as to perform a geometric correction of the Landsat data. The soil map and land use map were used as ancillary information about soil types in the area, and associated land use in the past. Field visits were conducted for training and verification.
Contour lines at 20 m intervals in hilly areas, and at 10 m intervals in the plains, were digitized from the topographic map, using ARC/INFO software. These data were imported for interpolation by ERDAS software.
ERDAS employs a 3 x 3 weighted moving average window to process group of locations as blocks within a specified search radius. The search radius specifies the distance around each pixel to be interpolated to search for terrain data points, the elevation values of which are used to interpolate the elevation of the former pixel (ERDAS 1991). The weighting function determines how the surface values will be interpolated from the original elevation data points.
The Landsat-TM data were radiometrically calibrated by comparing the gray-level histograms of each image, and the spectral value of sample pixels for identified features, so that the same feature in two different scenes had similar spectral values.
A geometric correction of the image was carried out, registering the image to the 1:50,000 scale topographic maps by selecting 45 Ground Control Points (GCPs). The Root Mean Square error accepted was less than 1 pixel (30 m) at the first order and the nearest neighborhood transformation. The two corrected scenes were stitched together, and only the study area was clipped with a vector boundary layer.
Apart from band 6, the standard image processing techniques employed in this study, were NDV1, 6, band ratioing and PCA. Twelve Classification Schemes (CS 1 - 12) were decided prior to the classification (Table 1). As there could be a large number of band combinations, the selected band combination in each CS was decided according to the relative use of the individual band, and the band combinations which were most often used to distinguish land use and ground cover.
An unsupervised classification was run for each of the first 12 CS to automatically cluster the class signatures. Of the unsupervised classification techniques, Iterative Self-Organizing Data Analysis Technique (ISODATA) was used. The clusters formed were regrouped, using Ward's hierarchical clustering technique, which is designed to optimize the minimum variance within clusters (Ward 1963). The resulting classes were identified on the basis of information drawn from field surveys and earlier land use maps. The results of each CS were compared creating error matrics and overall classification accuracy. One CS was selected which showed a superior result during unsupervised classification.
Supervised classification was performed on the selected CS 11, employing the Bayesian Maximum Likelihood Classifier (MLC). MLC, a parametric decision rule, is a well-developed method from statistical decision theory that has been applied to the problem of classifying image data (Niblack 1985, Settle and Briggs 1987). At first, training signatures for identifiable classes were established on the basis of field information. They were then evaluated to make sure there was suitable discrimination of individual classes. After obtaining a suitable grouping for satisfactory discrimination between the classes during signature evaluation, the final classification was carried out (Table 1, Class 13).
The classification accuracy was calculated using a confusion or error matrix (Lilesand and Kiefer 1994), which showed the accuracy of both the producer and the user. The producer's accuracy is obtained by dividing the number of correctly classified pixels in a class by the total number of pixels of that class as derived from the ground data. The user's accuracy is obtained by dividing the total number of pixels that were classified in that class.
The output digital elevation model (DEM) was a 16-bit single band raster image as the interpolated elevation ranged from 15 m to 1,671 m. Each pixel's value represented the elevation at that point. The output was evaluated by randomly selecting 100 sample points from the interpolated DEM, and comparing the interpolated values with corresponding points in the topographic map, using the Pearson product moment correlation coefficient. This correlation coefficient provides a measure of spatial correlation and spatial association between two variables (Haining 1990). A linear association was observed between interpolated elevation and real topographic elevation values with a correlation coefficient (R2) of 0.97.
The 16-bit DEM was rescaled to 8-bit raster so that it was compatible with that used for Landsat-TM. A regrouped DEM divided the terrain into 100 m intervals.
Both unsupervised and supervised classification were used to classify the Landsat TM data, with combinations of various original bands, single ratio bands, NDVI, principal component bands and DEM to test the classification results. A total of 12 CS for unsupervised, and one CS for supervised, classification were tested. Two descriptive statistics, namely producer's accuracy and user's accuracy, were calculated to estimate the Overall Classification Accuracy (OCA) for each band combination. For the unsupervised technique, the OCA ranged from about 58 to 78 % (Table 1). Among 12 band combinations, combinations which included DEM gave better results most of the time. For supervised classification, an OCA of 82.3% was obtained.
The Universal Soil Loss Equation (USLE) is an empirical model that is widely used all over the world for the assessment and prediction of soil erosion due to water runoff. When the equation was originally developed, it was not intended to be valid for a large area. However, satisfactory results were reported by various researchers who used it on a large scale for watersheds. One was Mellerowicz et al. (1994), who comments that it is still by far the most widely used method, but it is necessary to adjust the USLE factors to a specific location for reliable results.
Fig. 2 shows the GIS methodology for estimating the amount of soil loss.
Wischmeier and Smith (1978) gave the modified equation of USLE as follows:
E = R.K.L.S.C.P
Where, E is the mean annual soil loss (mt/ha/yr), R is the rainfall erosivity factor, K is soil erodibility, L is slope length, S is the steepness of the slope, C is crop management and P is the erosion control practices.
The R factor (rainfall erosivity index) is the principle function of USLE. The erosivity of a rainstorm is the combined function of its intensity and duration, and of the mass, diameter and velocity of raindrops. Previous studies have shown that maximum rainfall intensity for 30 minutes (El30), expressed as kinetic energy, has the best correlation with erosivity. Several location-specific equations have been proposed by various researchers to compute the R-factor, as cited in Morgan (1996), such as:
where P is the mean annual preci-pitation in mm.
The first equation gave very high erosivity, whereas the third gave a very low one. Many previous studies made in Thailand have used Equation 4. Morgan suggests that the average of Equations 1 and 2 gives a reliable result when computing the R-factor.
Soil erodibility is the function of physical characteristics of soil and its management, including both land and crop management. The erodibility index (K) can be obtained from an already developed nomograph*. The erodibility index was calculated for each soil series in this study, using the following equation:
All M, a, b and c values were extracted from soil reports of individual soil series.
According to this equation, the K value is calculated on the basis of the organic matter content, the type of soil structure, the soil permeability and the soil texture.
In this study, the slope length for different types of land use was determined, using information from the field survey. Data on aspect, land use and slope maps were overlaid to map the slope length. L and S were treated as a combined factor to find the LS index. The LS factor was computed for slope gradients of less than 8%, using the following equation.
C is the crop management factor, while P is the erosion control practice or conservation factor. Interviews conducted during the field survey indicated that little or no land management is practiced in the area. Hence, the C index value used by Funnpheng et al. (1991) was also used in this study.
The land use types were derived from the digital image classification of Landsat TM data, as mentioned in the previous section.
The amount of soil loss, and division of an area into classes according to this amount, does not follow an absolute value, but is relative in each area. The Sakae Krang watershed was divided into eight classes, depending on their annual soil loss (Table 2 and Fig. 3).
The soil losses each year in the area ranged from only 0.04 mt/ha/year to as much as 39 mt/ha/yr. Most of the area fell within the slight (38%) and moderate category (24%) of erosion hazard severity. Some slight erosion was observed in rainfed lowland paddy fields and ratoon crops of sugarcane. Moderate erosion was mostly seen in upland areas. About 30% of the watershed had a very slight erosion hazard. Nearly 10% of the watershed area had severe to extremely severe erosion. Such areas were found in foothills covered in dry dipterocarp forest. This sheds its leaves during the dry season, allowing rain to splash directly on the soil surface and detach soil particles.
Soil losses were also computed for each type of land use type. The results will serve as a guide as to what kind of conservation measures and agronomic practices should be adopted in the area. The average rate of soil loss associated with different types of land use and vegetation cover is presented in Table 3.
Soil nutrients are indispensable for crop growth and yield. As land use has become more intensive in the area, the application of chemical fertilizers has become common. Presumably, the rate of fertilizer application follows the recommendation of extension staff, who are perhaps striving to increase the yields of low-income farmers. Thus, the examination of nutrients lost in soil is a very important aspect of realizing the value of soil resources.
The soil nutrient balance was calculated, based on the spatial distribution of various crops, the nutrient uptake by the different crops in each cropping season, the nutrients lost by erosion, and the nutrients available in the soil system. Nitrogen (N), Phosphorus (P) and Potassium (K) are the main plant nutrients. The results indicated that average rates of nutrient loss in the area were 10.9 kg N, 0.16 kg P and 1.02 kg K per hectare per year.
In the case of annual N depletion, 54% of the area was found to have losses of less than 5 kg N/ha/year. Thirteen percent of the area had annual loses of 6 - 10 kg N, 4% of the area had losses of 11 - 15 kg N and about 5% area of the area had annual losses of 15 - 150 kg N (Fig. 4). In the case of P depletion, about two thirds of the watershed area was found to have losses of less than 0.1 kg/ha/year. About 4% of the area had losses of 0.1 to 0.2 kg/ha P per year, and 6% had losses of 0.2 to 0.6 kg/ha P per year. Similarly, in the case of K loss, 49% of the area was found to have losses of less than 0.2 kg/ha each year. About 16% of area had annual losses of 0.2 to 0.6 kg/ha K, and about 13% of the area lost 0.7 to 1.5 kg.
This Bulletin illustrates the application of SRS, DEM and GIS for monitoring a watershed. The case of the Sakae Krang watershed illustrates the use of high-resolution satellite data, Landsat TM, through digital image processing techniques. This can be combined using DEM, derived from reliable data and maps on land use. GIS technology can be employed to create DEM, and to map soil erosion by using the USLE. The soil nutrient balance can be determined, if data about nutrient loss by erosion, crop uptake and the nutrient supply in the soil system are available. Furthermore, the integration of SRS and GIS is a very useful technique for land capability and suitability classification, and in planning new alternative land use options.
Figure 1 Location of the Study Area
Figure 2 Gis Methodology to Estimate Soil Losses
Table 2 Annual Soil Loss Ratings
Figure 3 Soil Loss Map of Sakae Krang Watershed, Thailand
Table 1 Classification Scheme and Overall Classification Accuracy
Figure 4 Nitrogen Loss Map
Table 3 Rate of Soil Loss for Different Types of Land Use
Download the PDF. of this document, 246,900 bytes (241 KB).