Numerical Modelling of the Topographic Wetness Index : An Analysis at Different Scales

A variety of landscape properties have been modeled successfully using topographic indices such as topographic wetness index (TWI), defined as ln(a/tanβ), where a is the specific upslope area and β is the surface slope. In this study, 25 m spatial resolution from digital elevation models (DEM) data were used to investigate the scale-dependency of TWI values when converting DEMs to 50 and 100 m. To investigate the impact of different spatial resolution, the two lower resolution DEMs were interpolated to the original 25 m grid size. In addition, to compare different flow-direction algorithms, a second objective was to evaluate differences in spatial patterns. Thus the values of TWI were compared in two different ways: 1) distribution functions and their statistics; and 2) cell by cell comparison of DEMs with the same spatial resolution but different flowdirections. As in previous TWI studies, the computed specific upstream is smaller, on average, at higher resolution. TWI variation decreased with increasing grid size. A cell by cell comparison of the TWI values of the 50 and 100 m DEMs showed a low correlation with the TWI based on the 25 m DEM. The results showed significant differences between different flow-diretction algorithms computed for DEMs with 25, 50 and 100 m spatial resolution.


Introduction
Widely available digital elevation models (DEMs), together with corresponding tools for spatial analysis, have found extensive use in the development of research in many areas of the environmental sciences, including agriculture, hydrology, ecology, geography and branches of engineering using topographically-dependent variables.The tools for digital elevation modeling [1-4, among others] contain various options for the analysis of topographical attributes, such as algorithms for extracting drainage networks [5][6][7], topographic indices [8], indices of sediment transport capacity [9] and erosive potential [9].
The most common sources of data for digital elevation modeling are [10]: topographic bases with isolines and points denoting maxima and minima, and transformed field surveys undertaken, for example, with the use of GPS, photogrammetry, stereoscopy from remote sensing and radar interferometry, given in regular grid structures, triangular grids and isoline nets [9].
Data are generated at different scales and spatial resolutions.Quantities derived from topographic data need to take account of the spatial limits imposed by the respective scales.At present, a number of spatial resolutions are available for public use.Higher spatial resolutions are obtained from photogrammetry and laser profiling, while coarser spatial resolutions are freely available through the Shuttle Radar Topographic Mission (SRTM) [11].The accuracy of data generated at different spatial resolutions depends exclusively on how closely a calculated measure lies to its true value: the reliability of the data is therefore linked to its accuracy.In this context, the higher the spatial resolution of a DEM, the greater is the tendency to accuracy.When analyzing data from a DEM, a fundamental step is the evaluation of the precision of the derived data; the user must know how trustworthy his derived data is.
A number of authors have studied the effects of spatial resolution on topographic data derived from DEMs, including topographic derivatives of the first, second and third orders [10,11], effects of spatial resolution in mod-els for the direction of flow and accumulated overland flow [12] and variations in the topographic index [13,14].
A diagram relating the information in DEMs to hydrological information shows a qualitative and quantitative decrease in hydrological information with reduction in spatial resolution of DEMs [10].Two situations can determine the quality of hydrological information in high resolution models: 1) the topography controls the hydrological flows, giving an increase in hydrological information; 2) soil characteristics determine the preferred hydrological pathways, and in this case there is no increase in hydrological information.
In this context, we attempted to evaluate the effects of different flow-direction algorithms and spatial resolutions on the TWI.The analysis of the scale effect on topographic parameters in the drainage basin of the Potiribu River [15] helps develop an understanding of their magnitude and the way in which indices, such as soil saturation, are affected.

