HYDROLOGICAL PROCESSES
Hydrol. Process. 20, 1385– 1405 (2006)
Published online in Wiley InterScience (www.interscience.wiley.com). DOI: 10.1002/hyp.6098

A downscaling method of topographic index distribution
for matching the scales of model application andparameter identification
N. R. Pradhan, Y. Tachikawa* and K. Takara
Disaster Prevention Research Institute, Kyoto University, Gokasho, Uji, Kyoto 611-0011, Japan

Abstract:
Higher resolution topographic information contained in the topographic index of TOPMODEL is lost when digital
elevation models (DEMs) with a coarse grid resolution are used; thus, the topographic index is scale dependent,
demonstrating identified model parameter values to be dependent on DEM resolution. This makes it difficult to use
model parameter values identified through a different resolution of TOPMODEL. The inconsistency is the result of thedifference between the scale at which the model parameters are identified and the scale at which the model is applied.
To overcome this problem, scale laws that govern the relationship between the resolution of digital elevation data andgeomorphometric parameters of the topographic index were analysed and a method to downscale the topographic index
distribution developed to account for the difference in scales between model application and parameter identification.
The method to downscale the topographic index is composed of two ideas: one involves introducing a resolution
factor to account for the scale effect in upslope catchment area per unit contour length in the topographic index; theother utilizes a fractal method through steepest slope scaling to account for the scale effect on slopes. This method
successfully derived a topographic index distribution of a fine-resolution DEM by using only a coarse-resolution DEM.
The method has been applied successfully to the Kamishiiba catchment (210 km2 ) in Japan and has demonstrated thatthe downscaled topographic index distribution derived using a 1000 m grid DEM is very similar to the topographic
index distribution derived via fine-target-resolution DEMs. The method is then coupled with a TOPMODEL simulation
to match the scales of model application and parameter identification. It is shown that the simulated runoff from thedownscaled TOPMODEL applied at 1000 m resolution of the Kamishiiba catchment, with the same set of effective
parameter values derived from 50 m resolution DEM, matched the simulated runoff in the 50 m DEM resolution
TOPMODEL. It was also shown that TOPMODEL coupled with the downscaling method of the topographic index
accurately simulated runoff for different rainfall events in the catchment without recalibration. Copyright  2006 John
Wiley & Sons, Ltd.
KEY WORDS

scale invariance; downscale; fractal method; ungauged basin; topographic index; TOPMODELINTRODUCTION
Despite the enormous capacity of today’s (and tomorrow’s) information technologies, the complexity of the
Earth’s surface is such that the most voluminous descriptions are still only coarse generalizations of what is
actually present (Goodchild, 2001). This implies that the need for continued and sustained research on scaleissues is required. Since the introduction of the first blueprint of a distributed hydrological model (Freeze andHarlan, 1969), the desire to develop increasingly realistic distributed models has been spurred by requirements
to improve forecasting changes in hydrological behaviour due to a variety of land use and climate change andfor hydrological predictions in ungauged basins. However, in macroscale/mesoscale hydrological modelling,
the concept of the aggregation approach, which assumes governing equations valid at the small scale can* Correspondence to: Y. Tachikawa, Disaster Prevention Research Insitute, Kyoto University, Gokasho, Uji, Kyoto, Japan.
E-mail: tachikawa@flood.dpri.kyoto-u.ac.jp

Copyright  2006 John Wiley & Sons, Ltd.

Received 16 June 2004
Accepted 12 September 2005

1386

N. R. PRADHAN, Y. TACHIKAWA AND K. TAKARA

be applied at larger scales using effective parameter values, fails because effective parameter values are
scale dependent. To resolve this issue, a method that enables use of model parameters at different scales ofapplications is required. One method to solve the problem is to develop a scaling theory for effective model
parameters with regard to modelling scales, and another is to match the scales of model application andparameter identification through introducing a scaling method for a hydrological model. This paper adopts
the latter idea and discusses the method in the TOPMODEL framework.
As the spatial extent is expanded from point experiments to larger watershed regions, the direct extension of
the point models requires an estimation of the distribution of model parameters and process computations overthe heterogeneous land surface. If a distribution of a set of spatial variables required for a given hydrologicalmodel (e.g. surface slope, soil hydraulic conductivity) can be described using a joint density function, thendigital elevation models (DEMs) and geographical information systems (GISs) may be evaluated as a toolfor estimating the joint density function. The question still remains as to whether current GISs and available
spatial data sets are sufficient to estimate these density functions adequately.
Several researchers (Quinn et al., 1991; Iorgulescu and Jordan, 1994; Wolock and Price, 1994; Zhang andMontgomery, 1994; Bruneau et al., 1995; Saulnier et al., 1997; Mendicino and Sole, 1997) have discussed
the effects of DEM map scale and data resolution on the distribution of the topographic index, concludingthat there is interdependence between DEM scale and topographic index distribution. Lack of methods forthe translation of the scale dependence relations into effective hydrological models poses a serious problem
for the ungauged basins of developing countries, where either only coarse-resolution DEM data are available
or where the information gained at one scale is to be utilized in making predictions at other (usually lower)
scales. The problem of transferring information gained at one scale for making predictions at a differenthydrological scale is a scaling problem (Beven, 1995).
Band and Moore (1995) point out that higher frequency topographic information is lost because the larger
sampling dimensions of the grids act as a filter. This is one of the fundamentals of the scale problem,
and the prediction in ungauged basins becomes increasingly difficult and inaccurate without a method toresolve this issue. Thus, a downscaling/disaggregation method is needed to acquire a more realistic subgridscale parameterization in hydrological modelling to resolve the limitations of the aggregation method in ascale-dependent hydrological model while scaling small-scale hydrological processes and parameters to larger
scales.
Analysis of physically based models such as SHE (Abbott et al., 1986a,b), IHDM (Rogers et al., 1985),
and TOPMODEL (Beven and Kirkby, 1979) demonstrated that these models have, in theory, no need ofpreliminary calibration, since the model parameters offer clear physical meaning, which makes estimation of
their values possible with regard to knowledge of the basins characteristics. However, it seems to be necessary
to differentiate between ‘physically based’ in the sense of being based on defined assumptions and theories, and‘physically based’ in the sense of being consistent with observations (Beven, 2002). A hydrological modelthat is physically based in theory but not consistent with observations often results from the discrepancy
between the scale at which the model parameters are identified and the scale of model applications.
Even in a fully distributed physically based hydrological model, the differential equations concerning thevarious hydrological processes (overland flow, infiltration, percolation, etc.) are solved for single cells into
which the basin is subdivided, where, by introducing a conceptualization of the phenomenon for hydrologicalprocesses, the heterogeneity of the hydrological quantities inside the cell are ignored. Conceptualizations
introduced in this way result in different performances in models with variations in the assumed scale (Woodet al., 1988).
TOPMODEL, in practice, represents an attempt to combine the computational and parametric efficiency ofa distribution function approach with the link to a physically based process description . Though it is used fora wide variety of applications, its dominating geomorphometric parameters that account for the hydrologicalsimilarity condition are also strongly influenced by the resolution of the DEM used. This results in parameterinconsistency and predictive uncertainty across scales.
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

DOWNSCALING METHOD OF TOPOGRAPHIC INDEX DISTRIBUTION

1387

In this study we focus on the influence of DEM resolution on dominating geomorphometric parameters, slope
angle and upslope contributing area, which are considered as the main controls in a number of hydrologicalprocesses, and by incorporating scaling laws we developed a method to downscale the topographic index
of TOPMODEL. By using this method, the topographic index distribution of a fine-resolution DEM is
successfully derived by using only a coarse-resolution DEM (Pradhan et al., 2004). We then coupled themethod to downscale the topographic index distribution with TOPMODEL to match the scales of model
application and parameter identification. This allowed the use of topographic information gained at a coarser
DEM scale in making hydrological predictions at a finer DEM scale.

DEPENDENCE OF TOPOGRAPHIC INDEX DISTRIBUTION ON DEM RESOLUTION
Theory of TOPMODEL
The topographic index (Kirkby, 1975) of TOPMODEL is defined asa
TI D lntan ˇ

1where a is local upslope catchment area per unit contour length and ˇ is slope angle of ground surface. In
the TOPMODEL framework, the topographic index is used to distribute the local saturation deficits, givenknowledge of the mean storage deficit (Beven and Kirkby, 1979; Beven, 1986, 2000a). Equation (2) showsa relationship between St [L], the spatial mean storage deficit and S(i,t) [L], local saturation deficit at eachlocation i in a catchment asaiSi, t D St C m  m ln
2
T0 tan ˇiwhere tan ˇi is slope angle, T0 [L2 T1 ] is the lateral (horizontal) transmissivity at zero storage deficit (i.e.
when the soil is just saturated to the surface), ai is the area draining through i per unit contour length [L2 ],m is a decay factor of saturated transmissivity of soil with respect to saturation deficit with dimensions oflength [L], and  is a constant for the basin given by



1
ai
D
ln
dA
3
A A
T0 tan ˇi
According to Beven (1986, 2000a), subsurface contributions to streamflow Qb t [LT1 ] can be derived
from Equation (2) asQb t D T0 e eSt/m
4
Considering lateral transitivity to be constant in a catchment or subcatchment, the key factor forhydrologically similar conditions is the distribution function of the topographic index. Higher frequency
topographic information contained in topographic indexes is lost, as the larger sampling dimensions of
the grids act as a filter. This results in the hydrological similarity condition for computing the combinedsoil–topographic index in Equation (2) to vary accordingly with the variation in the DEM resolution used. To
overcome this problem, this research has developed a method to downscale topographic index distribution.
Influence of DEM resolution on topographic index
Figure 1 shows the density function of the topographic index at four different DEM resolutions in the
Kamishiiba catchment (210 km2 ) in Japan without taking into account the scale effect. A distinct shift of the
topographic index density functions towards higher values is seen in Figure 1 as the resolution of the DEM
becomes coarser, using the steepest slope algorithm. Saulnier et al. (1997) have also demonstrated a similar
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1388

N. R. PRADHAN, Y. TACHIKAWA AND K. TAKARA

0.35

Topographic index from 50m DEM resolution

0.3
Density

0.25

Topographic index from 150m DEM resolution

0.2
Topographic index from 450m DEM resolution

0.15
0.1

Topographic index from 1000m DEM resolution

0.05

35

32.8

30.6

28.4

26.2

24

21.8

19.6

17.4

15.2

13

10.8

8.6

6.4

4.2

2

0

Topographic index value [ln(m2)]
Figure 1. Effect of DEM resolution on density distribution of topographic index

Table I. DEM resolution effect on spatial mean of topographic index  in the Kamishiiba
catchment
DEM resolution (m)Spatial mean of topographic
index  lnm2 

50
6Ð076

150
7Ð423

450
9Ð222

600
9Ð622

1000
10Ð353

shift of topographic index density functions towards higher values as the resolution of the DEM becomescoarser when using a multi-directional algorithm (Quinn et al., 1995). This is a clear indication of the lossof higher frequency topographic information, as the larger sampling dimensions of the grids act as a filter.
Table I shows the distinct effect of DEM resolution on the spatial mean value of the topographic index  in
Equation (5) as



1
ai
D
ln
dA
5
A A
tan ˇiAnalysing Figure 1 and Table I, we can readily conceive the error in hydrological prediction in ungauged
basins using a DEM at a lower resolution when using parameter values identified at different DEM resolutions.
DEVELOPMENT OF A METHOD TO DOWNSCALE UPSLOPE CATCHMENT AREA FOR THE
SOLUTION OF DEM RESOLUTION EFFECT
In a DEM-based distributed hydrological model, upslope catchment area at a point is the number of pixels
draining through that point (Jenson and Domingue, 1988; Rodriguez-Iturbe and Rinaldo, 1997). In hydrologicalgeomorphology, upslope catchment area is a key variable because of its intrinsic capability to describe
the nested aggregation structure embedded in the fluvial landforms and its important physical implications(e.g. Leopold and Maddock, 1953; Rodriguez-Iturbe and Rinaldo, 1997). Most physically based models ofhydrological and geomorphic processes rely on spatially distributed characterization of drainage area (e.g.
Beven and Kirkby, 1979; O’Loughlin, 1986; Woods et al., 1997). Any upslope catchment areas smaller than
the grid area of the DEM resolution used are completely lost. Thus, a portion of the upslope catchment areathat can be well defined by a finer DEM resolution gets completely filtered out if that portion of the upslope
catchment area is less than the grid area of the coarse-resolution DEM used.
Figure 2 shows the frequency distribution of the upslope catchment area at four different DEM resolutionsin the Kamishiiba catchment (210 km2 ) in Japan without taking into account the scale effect. Distinct decreases
in the finer values of upslope catchment area are observed in Figure 2 when the resolution of DEM becomescoarser. In Figure 2, the smaller upslope catchment areas less than 1 km2 that appear in over 97% of thecatchment in 50 m DEM resolution are completely lost when switched to a 1000 m DEM resolution. Figure 3
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1389

Fraction of catchment area

DOWNSCALING METHOD OF TOPOGRAPHIC INDEX DISTRIBUTION

1.1
1
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0

At 50m DEM resolution
At 150m DEM resolution
At 450m DEM resolution
At 1000m DEM resolution

1

15

29

43

57

71

85

99

113 127 141 155 169 183 197 211

Upslope contributing area (km2)
Figure 2. Comparison of distribution functions of upslope catchment area obtained from different DEM resolutions in the Kamishiiba
catchment (210 km2 )LegendUpslope contributing area valueLess than 1 km2
Equal or greater than 1 km2

(a)

(b)

Figure 3. Comparison of upslope catchment area from (a) 50 m DEM resolution and (b) 1000 m DEM resolution in the Kamishiiba catchment(210 km2 )further clarifies this fact. Figure 3a shows that, at 50 m DEM resolution, most of the catchment is displayed
as possessing upslope catchment area of less than 1 km2 , whereas Figure 3b shows that, at 1000 m DEM
resolution, there is no grid having an upslope catchment area less than 1 km2 .
In fact, the smallest upslope catchment area derived from a DEM resolution is a single grid of the DEM
at that resolution. From this point of view, we introduced the total number of sub-grids Ns within a coarseresolution grid (see Figure 4) to derive scaled upslope catchment area asCi, scaled D

Ci
Ns If

6

where the suffix i is a location in a catchment, Ci, scaled is the scaled upslope catchment area at a point I, and
If is introduced as an influence factor.
Figure 2 shows that, as the upslope catchment area increases, the distributions of the upslope catchment
area values given by coarse- and fine-resolution DEMs converge; thus, the influence of Ns on Ci must
gradually decrease in Equation (6) as Ci becomes larger. For this reason we introduced influence factor If in
Equation (6). From the discussion of influence of Ns on Ci , the following three points are proposed:
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1390

N. R. PRADHAN, Y. TACHIKAWA AND K. TAKARA

W*i

1

2

3

4

5

6

7

8

9

Flow Direction

Wi

Coarse grid resolution

Sub-grid

2

Ns =

Wi
= Rf2 (=9)
W*iFigure 4. Concept of Ns as number of sub-grids within a coarse grid resolution

1. At the catchment divide portion where the upslope catchment area in a coarse-resolution DEM is asingle coarse-resolution grid area, the value of the influence factor If in Equation (6) is equal to unity,
demonstrating that Ns has complete influence on Ci .
2. Considering the upslope catchment area given by a coarse-resolution DEM and target fine-resolution DEM
are equal at the outlet of the catchment, the value of influence factor If in Equation (6) is equal to 1/Ns ,showing that Ns has no influence on Ci .
3. Exponential decay of the influence factor is taken from the value 1 defined in point 1. to 1/Ns defined inpoint 2. as the upslope catchment area increases. Thus, If is defined as
If D e[1Ni H]/No

7

where Ni is the number of coarse-resolution grids contained in the upslope catchment area at a location iin the catchment and No is the number of coarse-resolution grids contained in the upslope catchment areaat the outlet of the catchment. H in Equation (7) is introduced as a parameter to define the variation in theinfluence factor as a function of No and Ns . Considering the influence of Ns on Ci in Equation (6) is almost
negligible at the outlet of the catchment, where Ni D No and If D 1/Ns , the value of H can be obtained from
Equation (7) by substituting Ni D No and If D 1/Ns as
HD

No
ln Ns
No  1

8

Finally, we developed a method to downscale the upslope catchment area from Equations (6) and (7) as
Ci, scaled D