The Potiribu River Basin
The study area lies within the Taboão River basin which drains to the Potiribu River, lying to the northwest of the Brazilian State of Rio Grande do Sul.The area is representative of the basaltic plateau of the River Paraná basin (Figure 1) where basaltic outflows extend over an area of about 230,000 km 2 [16].The Potiribu River is a tributary of the Rio Ijuí, itself a left-bank tributary of the River Uruguay.The total area of the Potiribu River basin is 563 km² and consists of a number of nested sub-basins monitored as part of the Potiribu Project.Physical characteristics of these nested sub-basins are given in Table 1.
Altitude varies between 320 and 700 m.The relief is undulating with gentle slopes between 3% and 15%.Soils are fairly wet and well-developed as a consequence of the abundant rainfall regime.Soils can be classified as purple and dark-red latosol.Despite high clay content (with percentages of clay in excess of 60%), there is a well-developed drainage network as a result of soil macropores [17].
The climate of the Potiribu River basin is classified as Cfa: a mild, mesothermic and very damp climate with no well-defined dry season.In this classification, mean maximum temperatures are greater than 22˚C and mean minimum temperatures are between -3˚C and 18˚C.Mean temperatures in the region fluctuate between 18˚C and 19˚C, with July the coldest month (13˚C -14˚C) and January the hottest (24˚C).The mean maximum temperature is about 32˚C and the mean minimum about 8˚C.
Rainfall is determined by the interaction of large air masses [15] resulting from: 1) turbulent southerly air currents associated with polar anticyclone activity, termed the polar front.This polar front occurs fairly regularly with periodicity between 4 and 10 days, and is accompanied by rainfall events that are generally long and of moderate intensity, lingering for several consecutive days; 2) turbulent westerly currents associated with lines of tropical that cross the Southern Region from the middle of spring to the mid-October.Depressions are induced by small tropical anticyclones of Amazonia.These currents give rise to thunderstorms and convective rainfall that is generally intense but of short duration.
Researchers at the Institute of Hydraulic Research [15,

Methodology
The effects of different spatial resolutions on TWI have been assessed by means of histograms, descriptive statistics and regression analysis.The DEM with 25 m spatial resolution was obtained by aerial photogrammetry; that at 50 m resolution by digitizing contours of topographic maps [19], and that at 100 m by resampling of data from the Shuttle Radar Topographic Mission (SRTM) [22].Indices derived from the DEMs were obtained using different algorithms to determine flow direction: 1) D8 [21], 2) D MF [22] and 3) D [23].
The TWI, given by Equation ( 1), defines areas of saturated soil typically found in geomorphologically convergent segments, where a is the area of specific contribution based on flow-direction (summed over the area upstream of the cell) and β is surface slope.It is determined from the flow direction and accumulated runoff.
The specific contributing area is related to the concept of accumulated runoff and takes into account the complexities of hill-slope form.The curvature of slopes both in the plane and in profile effectively determines the hydro-sedimentological behavior of erosive processes.It is defined as the accumulated drainage area divided by the width of the downstream contributing area (m 2 •m -1 ).In the case of convex slopes for which accumulated flow is divergent, the specific contributing area tends to diminish.For concave slopes, the specific contributing area tends to increase, giving rise to a rapid increase in accumulated downstream flow, defined as a specific contributing area [24].
In the D8, the flow directions are in the form of a regular grid in which each cell's value corresponds to one of the eight possible directions.Obtaining the direction of flow for each pixel to one of its eight neighbors can be done automatically, based on the difference in level between them, weighted by distance.Thus a value is attached to each pixel indicating one of the eight flow directions.Following the flow directions, the number of upslope cells is calculated which drain to a downslope cell, so that the accumulated runoff is calculated.
However, the flow patterns from cell to cell are almost never as simple, since a cell can receive flows from several cells, so that the accumulated flow may be transferred to several other cells [25].Therefore a proposed multiple flow-direction algorithm (D MF ) [22], in which the flow transferred to each downslope cell is proportional to the product of the distance over which flow is accumulated and a geometrical weighting factor which depends on the flow direction.
Since the D8 deals best with surface flow largely concentrated in the main channel, and the D MF deals with dispersed surface flow, an algorithm with an infinite number of possible flow directions was proposed [23], to overcome the limitations of both D8 and D MF .The D gives an unlimited number of flow directions using triangular facets in a 3 × 3 pixel window.Thus the flow distribution is proportioned amongst downslope pixels according to the slope of each triangular facet.For a hypothetical conical surface, there is a concentration of surface flow in the D8, and a substantial dispersion of surface flow in the D MF .The D gives no substantial spreading of surface flow nor a concentration of flow in a single main channel, such as that occurring in the D8 and D MF .