Ci[1Ni H]/No

Ns e

9

The analysis and derivation of downscaling method of upslope catchment area was developed from theassumption that upslope catchment area at a point is the number of pixels draining through that point. Although
the multiple flow-direction approach is an approximate form of subgrid-scale flow pathway interpolation(Quinn et al., 1991), it is also affected by the loss of information of finer resolution upslope catchment area
when using coarse-resolution DEMs. This is because the values of distributed upslope catchment area at thecatchment divide, according to weighted percentages relative to the value assumed by the slope, are higher
in a coarse grid cell than in a fine grid cell.
The method to downscale the upslope catchment area is applied to the Kamishiiba catchment (210 km2 ).
Using Equation (9), the upslope catchment area is downscaled from 1000 m DEM resolution to various DEM
grid resolutions as shown in Figure 5. In contrast to Figure 2, Figure 5a–d shows that the distributions of
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1391

DOWNSCALING METHOD OF TOPOGRAPHIC INDEX DISTRIBUTION

At 50m DEM resolution
0.9

Downscaled from 1000m DEM resolution
to 50m DEM resolution

At 150m DEM resolution

0.9

Downscaled from 1000m DEM resolution
to 150m DEM resolution

0.8

209

196

183

170

157

144

131

118

92

105

79

66

53

40

27

1

14

209

196

183

170

157

144

131

118

92

105

79

66

53

40

27

1

Upslope catchment area (km2)

Upslope catchment area (km2)

(b)

1
0.9

At 450m DEM resolution

0.8
Downscaled from 1000m DEM resolution
to 450m DEM resolution

0.7

Fraction of catchment area

1.1

1.1

1
0.9
0.8

At 600m DEM resolution

0.7

Downscaled from 1000m DEM resolution
to 600m DEM resolution

0.6

(d)

209

196

183

170

157

144

131

118

105

92

79

66

53

40

27

1

209

196

183

170

157

144

131

118

105

92

79

66

40

27

1

53

Upslope catchment area (km2)

14

0.5

0.6
14

Fraction of catchment area

1

0.7

0.8

(a)

(c)

Fraction of catchment area

1

14

Fraction of catchment area

1.1
1.1

Upslope catchment area (km2)

Figure 5. Comparison of distribution function of upslope catchment area from a 1000 m grid-resolution DEM to a finer grid-resolution DEM
and the distribution function of the upslope catchment area at the fine scale in the Kamishiiba catchment (210 km2 ). (a), (b), (c) and (d) are
the comparisons for 50 m, 150 m, 450 m and 600 m grid resolutions respectively

upslope catchment area from 50, 150, 450 and 600 m DEM resolutions match the distribution of downscaled
upslope catchment area from 1000 m DEM resolution to the respective DEM resolution. This shows that the
proposed method to downscale the upslope catchment area given by Equation (9) can be successfully used
to obtain higher resolution assessments of upslope catchment area at finer grid sizes by using only a coarse
DEM resolution.

DEVELOPMENT OF THE METHOD TO DOWNSCALE TOPOGRAPHIC INDEX OF TOPMODEL
As topographic index is scale dependent, this results in identified parameter values being dependent on the
DEM resolution. This makes it difficult to use model parameter values identified in a model of different
resolution. To overcome the problem, a scale-invariant model of topographic index is proposed. A resolution
factor and a scaled slope with a fractal method is introduced in the scale-invariant model of the topographic
index distribution to scale upslope catchment area per unit contour length a and slope angle of the groundsurface ˇ.
Resolution factor in topographic index
Figure 1 shows that higher frequency topographic information contained in the topographic index distribution is lost. One of the main reasons for the movement of the topographic index distribution towards higher
values while using coarse-resolution DEMs is due to the filtering of the upslope catchment area defined by afine DEM resolution, as explained previously. To obtain the lost portion of the upslope catchment area definedby a target fine-resolution DEM, we introduced scaled upslope catchment area, as explained in Equation (6), to
the topographic index, as shown in Equation (10). In the TOPMODEL framework, the conceptual derivationof contour length by a single flow-direction approach or multiple flow-direction approach overestimates thecontour length at the catchment divide portion when a coarse DEM is used. Thus, to achieve the unit contour
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1392

N. R. PRADHAN, Y. TACHIKAWA AND K. TAKARA

length at the target subgrid scale, the unit contour length at the available coarse-resolution DEM is divided
by the resolution factor Rf in Equation (10) as


Ci


Ns If

10
TI D ln   Wi tan ˇiRfwhere TI is the topographic index, Ci is the upslope catchment area of the coarse-resolution DEM, Wi is theunit contour length at the available coarse-resolution DEM (see Figure 4), Ns is the total number of subgrids
within a coarse-resolution grid (see Figure 4), and i is a location in a catchment. Resolution factor Rf in
Equation (10) is defined as
Coarse DEM Resolution
Wi
Rf D
11
D Ł
Target DEM ResolutionWi
where WŁi is the unit contour length of the target DEM resolution (see Figure 4). It is clear from Figure 4
that
Ns D Rf2
12
From Equations (10) and (12), the resolution factor is introduced in the topographic index asCi ln If
TI D ln
Wi Rf tan ˇi

13Fractal method for scaled steepest slope
Figure 6 shows the frequency distribution of the steepest slope at four different DEM resolutions in the
Kamishiiba catchment (210 km2 ) in Japan without taking into account the scale effect. Figure 6 clearly showsthat slopes derived from DEMs vary with the spatial resolution, with values decreasing at larger pixel sizes.
The underestimation of slopes when using the coarse-resolution DEMs can seriously affect the accuracy ofhydrological and geomorphological models (Zhang et al., 1999). To scale the local slope, we followed thefractal theory in topography and slope proposed by Klinkenberg and Goodchild (1992) and Zhang et al. (1999)
and developed a modified fractal method for steepest slope.
Fractal method for average slope estimation proposed by Zhang et al. (1999). The variogram technique
is used to calculate the fractal dimension in a region when the logarithm of the distance between samples

Fraction of catchment area

1.2
1
0.8
From 50m DEM resolution

0.6

From 150m DEM resolution

0.4

From 450m DEM resolution
From 600m DEM resolution

0.2
From 1000m DEM resolution

2.3

2.2

2.1

2

1.8

1.7

1.6

1.5

1.3

1.2

1.1

1

0.8

0.7

0.6

0.5

0.3

0.2

0.1

0

0
Slope value (m/m)

Figure 6. Comparison of steepest slope distribution function obtained from different DEM resolutions in the Kamishiiba catchment (210 km2 )Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

DOWNSCALING METHOD OF TOPOGRAPHIC INDEX DISTRIBUTION

1393

is regressed against the logarithm of the mean squared difference in the elevations for that distance. The
variogram equation used by Klinkenberg and Goodchild (1992) to calculate the fractal dimension of topographyis14Zp  Zq 2 D kd42D
pq
where Zp and Zq are the elevations at points p and q, dpq is the distance between p and q, k is a constant,and D is the fractal dimension. The topographic fractal properties of Equation (14) can be used to scale slope
as follows:
Zp  Zq
D ˛d1D15
pq
dpq
where ˛ D šk 0Ð5 is a constant. Because the left part of Equation (15) represents the surface slope, it can beassumed that the slope value  is associated with its corresponding grid size d by the equation
 D ˛d1D16

This implies that slope will be a function of the measurement scale if the topography is unifractal in a
specified range of measurement scale (Zhang et al., 1999). However, it is impossible to predict the spatial
patterns of slopes due to the single value of the fractal dimension and the coefficient in the fractal slope
equation for the whole DEM. To overcome this problem, Zhang et al. (1999) proposed that the coefficient ˛and fractal dimension D of Equation (16) are mainly controlled by the standard deviation of the elevation of
the sub-regions in a DEM and derived the regressed relations between ˛ and D separately with the standard
deviation of the elevation. In deriving the regressed relation, Zhang et al. (1999) considered the smallest
sub-area (window) to be composed of 3 ð 3 pixels. Hence, elevations of nine neighbouring grids in the DEM
were taken to obtain the standard deviation of the elevation for a sub-area.
We found that the slope derived from the method by Zhang et al. (1999) tends to match only with theaverage slope within the 3 ð 3 moving-window pixels of the coarse-resolution DEM, but completely failed
to take into consideration the steepest slope defined as the direction of the maximum drop from centre pixel
to its eight nearest neighbours, known as the D8 method. Thus, we propose a modified fractal method for the
steepest slope.
Fractal method for steepest slope. A modified model for the fractal method to account for the steepest slopechange due to change in DEM resolution is described as follows:
1. Unlike the distance d in Equation (16) being represented by constant grid size, this distance d in every grid
point in a 3 ð 3 pixels moving window is provided as the steepest slope distance dsteepest . Figure 7a shows
the steepest slope distance dsteepest to be dx, dy and dx 2 C dy 2 according to the direction of steepest
descent of the slope in the X-axis, Y-axis and diagonal axis DD respectively.
2. It is found that there is not much variation in standard deviation of elevation from a high-resolution DEM
to a low-resolution DEM in the same sub-area. Fractal dimension D is related to the standard deviation ofelevation  (m) in a 3 ð 3 pixels moving window as per Zhang et al. (1999).
D D 1Ð135 89 C 0Ð084 52 ln 

17

3. The fluctuations of the coefficient ˛ values in Equation (16) over different locations were found to behigher in comparison with the D value. Unlike the method by Zhang et al. (1999), in which ˛ values
are derived from standard deviation  of the elevation in 3 ð 3 pixels moving window, we developed a
new method where the coefficient ˛ values are derived directly from the steepest slope of the available
coarse-resolution DEM, while maintaining that the steepest slope itself represents an extreme fluctuation.
The modified equation issteepest D ˛steepest d1D
18
steepest
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1394

N. R. PRADHAN, Y. TACHIKAWA AND K. TAKARA

qscaledYD

dx2 + dy221asteepest=(1 − D)
dscaledD
3

2
2
Y ∆x + ∆yD

D

8

1

4

i

X

X

X 8
∆x

3
i

7

dx
7

6

5

D

6
D

2

dy

Y

qsteepest

4 X
5∆y Y qscaled

D

D

(a)

(b)

Figure 7. Fractal method for scaled steepest slope at a location i of the 3 ð 3 pixels moving window: (a) steepest slope structure of theavailable coarse-resolution DEM; (b) scaled steepest slope structure of the target resolution DEM

In Figure 7a for example, where the steepest slope is shown in diagonal direction, ˛steepest at that locationi from Equation (18) is given bysteepest
˛steepest D19 dx 2 C dy 2 1D
where dx and dy are the grid sizes in the X-axis and Y-axis respectively and steepest is the steepest slope
using the coarse-resolution DEM.
4. While downscaling, the distance variation in the target resolution DEM is defined according to the directionof the steepest slope in the coarse-resolution DEM. Hence, in Figure 7b the downscaled steepest slope scaled
is shown to be in the same direction as that of the coarse-resolution DEM steepest slope (Figure 7a). Thus,
in Figure 7b, the downscaled steepest slope scaled is given as
scaled D ˛steepest d1D
scaled20

where dscaled D x 2 C y 2 in Figure 7b and x, y are the grid sizes of the target resolution DEM in the
X-axis and Y-axis respectively. Thus, using the coarse-grid DEM for each grid, the value of D is obtained
from Equation (17) and ˛steepest is obtained from Equation (19). Then, by using Equation (20), the downscaled
steepest slope scaled is obtained by defining the distance dscaled at target DEMs.
Figure 8a–d shows the comparison of scaled slope distribution function from the 1000 m grid-resolution
DEM to the 600 m, 450 m, 150 m and 50 m grid-resolution DEMs by using Equation (20) and the distributionfunction of the slope at those respective fine-scale DEMs. A close fit of the frequency distribution of thescaled slope from the 1000 m grid-resolution DEM to 600 m, 450 m, 150 m grid-resolution DEMs is shown in
Figure 8a–c, and Figure 8d demonstrates that Equation (20) overestimated the slope when downscaling from
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1395

(c)

Slope (m/m)

(d)

2.45

2.1

2.28

1.92

1.75

1.4

1.57

1.22

1.05

0.7

1.2
1
0.8

Downscaled from1000m DEM
resolution to 50m
DEM resolution

0.6
0.4

2.45

2.28

2.1

1.75

1.57

1.4

1.22

1.05

0.88

0.7

0.53

0

1.92

From 50m DEM
resolution

0.2
0.35

2.45

2.28

2.1

1.92

1.75

1.57

1.22

1.05

0.88

0.7

0.53

0.35

0

0

1.4

From 150m DEM
resolution

0.2

Slope (m/m)

0

0.4

0

0.18

Downscaled from 1000m
DEM resolution to 150m
DEM resolution

Fraction of catchment area

1

0.6

From 450m DEM
resolution

0.2
0.88

2.45

2.1

2.28

1.92

1.75

1.4

1.57

1.2

0.8

0.4

(b)

Slope (m/m)

0.18

Fraction of catchment area

(a)

1.22

0.7

0.88

0.53

0.35

0

0

1.05

From 600m DEM resolution

0.2

0.53

0.4

Downscaled from 1000m
DEM resolution to 450m
DEM resolution

0.6

0.35

Downscaled from 1000m DEM
resolution to 600m DEM resolution

0.6

1
0.8

0

0.8

1.2

0.18

1

Fraction of catchment area

1.2

0.18

Fraction of catchment area

DOWNSCALING METHOD OF TOPOGRAPHIC INDEX DISTRIBUTION

Slope (m/m)

Figure 8. Comparison of distribution function of scaled steepest slope from a 1000 m grid-resolution DEM to a finer grid-resolution DEM
and the frequency distribution of the steepest slope at the fine scale in the Kamishiiba catchment (210 km2 ). (a), (b), (c) and (d) are thecomparisons for 600 m, 450 m, 150 m and 50 m grid resolutions respectively

1000 m to 50 m DEM resolution. This causes the downscaled topographic index distribution from 1000 m
to 50 m DEM resolution to shift slightly below that of the original 50 m DEM resolution topographic index
distribution. In Equation (20), slope is a function of the measurement scale by assuming that topography
is unifractal in a specified range of measurement scale. This unifractal concept can break down at veryfine scales (Andrle and Abrahams, 1989; Klinkenberg and Goodchild, 1992). This possible break in theunifractal condition at fine scales and its possible solution method is for further research work. In this study
we assume that the lower applicable limit of Equation (20) is 150 m when downscaling the slope distribution
from a 1000 m grid resolution to lower grid-resolution DEMs. Thus, when downscaling slope distribution
from 1000 m to 50 m DEM resolution, the slope is downscaled to 150 m, and the distribution is used fordownscaling to 50 m.
Scaled topographic index distribution
By combining Equations (13) and (20), the method to downscale the topographic index is defined asCi ln If
21
TIscaled D ln
Wi Rf scaled
The scaled topographic index TIscaled includes the resolution factor to account for the effect of scale on
upslope catchment area per unit contour length and the fractal method for scaled steepest slope as an approachto account for the effect of scale on slope.
COUPLING SCALE-INVARIANT TOPOGRAPHIC INDEX DISTRIBUTION IN THE TOPMODEL
FRAMEWORK
Total runoff is calculated as the sum of two flow components: saturation excess overland flow from variable
contributing areas (Dunne and Black, 1970) and subsurface flow from the saturated zone of the soil, as shown
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1396

N. R. PRADHAN, Y. TACHIKAWA AND K. TAKARA

in Figure 9. Areas with Si, t  0 (in Equation (2)) are contributing areas for saturation excess overland flow.
Following Equations (2), (3) and (5), the dependency of S(i,t) on other variables is shown as
Si, t D fSt, T0 , , m, TI

22