Results and Discussion
The analyses reported here have two distinct aspects: 1) the effect of the flow-direction algorithm and 2) the effects of different spatial resolutions.At 25 m resolution, along the main channel and its tributaries and in these areas the highest values were found.Values found on the slopes were much smaller than those in the main channels.When the D MF was used, TWI showed a much smoother behavior than that found using the D8.The main channels could be seen but the tributaries could no longer be distinguished, having values very similar to those on the hill-slopes.The largest values found in the main channels when the D MF was used were smaller than those given in the main channels by the D8 algorithm.The TWI generated by the D gave results very similar to those given by the D8, although tributaries to the main channel were less in evidence as a consequence of the considerably greater dispersion of runoff.Descriptive statistics of indices generated by the three flow-direction algorithms were very similar, although the summary statistics increased with the pixel size (Table 2).The frequency distribution shows significant differences between the flow-direction algorithms (p-value < 0.05) and between spatial resolutions (p-value < 0.05) (Figure 2).
For indices generated by the D8, the mean gradually increases from 7.45 to 8.35 for resolutions 25 and 100 m  Taking the 25 m spatial resolution as the most representative of TWI, the D8 gives increases of 5.23% and 12.08% in the mean value of indices generated at resolutions 50 and 100 m respectively.For the D MF , increases are 3.68% and 8.65% at resolutions 50 and 100 m, while for the D the increases are 6.71% and 12.88%.
The greatest maximum values were found in indices generated at 50 m spatial resolution, whilst the smallest maxima values were found at 100 m resolution.For the D8, variations from 10.20% to -10.36% were found when comparing resolutions 50 and 100 m respectively, with the 25 m resolution.For the D MF , these changes were between 4.01% and -9.93%, respectively, and for the D, between -3.32% and -12.33%.The magnitude of these variations in maximum values of TWI is probably a consequence of smoothing the DEM and the slope indices, as for each slope segment at 100m resolution there are 2 and 4 slope segments at resolutions 50 and 25 m respectively.
A pixel-by-pixel comparison of root mean square error (RMSE) at the different spatial resolutions shows where the greatest variations in scale effects occur.A bilinear interpolator was used to resample the topographic indices from 50 and 100 m to 25 m.The greatest differences are found along the main channel and its tributaries, whilst differences on hill-slopes are smaller.In the case of D8, the RMSE between 25 and 50 m was 0.96 ± 1.42 and between 25 and 100 m, the RMSE was 1.44 ± 1.60.For the D, the RMSE was 0.92 ± 1.30 and 1.43 ± 2.51 at spatial resolutions from 25 to 50 m and 25 to 100 m, respectively.Considering D MF , the RMSE between 25 and 50 m, and between 25 and 100 m, were smaller than the differences generated by the D8 and D.
It can be concluded, from these analyses at different resolutions, that the greatest differences are found between indices generated at 25 and 100 m resolutions using the D8 and D, since the greatest flows are propagated in the main channels.In this context, indices generated by the multiple-flow D MF show differences that are substantially smaller than those given by the D8 and D algorithms, mainly as a result of the way in which flow is distributed.Figure 3 shows a box-plot diagram

Concluding Remarks
The spatial resolution of DEMs is a major influence on the values of indices such as topographic wetness index, derived from numerical models.Equally, determination of direction of flow and of specific contributing area depends on the flow-direction algorithm used.Algorithms such as the simple flow D8, the multiple-flow D MF and the algorithm based on triangular facets such as D, can give rise to quite different results.
From the analyses using different algorithms for flow direction, it can be concluded that 1) the D8 concentrated runoff along the main channel, whilst the D MF algorithm gave dispersal of runoff, and the D gave intermediate results; 2) descriptive statistics showed subtle differences, mainly in terms of minimum values, means, lower quartile and the median.
Analyses of the effects of different spatial resolutions showed significant variation between descriptive statistics, with quite different results.Little correlation was found between indices at resolutions 25, 50 and 100 m, in agreement with other results [14].However in more than 75% of cases the differences were quite subtle, with differences in TWI lower than two points.A pixel-bypixel analysis showed that the greatest differences occurred mostly along the main channel, independently of the algorithm used to determine flow direction.There was also greater concentration of differences along tributaries of the main channel, when the D8 and D were used.Differences in the topographic index tended to be smaller in the headwater areas of the drainage network and on hillslopes, where runoff is smaller.
It can be concluded from the results that there are many differences between values of TWI generated at different spatial resolutions and by different flow-direction algorithms, but it is not clear that results obtained at more detailed spatial resolutions are necessarily better.The effects of scale can become more accentuated or less so, depending on the topographic influences on hydrologic processes.

Figure 1 .
Figure 1.Position of the Potiribu River basin in Brazil.The broken line marks the limit of the basaltic plateau of the Parana River basin.

Figure 2 .
Figure 2. Frequency distribution of topographic wetness index (TWI) for different spatial resolutions (25, 50 and 100 m) using 8, D MF and D∞ flow-direction algorithms.D

Figure 3 .
Figure 3. Ranges of root mean square error (RMSE) between spatial resolutions from 25 to 50 m (1) and from 25 to 100 m (2) for TWI generated by the D8 algorithm; from 25 to 50 m (3) and from 25 to 100 m (4) using D MF algorithm; from 25 to 50 m 5) and from 25 to 100 m (6) using D∞ algorithm.(

Figure 4 .
Figure 4. Correlation of the topographic wetness index (TWI) between different spatial resolutions.