Again, from Equation (4), the dependency of subsurface flow from the saturated zone (SZ [L] in Figure 9),
Qb t is written as23Qb t D gSt, T0 , , m
In Equations (22) and (23), the only observable independent variable accounting for heterogeneity is
the topographic index TI. The parameters T0 , , m and St in the catchment or subcatchment are directly
influenced by the topographic index value, which changes with the resolution of the DEM used. This makes
recalibration in the model compulsory when the scale of application of the model differs from the scale atwhich the model parameters are identified. Thus, in the TOPMODEL concept, to make the topographic index
value scale invariant, Equation (1) has been replaced by Equation (21). The scale invariant-function definedby Equation (21) is S0 based on scale laws and does not add any extra parameter burden when coupledwith TOPMODEL.
Root zone store RZ [L] in Figure 9 for each topographic index value is depleted only by evapotranspiration,and that water is added to the unsaturated zone drainage storage UZ [L] only after the root zone reaches field
capacity or maximum root zone storage RZmax [L]. The drainage from the unsaturated zone is assumed to be
essentially vertical and drainage flux per unit area qv [LT1 ] is calculated for each topographic index class(Beven, 1986) asqv i, tt D minfT0 eSi,t/m t, UZ i, tg
24
where t [T] is calculation time step. The initial condition for average saturation deficit is derived from
Equation (4) taking Qb 0 as the initial observed discharge (Beven, 1987). S in the successive time step iscalculated bySt C 1 D St  Qv tt C Qb tt
25

Zero saturation deficit area

Downscaling of topographic index distribution

S(i,t)

RZmax

RZ
UZ
SZ

Routed overland flow Qof contributed from saturated area
Baseflow, Qb
Total simulated discharge, Qs

Figure 9. Coupling the method to downscale topographic index distribution in the TOPMODEL framework

Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1397

DOWNSCALING METHOD OF TOPOGRAPHIC INDEX DISTRIBUTION

where Qv t, the total input to groundwater from the unsaturated zone, is the sum of qv t over all grids inthe catchment, and Qb t is the groundwater discharge to the stream. The Muskingum-Cunge routing method
is used for calculating hill slope channel routing (Cunge, 1969).

RESULTS AND DISCUSSION
Application of the method to downscale the topographic index distribution
The method to downscale the topographic index of TOPMODEL is applied to the Kamishiiba catchment(210 km2 ) in Japan. Using Equation (21), the topographic index is then downscaled. Figure 10 shows thecomparison of the downscaled topographic index distribution function using Equation (21) from 1000 m to
50 m DEM resolution with variable If given by Equation (7) and with constant If D 1. In Figure 10, no
significant change is seen in the topographic index distribution, except for the highest topographic index
values in the lowest density range in a catchment. This is because If values from Equation (7) are also almost
equal to unity at the higher density distribution of topographic index values. This is why, for the rest of theanalysis in this discussion, If D 1 is assumed in Equation (21). This is the same as taking only the first termin Equation (13). Moreover, in the TOPMODEL framework, the conceptual derivation of contour length by
the single flow-direction approach or the multiple flow-direction approach underestimates the contour length
at the catchment outlet when a fine DEM resolution is used. Also, the introduction of If with Rf and Wi in
Equation (10) can correct this error by varying the contour length derivation from the target fine grid resolution
at the catchment divide portion to the available coarse grid resolution at the catchment outlet portion. The
introduction of If with Rf in Equation (10) is also the same as taking only the first term in the right hand sideof Equation (13).
Table II shows the scaled topographic constant  from a 1000 m grid-resolution DEM to various target DEM
resolutions by applying the downscaling method. The downscaled values of  from 1000 m grid resolution
to finer DEM resolutions in Table II are almost equal to the values of  in Table I derived from those finegrid-resolution DEMs.
Figure 11a is the topographic index distribution using a 1000 m DEM. Figure 11b–e shows the scaled
topographic index distribution obtained by using the downscaling method with the same 1000 m gridresolution DEM to the 600 m, 450 m, 150 m and 50 m grid-resolution DEMs respectively. Figure 11f is
the topographic index distribution using the 50 m DEM. Distinct differences can be seen between the spatial
distribution of the topographic index in Figure 11a and Figure 11f for the 1000 m grid-resolution DEM and the
0.3
Scaled topographic index from 1000m grid resolution
DEM to 50m grid resolution DEM with varying influence
factorDensity

0.25
0.2

Scaled topographic index from 1000m grid resolution
DEM to 50m grid resolution DEM with constant influence
factor equal to 1

0.15
0.1
0.05

35

33.4

31.7

30.1

28.4

26.8

25.1

23.5

21.8

20.2

18.5

16.9

15.2

13.6

11.9

10.3

8.6

6.95

5.3

2

3.65

0

Topographic index value [ln(m2)]
Figure 10. Effect of influence factor on scaled topographic index distribution. Comparison of density function of scaled topographic index
using Equation (21) with varying influence factor given by Equation (7) and with constant influence factor equal to unity. In both cases thedownscaling is made from a 1000 m to a 50 m grid-resolution DEM

Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1398

N. R. PRADHAN, Y. TACHIKAWA AND K. TAKARA

Table II. Spatial mean of topographic index  using downscaled topographic
index for target resolutions derived from 1000 m grid DEM in the Kamishiiba catchment
Target DEM resolution (m)Spatial mean of topographic
index  lnm2 

50
6Ð474

150
7Ð573

450
9Ð110

600
9Ð604

1000m grid resolution DEM
topographic index distribution(available DEM)
Coarse resolution DEM
Down scaled tographic index
distribution from 1000m grid
resolution DEM to 600m grid
resolution DEM

Scale Invariant Model for topographic
index of TOPMODEL

Down scaled topographic index
distribution from 1000m grid
resolution DEM to 450m grid
resolution DEM

Resolution factor in topographic index
Fractal method for scaled steepest slope

Down scaled topographic index
distribution from 1000m grid
resolution DEM to 150m grid
resolution DEM(a)Down scaled topographic index
distribution from 1000m grid
resolution DEM to 50m grid
resolution DEM

Scaled topographic index attarget DEM resolution(b)

Topographic index at
50m grid resolution DEM

LegendTopographic index range
3.1 to 4.2
4.3 to 5.3
5.4 to 6.4
6.5 to 7.5
7.6 to 8.6
8.7 to 9.7
9.8 to 10.8
> 10.8(c)(d)
Similar spatial distribution of topographic index

(e)

(f)

Figure 11. Spatial distribution of scaled topographic index applied to the Kamishiiba catchment (210 km2 ). (a) Topographic index distribution
using a 1000 m DEM resolution. (b), (c), (d) and (e) Scaled topographic index distributions obtained from a 1000 m DEM resolution to
600 m, 450 m, 150 m and 50 m DEM resolutions respectively. (f) Topographic index distribution using a 50 m DEM resolution

50 m grid resolution DEM respectively. The spatial distribution of topographic index displayed in Figure 11ematches with the existing reality displayed in Figure 11f.
Figure 12a–d shows the comparison of scaled topographic index distribution from the 1000 m gridresolution DEM to the 50 m, 150 m, 450 m and 600 m grid-resolution DEMs by using the downscaling
method and the density function of the topographic index at a different finer scale. A close fit of density
functions of scaled topographic indices from the 1000 m grid-resolution DEM to various grid-resolution
DEMs is shown in Figure 12. The proposed method pursues a framework for downscaling the topographic
index distribution from a coarse resolution to a specific fine-resolution DEM to obtain higher resolution
topographic information. Obtaining the specific fine-resolution DEM that is most suitable for a hydrologicalmodel is not within the scope of this paper. Zhang and Montgomery (1994) suggest from their analysis that,for many landscapes, a 10 m grid size presents a rational compromise between increasing resolution and data
volume required for simulating geomorphic and hydrological processes. Quinn et al. (1991) demonstrated
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1399

0.3
0.25
0.2
0.15
0.1
0.05
0

(d)

35

32.8

30.6

28.4

24

26.2

21.8

19.6

17.4

13

15.2

8.6

10.8

6.4

4.2

Topographic index from 150m grid resolution DEM

2
Density

35

32.8

30.6

28.4

26.2

24

21.8

19.6

17.4

13

15.2

10.8

8.6

6.4

2

Topographic index from 450m grid resolution
DEM

Scaled topographic index from 1000m grid
resolution DEM to 150m grid resolution DEM

Topographic index value [ln(m2)](b)

Scaled topographic index from 1000m grid
resolution DEM to 450m grid resolution DEM

Topographic index value [ln(m)2]

0.3
0.25
0.2
0.15
0.1
0.05
0

35

30.6

32.8

28.4

26.2

24

21.8

19.6

17.4

13

15.2

8.6

Topographic index value [ln(m2)]

4.2

Density

10.8

6.4

2

Topographic index from 50m grid
resolution DEM

(a)

(c)

Density

Scaled topographic index from 1000m grid
resolution DEM to 50m grid resolution
DEM

0.35
0.3
0.25
0.2
0.15
0.1
0.05
0

Scaled topographic index from 1000m grid
resolution DEM to 600m grid resolution DEM
Topographic index from 600m grid resolution
DEM

2
3.65
5.3
6.95
8.6
10.3
11.9
13.6
15.2
16.9
18.5
20.2
21.8
23.5
25.1
26.8
28.4
30.1
31.7
33.4
35

0.3
0.25
0.2
0.15
0.1
0.05
0
4.2

Density

DOWNSCALING METHOD OF TOPOGRAPHIC INDEX DISTRIBUTION

Topographic index value [ln(m2)]

Figure 12. Comparison of density function of scaled topographic index from a 1000 m grid-resolution DEM to a finer grid-resolution DEM
and the density function of the topographic index at the fine scale in the Kamishiiba catchment (210 km2 ). (a), (b), (c) and (d) are thecomparisons for 50 m, 150 m, 450 m and 600 m grid resolutions respectively

that features of much less than 50 m resolution are significant in hydrological routing. It is hoped that
this study can be further extended to develop the downscaling method to 10 m DEM resolution andfiner.
Analysis of prediction error when the parameter identification scale and model application scale aremismatched
Figure 13 shows three different rainfall events, Event(1), Event(2) and Event(3), in the Kamishiiba
catchment (210 km2 ). In Figure 13, A(1), A(2) and A(3) are the simulation results using a 50 m DEM
resolution TOPMODEL; B(1), B(2) and B(3) are the simulation results using a 1000 m DEM resolution
TOPMODEL; C(1), C(2) and C(3) are the simulation results from TOPMODEL by coupling the downscaling
method of topographic index distribution, applied at 1000 m DEM resolution and the topographic index
downscaled to 50 m DEM resolution. The labels (1), (2) and (3) denote the simulation results that belong torainfall Event(1), Event(2) and Event(3) respectively.
For all the simulation results shown in Figure 13, the effective parameters of TOPMODEL are identified by
the 50 m DEM resolution TOPMODEL, as shown in Table III. Although the effective parameters identified inTable III are calibrated by the 50 m DEM resolution TOPMODEL in Event(2) as shown for A(2) in Figure 13,
resulting in a Nash efficiency of 94%, consideration is also given to the model consistency when applying the
50 m DEM resolution TOPMODEL with the same set of parameters to other events. In Event(1) and Event(3) the Nash efficiencies are 90% and 78% respectively, as shown in Figure 13. In accordance with the workof Ambroise et al. (1996) and Güntner et al. (1999), a parameter value for m was derived from the first-order
hyperbolic function that fitted the recession curve obtained from the observed discharge of rainfall Event(1).
The value of the parameter m derived by the recession analysis is 0Ð07 m, which is found to be equal to thecalibrated parameter value of m by the 50 m DEM resolution TOPMODEL.
Comparing B(1), B(2), and B(3) with A(1), A(2), and A(3) respectively in Figure 13, it is seen that themagnitudes of the differences in TOPMODEL predictions based on different DEM resolutions are significant
when the same set of parameter values are used. The simulated discharges for B(1), B(2) and B(3) in Figure 13
are very sensitive to rainfall, showing overestimation of discharge in the rainfall duration and underestimationof discharge as soon as the rainfall stopped or diminished. Because of this inconsistency in the 1000 m DEM
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

Precipitation (m/Hour)

Discharge (m3/s)

Discharge (m3/s)

Observed

CalculatedObserved

CalculatedTime (Hour)
A(1), Nash Efficiency=90%
1400
1200
1000
800
600
400
200
0

CalculatedObserved

Time (Hour)
A(2), Nash Efficiency=94%

Observed

CalculatedObserved

CalculatedObserved

CalculatedTime (Hour)
B(2), Nash Efficiency=-45%

Time (Hour)

B(1), Nash Efficiency=50%
800
700
600
500
400
300
200
100
0

700
600
500
400
300
200
100
0

0
0.005
0.01
0.015
0.02
0.025
0.03
0.035

Event(2): 1999[9/22, 6:00~9/27, 6:00]

Time (Hour)
C(2), Nash Efficiency=89%

Time (Hour)

C(1), Nash Efficiency=92%

0

50

100

150

200

250

0

100

200

300

400

500

600

0

50

100

150

200

0
0.005
0.01
0.015
0.02
0.025
0.03
0.035
250

Event(3): 1999[6/24, 0:00~7/3, 21:00]

CalculatedObserved

CalculatedObserved

A(3), Nash Efficiency=78%

Time (Hour)

Time (Hour)

CalculatedObserved

B(3), Nash Efficiency=-180%

Time (Hour)
C(3), Nash Efficiency=71%

Figure 13. Simulation results of TOPMODEL (with and without coupling the downscaling method of topographic index) applied to the Kamishiiba catchment (210 km2 ). A(1),A(2), A(3) are the simulation results from a 50 m DEM TOPMODEL. B(1), B(2), B(3) are the simulation results from a 1000 m DEM TOPMODEL without scale-invariant
model for topographic index. C(1), C(2), C(3) are the simulation results from TOPMODEL at a 1000 m DEM resolution and topographic index downscaled to a 50 m DEM
resolution. The same set of effective parameter values, identified by a 50 m DEM resolution TOPMODEL using Event(2), is used for all the simulation results in A(1), A(3),B(1), B(2), B(3), C(1), C(2) and C(3)

1600
1400
1200
1000
800
600
400
200
0

0

500

1000

1500

2000

1600
1400
1200
1000
800
600
400
200
0

0
0.005
0.01
0.015
0.02
0.025
0.03
0.035
0.04

Event(1): 1997[9/15, 1:10~9/19, 0:00]

1
6
11
16
21
26
31
36
41
46
51
56
61
66
71
76
81
86
91

Discharge (m3/s)

1
6
11
16
21
26
31
36
41
46
51
56
61
66
71
76
81
86
91

1
6
11
16
21
26
31
36
41
46
51
56
61
66
71
76
81
86
91

1
6
11
16
21
26
31
36
41
46
51
56
61
66
71
76
81
86
91

Precipitation (m/Hour)
Discharge (m3/s)
Discharge (m3/s)
Discharge (m3/s)

1
7
13
19
25
31
37
43
49
55
61
67
73
79
85
91
97
103
109
115
1
7
13
19
25
31
37
43
49
55
61
67
73
79
85
91
97
103
109
115
1
7
13
19
25
31
37
43
49
55
61
67
73
79
85
91
97
103
109
115
1
7
13
19
25
31
37
43
49
55
61
67
73
79
85
91
97
103
109
115

Precipitation (m/Hour)
Discharge (m3/s)
Discharge (m3/s)
Discharge (m3/s)

1
12
23
34
45
56
67
78
89
100
111
122
133
144
155
166
177
188
199
210
221
232
1
12
23
34
45
56
67
78
89
100
111
122
133
144
155
166
177
188
199
210
221
232
1
12
23
34
45
56
67
78
89
100
111
122
133
144
155
166
177
188
199
210
221
232

Copyright  2006 John Wiley & Sons, Ltd.
1
12
23
34
45
56
67
78
89
100
111
122
133
144
155
166
177
188
199
210
221
232

1400
N. R. PRADHAN, Y. TACHIKAWA AND K. TAKARA

Hydrol. Process. 20, 1385– 1405 (2006)

DOWNSCALING METHOD OF TOPOGRAPHIC INDEX DISTRIBUTION

1401

Table III. Effective parameter values identified
by a 50 m DEM resolution TOPMODEL for
Event(2) used for different events and different
DEM resolution TOPMODELs in the Kamishiiba catchment
T0 m2 h1 
9Ð8

m (m)
0Ð07

RZmax (m)0Ð001resolution TOPMODEL prediction, the Nash efficiencies for B(1), B(2), and B(3) in Figure 13 decreased
drastically to 50%, 45% and 180% respectively when retaining the same parameter values identified forthe 50 m DEM resolution TOPMODEL.
Earlier, in Table I, we showed that, as the resolution of DEMs is coarsened, the spatial mean value of
topographic index  in Equation (5) increases. As  increases, the initial average saturation deficit derived
from Equation (4) decreases. The increase in  value also results in less variability of subsurface flow, as
defined in Equation (4). Thus, the predicted mean depth of the water table approaches the surface at theinitial stage, and the average saturation deficit obtained from Equation (25) remains around zero. In themodel, zero saturation deficit means that the state of the vertical soil profile is completely saturated up to
the surface. Any further rainfall after the zero local saturation deficit state contributes directly to surface
runoff. Figure 14a demonstrates that the average saturation deficit predicted by the 1000 m DEM resolution
TOPMODEL, while retaining the effective parameter values identified at the 50 m DEM resolution, shows
negative values throughout the simulation period in Event(1). This results in the complete saturation of thelocal soil moisture deficit in a large part of the catchment from the very beginning of the simulation. This is
why the overestimations in simulated discharge during rainfall duration for B(1), B(2) and B(3) are observedin Figure 13.
As we have discussed, the higher  value results in less variability of subsurface flow rate given by
Equation (4). Figure 14b shows that the subsurface flow Qb t obtained from the 1000 m DEM resolution
TOPMODEL is almost constant and relatively low during rainfall hours when compared with the other models.
The difference of subsurface flow obtained from the 1000 m DEM resolution TOPMODEL and that fromthe 50 m DEM resolution TOPMODEL or the downscaled 50 m DEM resolution TOPMODEL is significant.
This is the reason for the underestimation of the simulated discharge during hours of no rainfall for B(1),B(2) and B(3) in Figure 13.
The overestimation of the simulated discharge in rainfall duration and underestimation of the simulated
discharge in the no-rainfall duration is the reason for the increase in the variance in simulated discharge
when using a coarse-resolution DEM. Zhang and Montgomery (1994) and Wolock and Price (1994) have also
shown that increasing the coarseness of DEM data resolution tended to decrease the mean depth to watertable and increase the peak flow.
TOPMODEL coupled with the downscaling method of the topographic index distribution
TOPMODEL coupled with the downscaling method of the topographic index distribution is then applied
to solve the prediction error arising from the discrepancy in scales between the scale at which the model
parameters are identified and the scale of model applications. The method proposed enables downscaling of
the topographic index from a coarse grid-resolution DEM to a fine grid-resolution DEM. Incorporating thedownscaled topographic index in TOPMODEL, we can see that the average saturation deficit simulated is quite
similar to average saturation deficit simulated by the 50 m DEM resolution TOPMODEL in Figure 14a. As theresults of the scaled topographic index distribution and the scaled average saturation deficit produced by the1000 m grid resolution DEM shown in Figures 11 and 14a match with the topographic index distribution andaverage saturation deficit produced by the 50 m grid-resolution DEM, simulated discharge of the downscaled
TOPMODEL displayed by C(1), C(2) and C(3) in Figure 13 matches the simulated discharge of the 50 m
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1402

N. R. PRADHAN, Y. TACHIKAWA AND K. TAKARA

Average saturation deficit from TOPMODEL applied at1000m DEM
resolution and topographic index downscaled to 50m DEM resolution

0.3

Average saturation deficit from 50m DEM resolution TOPMODEL

Average saturation deficit (m)0.25
Average saturation deficit from 1000m DEM resolution TOPMODEL

0.2
0.15
0.1
0.05
0
-0.05

1

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89 93

-0.1
-0.15

(a) -0.2

Time (Hour)

Subsurface flow (m3/s)

1500

Subsurface flow from TOPMODEL
applied at 1000m DEM resolution andtopographic index downscaled to50m
DEM resolution

1200

Subsurface flow from 50m DEM
resolution TOPMODEL

900
600

Subsurface flow from 1000m DEM
resolution TOPMODEL

300
0
1

5

9 13 17 21 25 29 33 37 41 45 49 53 57 61 65 69 73 77 81 85 89 93(b)Time (Hour)

Figure 14. Analysis of average saturation deficit and subsurface flow results of TOPMODEL (with and without coupling the downscaling
method of topographic index) applied to the Kamishiiba catchment (210 km2 ) for Event(1). (a) Comparison of average saturation deficit
obtained from a 50 m DEM resolution TOPMODEL, a 1000 m DEM resolution TOPMODEL and TOPMODEL applied at a 1000 m
DEM resolution with scaled topographic index to a 50 m DEM resolution. (b) Comparison of subsurface flow obtained from a 50 m DEM
resolution TOPMODEL, a 1000 m DEM resolution TOPMODEL and TOPMODEL applied at a 1000 m DEM resolution with scaled
topographic index to a 50 m DEM resolution. In the results of (a) and (b), the applied effective parameters are identified by a 50 m DEM
resolution TOPMODEL

DEM resolution TOPMODEL displayed by A(1), A(2) and A(3) in Figure 13. Thus, the Nash efficiencies forC(1), C(2) and C(3) are 92%, 89% and 71% respectively in Figure 13.
An increase in , the mean of lna/ tan ˇ, caused by using a coarser scale DEM, would be compensated byan increase in calibrated value of T0 , producing similar observed and simulated hydrographs. The parametervalue, however, would exceed physically acceptable ranges, leading to the problem of false assumptions thatare less restrictive than otherwise thought, as calibration can often compensate for such deficiencies. Table IVshows the DEM resolution dependence of the effective parameter value of T0 in TOPMODEL. The parametervalue of T0 in Table IV is the calibrated value for each DEM resolution to get the highest Nash efficiency in
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1403

DOWNSCALING METHOD OF TOPOGRAPHIC INDEX DISTRIBUTION

Table IV. Scale dependence of effective parameter values in TOPMODEL and downscaling methodof topographic index distribution for its solution
DEM resolution (m)parameter
Calibrated effective parametric value of T0 m2 h1  forthe highest Nash efficiency [%]Nash efficiency (%) at T0 D 6Ð0 m2 h1Nash efficiency (%) with downscaled topographic index to
50 m DEM resolution and using T0 D 6Ð0 m2 h1

50

500

1000

6Ð0 [96Ð1]

90 [95Ð4]200 [95Ð6]

96Ð1
96Ð1

70Ð9
94Ð3

20Ð5
94Ð1

rainfall Event 1 in the Kamishiiba catchment. In calibrating the effective parameter value at different DEM
resolutions, Table IV shows that the T0 value has to be increased from 6Ð0 m2 h1 at a 50 m DEM resolution
to 200 m2 h1 at a 1000 m DEM resolution to obtain the highest Nash efficiency. Ibbitt andWoods (2004)
and Franchini et al. (1996) point out that such a high soil parameter value may preserve subsurface flow
properties in a desirable way, but at the same time could affect the infiltration behaviour, if the verticaland lateral hydraulic conductivities are assumed to be equal in the TOPMODEL framework. Table IV shows
the drop in Nash efficiency when the same value of T0 identified at a 50 m DEM resolution is used for
different DEM resolutions. Table IV also demonstrates that the Nash efficiency converges to a maximum
Nash efficiency at each of the DEM resolutions when the topographic index is downscaled to a 50 m DEM
resolution from the respective DEM resolutions and the same value of T0 that is identified at a 50 m DEM
resolution is used. This proves that TOPMODEL coupled with the downscaling method of topographic index
distribution defined by Equation (21) is consistent with observations (observed discharge data), although thescale of the DEM at which the parameters are identified and the scale of the DEM at which the model isapplied are dissimilar.

CONCLUSIONS
There is a long tradition in geomorphology of seeking general rules for landscape evolution so as to provide
a topographic framework for hydrological, ecologic, or other landscape-based models. Geomorphologicalfeatures obtained from DEMs are influenced by the resolution of the DEM. In this study we developed
a method to downscale the upslope catchment area to solve the effect of DEM resolution on the upslope
catchment area. Through analysis of scale laws, the concept of a resolution factor has been developed in thisstudy to account for the effect of scale on upslope catchment area per unit contour length in topographic index,
and a fractal method for scaled steepest slope has been developed as an approach to account for the effect ofscale on slopes. These are combined to develop the method to downscale the topographic index distribution.
The method for downscaling the topographic index represents a scale-invariant function for topographic
index distribution and has made possible the utilization of parameter values calibrated using high-resolution
DEMs in making predictions with TOPMODEL using coarser resolution DEMs. This study can be extended
in analysis of scaling behaviour across catchments through developing a topography-driven regionalization
model, at least for regions that are ‘homogeneous’ in some sense (Pradhan et al. 2005). Although assessmentof uncertainty has not been conducted as part of this study, it is worth noting that a large part of the uncertainty
in hydrological prediction can be reduced if the phenomenon of the uncertainty is understood (Sivapalan et al.,2003). It is hoped that the downscaling method can reduce the parametric uncertainty by reducing the feasible
parameter value range in grouping of models into functional classes according to similarity of predictions forsome variable of interest, as defined by Beven (2000b). It is also hoped that the findings of this study will
demonstrate their applicability as a tool to solve a wider range of problems, especially with respect to scaleissues, when modelling hydrological processes.
Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

1404

N. R. PRADHAN, Y. TACHIKAWA AND K. TAKARAACKNOWLEDGEMENTS

Special thanks go to Professor K. Beven at the University of Lancaster for his valuable comments on the
manuscript. We also thank an anonymous referee and Dr T. Nakaegawa for their valuable advice.

REFERENCES
Abbott MB, Bathurst JC, Cunge JA, O’Connell PE, Rasmussen J. 1986b. An introduction to the European Hydrological System–Système
Hydrologique Europeèn, SHE, 2: Structure of a physically-based, distributed modelling system. Journal of Hydrology 87: 61–77.
Ambroise B, Beven K, Freer J. 1996. Towards a generalization of the TOPMODEL concepts: topographic indices of hydrologic similarity.
Water Resources Research 32: 2135– 2145.
Andrle R, Abrahams AD. 1989. Fractal techniques and the surface roughness of talus slopes. Earth Surface Processes and Landforms 14:
197–209.
Band LE, Moore ID. 1995. Scale: landscape attributes and geographical information systems. Scale Issues in Hydrological Modelling,
Kalma JD, Sivapalan M (eds). Wiley: Chichester; 159–179.
Beven K. 1986. Runoff production and flood frequency in catchments of order n: an alternative approach. Scale Problems in Hydrology,
Gupta VK, Rodrı́guez-Iturbe I, Wood EF (eds). Reidel: Dordrecht; 107–131.
Beven K. 1987. Towards the use of catchment geomorphology in flood frequency predictions. Earth Surface Processes and Landforms 12:
69–82.
Beven K. 1995. Linking parameters across scales: subgrid parameterizations and scale dependent hydrological models. Scale Issues inHydrological Modelling, Kalma JD, Sivapalan M (eds). Wiley: Chichester; 263– 281.
Beven KJ. 2000a. Rainfall-Runoff Modelling. Wiley: New York.
Beven KJ. 2000b. Uniqueness of place and process representations in hydrological modeling. Hydrology and Earth System Sciences 4:
203–214.
Beven K. 2002. Towards an alternative blueprint for a physically based digitally simulated hydrologic response modelling system.
Hydrological Processes 16: 189– 206.
Beven KJ, Kirkby MJ. 1979. A physically based, variable contributing area model of basin hydrology. Hydrological Sciences Bulletin 24:
43–69.
Bruneau P, Gascuel-Odoux C, Robin P, Merot Ph, Beven KJ. 1995. Sensitivity to space and time resolution of a hydrological model using
digital elevation data. Hydrological Processes 9: 69–82.
Cunge JA. 1969. On the subject of a flood propagation computation method (Muskingum method). Journal of Hydraulic Research 7:
205–230.
Dunne T, Black RD. 1970. An experimental investigation of runoff production in permeable soils. Water Resources Research 6: 478–490.
Franchini M, Wendling J, Obled C, Todini E. 1996. Physical interpretation and sensitivity analysis of the TOPMODEL. Journal of Hydrology
175: 293– 338.
Freeze RA, Harlan RL. 1969. Blueprint for a physically based, digitally simulated hydrological response model. Journal of Hydrology 9:
237–258.
Goodchild MF. 2001. Models of scale and scales of modelling. Modelling Scale in Geographical Information Science, Tate NJ, Atkinson PM(eds). Wiley: Chichester; 3–10.
Güntner A, Uhlenbrook S, Seibert J, Leibundgut Ch. 1999. Multi-criterial validation of TOPMODEL in a mountainous catchment.
Hydrological Processes 13: 1603– 1620.
Ibbitt R, Woods R. 2004. Re-scaling the topographic index to improve the representation of physical processes in catchment models. Journalof Hydrology 293: 205– 218.
Iorgulescu I, Jordan J-P. 1994. Validation of TOPMODEL on a small Swiss catchment. Journal of Hydrology 159: 255–273.
Jenson SK, Domingue JQ. 1988. Extracting topographic structures from digital elevation data for geographic information systems analysis.
Photogrammetric Engineering and Remote Sensing 54: 1593– 1600.
Kirkby MJ. 1975. Hydrograph modeling strategies. Process in Physical and Human Geography, Peel R, Chisholm M, Haggett P (eds).
Heinemann: London; 69–90.
Klinkenberg B, Goodchild MF. 1992. The fractal properties of topography: a comparison of methods. Earth Surface processes and Landforms
17: 217–234.
Leopold LB, Maddock T. 1953. The hydraulic geometry of stream channels and some physiographic implications, United States, Geological
Survey, Professional Papers 252.
Mendicino G, Sole A. 1997. The information content theory for the estimation of the topographic index distribution used in TOPMODEL.
Hydrological Processes 11: 1099– 1114.
O’Loughlin EM. 1986. Prediction of surface saturation zones in natural catchments by topographic analysis. Water Resources Research 22:
794–804.
Pradhan NR, Tachikawa Y, Takara K. 2004. A scale invariance model for spatial downscaling of topographic index in TOPMODEL. Annual
Journal of Hydraulic Engineering, JSCE 48: 109–114.
Pradhan NR, Tachikawa Y, Takara K. 2005. Development of a transferable hydrologic modeling in TOPMODEL framework across scaleand region. Annual Journal of Hydraulic Engineering, JSCE 49: 217–222.
Quinn P, Beven K, Chevallier P, Planchon O. 1991. The prediction of hillslope flow paths for distributed hydrological modelling using
digital terrain models. Hydrological Processes 5: 59–79.

Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)

DOWNSCALING METHOD OF TOPOGRAPHIC INDEX DISTRIBUTION

1405

Quinn P, Beven K, Lamb R. 1995. The lna/ tan ˇ index: how to calculate it and how to use it within the TOPMODEL framework.
Hydrological Processes 9: 161–182.
Rodriguez-Iturbe I, Rinaldo A. 1997. Fractal River Basins: Chance and Self-organization. Cambridge University Press.
Rogers CCM, Beven KJ, Morris EM, Anderson MG. 1985. Sensitivity analysis, calibration and predictive uncertainty of institute ofhydrology distributed model. Journal of Hydrology 81: 179–187.
Saulnier G-M, Obled Ch, Beven KJ. 1997. Analytical compensation between DTM grid resolution and effective values of saturated hydraulicconductivity within the TOPMODEL framework. Hydrological Processes 11: 1331– 1346.
Sivapalan M, Takeuchi K, Franks SW, Gupta VK, Karambiri H, Lakshmi V, Liang X, McDonnell JJ, Mendiondo EM, O’Connell PE,
Oki T, Pomeroy JW, Schertzer D, Uhlenbrook S, Zehe E. 2003. IAHS Decade on Predictions in Ungauged Basins (PUB), 2003– 2012:
shaping an exciting future for the hydrological sciences. Hydrological Sciences Journal–Journal des Sciences Hydrologiques 48: 857–880.
Wolock DM, Price CV. 1994. Effects of digital elevation model map scale and data resolution on a topography based watershed model.
Water Resources Research 30: 3041– 3052.
Woods EF, Beven K, Sivapalan M, Band L. 1988. Effects of spatial variability and scale with implication to hydrology modeling. Journalof Hydrology 102: 29–47.
Woods RA, Sivapalan M, Robinson JS. 1997. Modeling the spatial variability of surface runoff using a topographic index. Water Resources
Research 33: 1061– 1073.
Zhang W, Montgomery DR. 1994. Digital elevation model grid size, landscape representation, and hydrologic simulations. Water Resources
Research 30: 1019– 1028.
Zhang X, Drake NA, Wainwright J, Mulligan M. 1999. Comparison of slope estimates from low resolution DEMs: scaling issues and afractal method for their solution. Earth Surface Processes and Landforms 24: 763– 779.

Copyright  2006 John Wiley & Sons, Ltd.

Hydrol. Process. 20, 1385– 1405 (2006)