Quality assessment from a hydrological perspective of a digital elevation model derived from WorldView-2 remote sensing data

Abstract

A digital elevation model (DEM) derived from a stereo pair of WorldView-2 (WV-2) images was assessed against ground-truth GPS point datasets. Two assessment methods were used: (a) vertical accuracy assessment and (b) hydrological assessment of surface runoff variables. Three agricultural plots with different topographic slopes were selected to perform a vertical accuracy assessment, followed by a comparative assessment of a set of hydrological variables. The results show an overall vertical accuracy of 0.45 m, confirming the potential of WV-2 stereo images to extract elevation information at high spatial resolution. Concerning plot-scale micro-topographic features, the WV-2 DEM performed better on the plot with rolling slopes (5–10%), extracting variables such as the total length and drainage area of flow paths with relative errors lower than 20%. However, some limitations were detected in the extraction of variables such as terrain slope, drainage points of flow paths and terrain depressions in areas of flatter slopes (<5%).

Editor Z.W. Kundzewicz; Associate editor D. Gerten

Résumé

Nous avons évalué la qualité, d’un point de vue hydrologique, d’un modèle numérique de terrain (MNT) issu d’un couple d’images stéréoscopiques WorldView-2 (WV-2) en le comparant à des mesures de terrain obtenues par GPS. Deux méthodes d’évaluation ont été utilisées: (a) l’évaluation de la précision verticale et (b) l’évaluation hydrologique de variables de ruissellement. Trois parcelles agricoles de pentes différentes ont été choisies pour effectuer une évaluation de la précision verticale, suivie d’une évaluation comparative d’un ensemble de variables hydrologiques. Les résultats montrent une précision verticale globale de 0,45 m, ce qui confirme le potentiel des images stéréo WV-2 à extraire des informations altimétriques à haute résolution spatiale. En ce qui concerne les caractéristiques micro-topographiques à l’échelle de la parcelle, le modèle de terrain WV-2 obtient de meilleurs résultats sur les terrains aux pentes prononcées (5–10%), permettant l’extraction de variables telles que la longueur et la superficie totale de drainage des chenaux d’écoulement avec des erreurs relatives inférieures à 20%. L’extraction de caractéristiques telles que la pente, les exutoires des chenaux d’écoulement et les dépressions, a cependant montré qu’elle avait des limites dans les zones de faible pente (<5%).

Key words:

• surface hydrology
• surface runoff
• DEM
• plot scale hydrology
• WorldView-2
• ArcHydro

Mots clefs:

• hydrologie de surface
• ruissellement de surface
• MNT
• hydrologie à l’échelle de la parcelle
• WorldView-2
• ArcHydro

1 INTRODUCTION

Gridded digital elevation models (DEM) are used in many disciplines due to their ability to represent topographic characteristics required in landscape, hydrological and erosional modelling (Li et al. 2005, Longley 2005). The uses of DEMs in hydrological analysis include: (a) automatic extraction of terrain properties and drainage networks, (b) determination of flow direction and accumulation, (c) delineation of catchment areas, and (d) extraction of watercourses (Wechsler 2007, Nigel et al. 2010a).

Several studies have been done at the basin scale to analyse the capabilities and limitations of DEMs in hydrological modelling (Callow et al. 2007, Hosseinzadeh 2011, El Hage et al. 2012, Rezak et al. 2012). However, at more detailed (plot) scales there is a lack of knowledge to support the operational and practical use of DEMs (Lyon et al. 2003), particularly for agricultural plots, where it is anticipated that high-resolution DEMs can provide farmers with accurate topographic information to help them develop better, more productive, environmentally conscious practices to increase yield and decrease environmental impacts (van der Velde et al. 2010, Ge et al. 2011).

Lowland watersheds—generally characterized by fertile soils, water availability and high accessibility—are commonly exploited for urbanization or agricultural practices. This gives rise to environmental problems such as land loss, and water and land pollution resulting from these land-use changes and exploitation (e.g. high nutrient inputs, modifications of discharge water flow paths). Hydrological modelling is an important task in such watersheds (van der Velde et al. 2010) and thus it is essential to evaluate topographic representations of watersheds using DEMs to identify their capabilities and limitations.

Topography defines water movement across watersheds and, thus, it is a major factor in modelling hydrological processes. Among the DEM-extracted topographic variables most used in hydrological modelling are slope, flow length, drainage area, watershed area and combined variables such as the topographic index (Wu et al. 2008). Several studies have been done (Wolock and Price 1994, Zhang and Montgomery 1994, Chaubey et al. 2005) that confirm that at lower spatial resolution, slopes decreases and drainage area increases, thus reducing the topographic index values. The latter is a key factor in many distributed hydrological models because it is used for estimating the spatial distribution of groundwater levels, surface saturation and runoff discharge (Mukherjee et al. 2013). Accordingly, increasing coarseness of the DEM tends to decrease the mean depth of water table and to increase daily and total flow. Moreover, surface saturation and surface runoff peak discharges increase when the DEM spatial resolution is decreased. Using TOPMODEL, a topographic-based tool to simulate hydrological flux in a watershed, flow length and watershed area did not show a systematic pattern against DEM spatial resolution, a fact confirmed by Charrier and Li (2012) when analysing LiDAR-derived DEMs. Nevertheless, Chaubey et al. (2005) detected a pattern when evaluating outcomes from the Soil and Water Assessment Tool (SWAT), where watershed and sub-watershed delineation at coarser spatial resolutions presented smaller areas compared to ground truth. The impact of flow lengths in hydrological modelling has also been evaluated and results show there might be a time lag between precipitation and flow peak discharge, altering the predicted hydrographs (Wu et al. 2008). Indeed, water stream lengths tend to gradually decrease with a coarser DEM spatial resolution, but in the prediction of floodplain boundaries, flooded surface tends to increase (Charrier and Li 2012). Flow length is extremely important at plot scale because flow paths determine the hydrological connectivity of saturated areas in agricultural areas; thus, incorrect flow path extraction and micro-watershed delineation could lead to erroneous estimations (Shore et al. 2013).

Hydrological variables derived from a DEM are highly sensitive to the accuracy of the DEM. Consequently, high-precision DEMs are required in hydrological modelling to produce unbiased predictions and outcomes (Vaze et al. 2010, Xiao et al. 2010). Furthermore, accurate elevation information is an important component of hydrological modelling at high spatial resolution. DEM-extracted variables, such as water-flow direction and accumulation, micro-basins, flow paths, drainage points, depression zones and drainage area, need to be derived from highly accurate elevation models (preferably also at high resolution) to ensure good quality outcomes (Erskine et al. 2007).

Currently, several remote sensing techniques are used to generate elevation models. These techniques, most of which have been developed and improved over the last 30 years (Toutin and Cheng 2002), range from standard and GPS topography, aero- and space-photogrammetry, and radar interferometry to airborne LiDAR. Some drawbacks of these technologies are cost, processing time and lack of diversified spectral information (with the notable exception of space photogrammetry) (Hobi and Ginzler 2012).

Very high spatial resolution (VHSR) stereoscopic satellite images, imagery with pixel sizes of 1 metre or less, have been commercially available since the launch of the IKONOS satellite in 1999. VHSR satellite stereo images are being used as a viable and fast alternative for producing elevation models (Jacobsen 2003). The WorldView-2 satellite, launched in October 2009, is equipped with panchromatic and multispectral sensors and is capable of capturing stereoscopic images. These stereo pairs can be used in automatic elevation extraction processes to generate DEMs covering the areas of the images. In addition, the satellite’s revisit time of 1–4 days makes it an important data source for monitoring spatial phenomena over time. WorldView-2 is the first satellite with one panchromatic band and eight multispectral bands at spatial resolutions of 0.5 m and 2 m, respectively. The panchromatic sensor can capture stereo image pairs that are used to generate 1-m pixel-size DEMs. Previously, such high-resolution stereo-based DEMs were only achievable with aerial photogrammetry and land surveying (Büyüksalih et al. 2012).

The aim of this research is to evaluate, at the agricultural plot scale, the capabilities and limitations of a WorldView-2 digital elevation model (WV-2 DEM) for representing elevation and extracting hydrological variables commonly used in environmental modelling. The WV-2 DEM will be evaluated against high-density and high-precision field-captured GPS point datasets based on an assessment of (a) vertical accuracy and (b) extracted hydrological variables describing spatial surface runoff.

2 METHODS

2.1 Study area

The La Chevrotière River basin is located in the Portneuf Municipality, part of the Capitale-Nationale Administrative Region, in the province of Quebec, Canada (Fig. 1). The basin lies on the north shore of the St. Lawrence River and has an area of 107 km². It is located in a low-lying region that is mainly agricultural, characterized by pasture, forest patches and small urban settlements. The combined population of Saint-Marc-des-Carrières and Saint-Gilbert, the two main urban centres, is around 3100 inhabitants. The La Chevrotière River is the union of two main branches and several small tributaries. It covers a distance of 29 km in a flat-to-undulating landscape before reaching the St. Lawrence. The highest point in the river basin has an altitude of 150 m a.s.l. The mean slope of the basin is around 5%, but in riverine zones, slopes can be greater than 8% (CAPSA 2011).

Fig. 1 The study area, the La Chevrotière River basin (107 km²), located on the north shore of the St. Lawrence River in the province of Quebec, Canada.

Like most river basins in agricultural regions, the La Chevrotière basin faces a variety of environmental pressures, including nonpoint-source pollution, erosion and land loss (CAPSA 2011). Many of these environmental problems are being addressed from a landscape perspective, thereby creating an operational gap between land managers and planners—who work at regional scales—and the occupants of the basin—who work at local scales and who typically lack technical information to understand and participate in the territorial planning process of their watershed. Tools to foresee potential problems at smaller analysis scales (e.g. agricultural plots) would give landscape planners a way to evaluate the potential hydrological problems faced locally (on plots) by farmers. Such tools would also give farmers and other watershed inhabitants information that could help them prevent environmental problems and improve their income through better practices. Provision and use of such tools would represent a viable, empowering management practice that could help to improve water quality and reduce environmental impacts (Lyon et al. 2003, van der Velde et al. 2010).

2.2 Data description

This research is based on two datasets: (a) high-precision and high-density GPS point datasets from field surveys (referred to in tables and figures as GPS) and (b) a 1-m gridded DEM with a vertical accuracy of 1 m extracted from a WorldView-2 stereo pair (referred to in tables and figures as WV-2).

2.2.1 Ground-truth GPS survey

During summer 2011 we carried out ground-truth GPS surveys in three agricultural plots with different topographic characteristics (Fig. 2). Plot 1 has a nearly level longitudinal slope of 0.3% (the longitudinal slope covers the length of the plot and is indicated by the arrow in the plot in Fig. 2). The elevation range of the plot is 2.3 m and its maximum slope is 25%. Plot 2 has a gentle longitudinal slope (1.0%), maximum slopes of ~27% and an elevation range of 4.6 m. Plot 3 has a nearly rolling longitudinal slope (5.2%) and an elevation range of 17.4 m, with maximum slopes of 63%.

Fig. 2 The three agricultural plots used for ground-truth GPS-RTK survey. The arrow in each plot represents the longitudinal slope used for describing the topography of the plot.

For this project, we mounted a GPS receiver on an all-terrain vehicle and, using a differential real-time kinematic (RTK) method, collected 131 519 xyz coordinates from a well-distributed grid of points, evenly spaced on the ground approximately every metre. The resulting dataset had an average density of 1.5 GPS points per square metre and a vertical accuracy of 10 cm (Table 1). As discussed in Höhle and Höhle (2009), the quantity and accuracy of the control points surveyed must be at least three times as great as those of the DEM evaluated. In the present work, this requirement was surpassed to ensure a correct comparative statistical analysis.

Table 1 Ground-truth dataset statistics (SD: standard deviation).

Elevation (m)	Plot
	1	2	3
Minimum	66.77	83.89	58.66
Maximum	69.05	88.44	76.12
Range	2.28	4.56	17.44
Mean	68.08	87.29	71.54
Median	68.11	87.85	73.72
SD	0.35	1.15	4.95
Longitudinal slope (%)	0.3	1.0	5.2
GPS points collected	40 603	23 603	67 313

2.2.2 Satellite DEM

The satellite stereo pair was taken by the WV-2 VHSR satellite on 14 April 2012 (15:58 and 15:59 h GMT). The stereo pair was in UTM projection and WGS84 datum.

In order to geo-rectify the stereo image pair, we provided the image supplier (in the present case, Effigis Inc.) with five ground control points. The ground control points were selected from well-distributed and easily recognizable image features, such as utility poles and the intersections of paved streets. The image supplier then took the stereo image pairs and control points and generated an ortho- and geo-rectified DEM using the BAE Systems NGATE module for stereo pairs (DeVenecia et al. 2007). That module applies image correlation algorithms to automatically extract elevation. During the process, image features, such as roads, tracks, rivers, railways and freeways, are additionally extracted using image interpretation and integrated into the DEM. The output DEM was in the local coordinates of MTM7 projection and NAD83 datum. Pixel spacing and vertical accuracy were both 1 m.

The two datasets (GPS points and WV-2 DEM) were then compared using two complementary assessment techniques. The first analysed the ability of the NGATE module to extract accurate elevation values in the form of a gridded model and the second performed a comparative evaluation of the surface hydrological variables extracted from the satellite DEM and the GPS-based DEM over the three sample plots.

The three plots used in the GPS survey were untilled and bare of plant life, with topographic features clearly visible, during both the field survey (summer 2011) and subsequent satellite image capture (spring 2012). However, there are small concave areas on the satellite image—mainly in Plot 3—where terrain characteristics were not extracted properly (during DEM generation) because of the presence of snow during image acquisition (some comments on this issue are included in the Discussion, Section 3).

2.3 Vertical accuracy assessment

For the elevation accuracy assessment, we used widely applied error metrics such as mean absolute error (MAE) and root mean square error (RMSE) (Eckert et al. 2005). Additionally, and to make this study comparable to other research and to commercially available DEM accuracy reports, we also calculated metrics such as 90th percentile linear error (LE90) and vertical accuracy at 95% confidence level (NDEP 2004). To accurately calculate absolute vertical accuracy, each sample plot was analysed individually. For this, at the location of each ground-truth GPS point, the corresponding elevation provided by the WV-2 DEM was extracted. For the vertical accuracy assessment, the GPS points were compared directly to the WV2-DEM data. We used this approach rather than a GPS-DEM versus WV2-DEM comparison to avoid errors introduced by interpolation algorithms.

2.4 Hydrological assessment

It well known that distributed hydrological models make use of gridded DEMs to extract topographic variables in order to simulate several hydrological phenomena in a watershed. Those variables are used to simulate and predict hydrological processes such as surface runoff, or they are used as input datasets for further analyses, for example, in modelling sediment/nutrient/pollutant transportation, in erosion and natural hazard modelling, or assessing biological diversity (Moore et al. 1991).

The hydrological assessment method used in this work relies on surface runoff variables derived from a grid DEM. Therefore, for the hydrological assessment, the GPS points were interpolated to generate a DEM for each plot. The DEMs were interpolated at 1-m resolution using the Inverse Distance Weighted (IDW) algorithm. This interpolation algorithm performs well when the density and precision of the source elevation data are high (Chaplot et al. 2006, Aguilar et al. 2010, Guo et al. 2010, Xiao et al. 2010). Our GPS dataset met both of these conditions, with an average point density of 1.5 points/m² and a vertical accuracy of 10 cm.

One quarter of the GPS points were randomly selected for validation purposes (i.e. not used in the interpolation process). As expected, small interpolation errors (0.03 m) were found in the three plots. However, such errors are more sensitive to the DEM’s sources than to interpolation algorithms and uncertainties (Erskine et al. 2007, Li et al. 2012). We concluded, therefore, that the IDW algorithm was an appropriate choice, offering a good balance between accuracy and computing speed.

The hydrological variables retained for assessment were (a) terrain depressions, (b) terrain slopes, (c) and flow paths. These variables and their quality indicators were chosen to ensure an easily reproducible hydrological quality assessment using available software and the simplest variables (Table 2). In the present case, the quality indicators were extracted using the software ArcGIS® and its ArcHydro® module, which is extensively used in hydrological analyses as it is effective at extracting hydrological parameters for further modelling (Maidment 2003, Nigel et al. 2010b, Zhang et al. 2010).

Table 2 Hydrological variables.

Variable	Quality indicators	Output information
Terrain depressions	Total area	Vector polygons
	Non-overlapping area
Slope	Longitudinal slope	Raster layers
	Mean slope
	Maximum slope
Flow paths	Mean length	Vector points and polylines
	Maximum length
	Total length
	Drainage outlets
	Drainage area

Pits or depressions are low points completely surrounded by higher elevation points. These features pose the main difficulty in determining hydrological flow directions because the depressions must be filled before flow can continue downstream (Jenson and Domingue 1988). Thus, a DEM that can correctly identify and represent depressions or pits can provide important information at the agricultural plot scale; such information is directly related to nutrient control and crop yields and has wide potential for use in precision agriculture (Lyon et al. 2003). Post-processing a DEM for hydrological and erosion modelling produces a ‘filled’ DEM from which pits are removed. To identify the location of depressions, we extracted terrain depression polygons by subtracting the original DEM from this ‘filled’ version. In order to evaluate the performance of the model using this technique, two quality indicators were used: total area of depressions; and non-overlapping areas of depressions. The first expresses the quantity of land surface detected as a depression, while the second expresses the spatial coincidence of these depressions over the surface.

After terrain depressions, terrain slopes were extracted. In a gridded DEM, slope is a local function and it represents the maximum rate of elevation change of a plane formed by a 3 × 3-cell neighbourhood around the centre pixel (equation (1)). The slope value of this plane was calculated using the average maximum technique (Burrough 1986) over the distance between the target cell and its surrounding eight pixels and it is expressed in decimal degrees as:

(1)

where dz/dx is the rate of elevation change in the east–west direction and dz/dy is the rate of elevation change in the north–south direction of the plane. To evaluate the performance of this variable, three quality indicators were used: the longitudinal slope (Fig. 2), the mean slope and the maximum slope.

Finally, flow paths were extracted and to evaluate their performance five quality indicators were calculated: mean length, maximum length, total length, drainage outlets and drainage area. The first three indicators are directly derived from the geometric properties of flow paths, while the other two—drainage outlets and drainage area—are extracted indirectly. Drainage outlets give precise information of the position at which water flows out of a catchment area, while drainage area corresponds to the total area drained towards the drainage point; the latter is commonly used as a proxy of surface runoff intensity. To calculate these five indicators, intermediate datasets such as flow direction and flow accumulation were extracted. The ArcHydro software uses a pre-filling algorithm to calculate flow direction (Jenson and Domingue 1988), which ensures a continuous water flow and enables further hydrological modelling. The flow direction algorithm used by the software, called D8, uses the direction of the steepest descent from the central pixel of a 3 × 3 moving window to one of its eight neighbours to determine the direction of water flow in every cell of the elevation model. Subsequently, flow paths were extracted using an accumulation threshold of 100 cells (i.e. 100 m²), which was selected after several tests were performed to determine the number of cells needed to ensure good, continuous representation of the surface runoff flow paths. Due to the importance of knowing the position and amount of water flow towards natural water streams or artificial pathways, an individual analysis was made for each side of the sample plots (i.e. northeast, northwest, etc.).

Once the hydrological variables were calculated, and to ensure easy comparability between the variables and across the three sample plots, we used the relative error (e), expressed in percentage terms (equation (2))), as a relative measure of change, as follows:

(2)

where value corresponds to the measured value from the satellite DEM and value_ref is the ground-truth value extracted from the GPS-derived elevation model. This error metric is frequently used when absolute values are expressed in different units or different magnitude scales. It must be noted that neither of the two elevation models, GPS and WV-2, were post-treated with drainage enforcement procedures in order to impose linear or polygonal entities, particularly to enforce flow direction.

3 RESULTS AND DISCUSSION

3.1 Vertical accuracy

First, we analysed the elevation residuals using the three GPS point datasets. This enabled us to identify the presence of a systematic elevation error (bias), with a global mean value of 0.78 m, revealing a trend of the satellite DEM to systematically overestimate the elevation in bare, flat terrain (Fig. 3). There are two error categories in elevation models: random and systematic. The latter category originates mostly in the production chain, which is composed of processes that make use of image matching, elevation extraction and gridding algorithms (Wise 2000). Similar elevation biases have been found in other DEMs derived from imagery from WorldView-2 and other satellites (Reinartz et al. 2010, Büyüksalih et al. 2012, Li et al. 2012). As these examples show, systematic errors are commonly present in satellite-derived elevation products. It is preferable that such errors be detected and eliminated before further processing, particularly in studies requiring absolute DEM values or operating at very high spatial resolution.

Fig. 3 Histograms of elevation residuals of the sample plots.

Table 3 shows the accuracy metrics calculated using the biased elevation residuals. We found a global RMSE of 0.81 m, vertical accuracy (VA) of 1.60 m, MAE of 0.78 m, global maximum AE of 2.36 m and a vertical linear error (LE90) at the 90th percentile of 1.05 m.

Table 3 Biased vertical accuracy metrics (m) of the WorldView-2 DEM. RMSE: root mean squared error; VA: vertical accuracy; MAE: mean absolute error; Maximum AE: maximum absolute error; LE90: vertical linear error.

Plot	Bias	RMSE	VA	MAE	Maximum AE	LE90
1	0.80	0.83	1.63	0.80	1.67	1.08
2	0.67	0.69	1.34	0.67	1.53	0.87
3	0.81	0.84	1.66	0.81	2.36	1.08
Global	0.78	0.81	1.60	0.78	2.36	1.05

Before further analysis, we removed the 0.78-m bias in order to obtain adjusted error metrics of the elevation model (Vaaja et al. 2011). Once the bias was removed (Table 4), we found a global RMSE of 0.23 m, VA of 0.45 m, MAE of 0.17 m, global maximum AE of 1.55 m and LE90 at the 90th percentile of 0.36 m.

In brief, the WV-2 DEM achieved a vertical accuracy of 0.45 m, consistent with the quality standards of the supplier. The WV-2 DEM vertical errors obtained in this work are similar to those found in other studies using WorldView-2 images, particularly those focused on analysing DEM performance on bare land (Mitchel and MacNabb 2010, Hobi and Ginzler 2012). This confirms both the capabilities of the satellite-derived elevation model and the correctness of the present methodology and data used.

Table 4 Vertical accuracy metrics (m) of the WorldView-2 DEM.

Plot	RMSE	Vertical accuracy	MAE	Maximum AE	LE90
1	0.22	0.42	0.17	0.88	0.36
2	0.17	0.33	0.13	0.87	0.56
3	0.25	0.49	0.18	1.55	0.38
Global	0.23	0.45	0.17	1.55	0.36

In addition to our analysis of vertical accuracy metrics, we made a visual comparison of transversal elevation profiles of the sample plots to evaluate the capabilities of the satellite DEM to detect micro-topographic features (Fig. 4). In the GPS-based DEM (solid line), we can clearly observe the micro-relief defined by the longitudinal water flow paths (valleys) and the cultivable areas (crests). Conversely, the WV-2 DEM (dotted line) shows mostly a smoothing effect in elevation, producing a flattening in all the profiles. The NGATE elevation extraction algorithm, which is a hybrid algorithm that uses a variably sized search window to apply edge and area matching correlation algorithms to every pixel of the image (DeVenecia et al. 2007), was unable to properly depict all of the micro-topographic features. The magnitude of the resultant flattening is evident in the maximum absolute errors recorded in the sample plots (0.88, 0.87 and 1.55 m for plots 1, 2 and 3, respectively).

Fig. 4 Transverse elevation profiles of the sample plots.

It is important to note that at the time the stereo pair was acquired (during the second week of April 2012), there were still snow patches in the study area. The snow was present mostly in depressions and drainage channels and thus certainly influenced the extraction of elevation information. Figure 5 shows the snow present in Plot 3 (the white patches on the ground). The hatched regions (polygons) have LE90 errors greater than 0.36 m.

Fig. 5 Spatial distribution of snow patches and vertical errors on Plot 3.

Furthermore, although the vertical accuracy assessment provides an indication of the DEM’s accuracy, this metric cannot be taken to express the overall performance of the WV-2 DEM in the study area. Our ground-truth check was performed only on agricultural, untilled plots. As Reinartz et al. (2010) has suggested, further research using different land cover types and probably different methodologies is needed in order to fully assess the ability of the WorldView-2 satellite to depict elevation in complex landscapes, including forested and urbanized areas. As demonstrated by the metrics calculated for the present study, the vertical precision of the elevation data in the satellite-derived DEM is less than 1 m, which may be accurate enough for environmental modelling at regional scales. However, at local (plot) scales, the inability of the DEM to pick up micro-topographic details could lead to erroneous results. To overcome this problem, DEM reconditioning methods may need to be applied. For example, AGREE, an algorithm incorporated into ArcHydro, modifies DEMs by imposing vector features such as channels and streams to ensure better surface hydrological analysis, especially in the calculation of flow direction and flow accumulation (Hellweger 1997).

3.2 Hydrological variables

Table 5 summarizes the comparison statistics of the three surface hydrological variables: terrain depressions, terrain slope and runoff flow paths. By analysing the relative errors of each quality indicator, we were able to assess the performance of the satellite elevation model at extracting these variables. We selected the relative error e (equation (2))) and its mean absolute value |e| to express the change of each quality indicator in each sample plot; this approach makes it easy to compare the variables to one another. All the error metrics are expressed in percentage terms to ease interpretation.

Table 5 Comparison statistics (terrain depressions, slope, flow paths) between the WV-2 DEM and the GPS-derived DEM for the three plots. The longitudinal slope covers the length of the plot as shown by the arrow in Fig. 2.

Hydrological variables	Plot 1			Plot 2			Plot 3			|e|
	GPS	WV-2	e	GPS	WV-2	e	GPS	WV-2	e
Depressions
Total area (m^2)	5117	9514	86%	1486	2206	48%	838	1513	81%	72%
Non-overlapping area (m^2)	3331	65%	941	63%	726	87%	72%
Slope
Longitudinal slope (%)	0.3	0.3	0%	1.0	1.0	0%	5.2	5.1	−2%	1%
Mean slope (%)	3.3	2.5	−25%	5.23	6.1	−54%	12.9	13.3	3%	27%
Maximum slope (%)	13.9	12.1	−13%	26.52	33.9	−48%	62.8	65.5	4%	22%
Flow paths
Mean length (m)	40	13	−67%	23	12	−46%	20	14	−28%	47%
Maximum length (m)	173	69	−60%	164	69	−58%	139	66	−53%	57%
Total length (m)	2994	2629	−12%	1062	833	−21%	2846	2290	−20%	18%

First, terrain depressions were analysed. This variable was the most difficult one for the satellite-derived elevation model to detect. For this assessment, we compared both the total area detected as land depressions in the two DEMs and the non-overlapping area between them. While the total area expresses the quantity of terrain surface that was detected as a depression, the non-overlapping area represents the spatial disparity between the areas detected in the two DEMs, making this a quality indicator of the discrepancy between the DEMs. In both cases, these calculations produced mean absolute errors as high as 72%, demonstrating the difficulty of using the WV-2 DEM to detect land depressions. Despite these unpromising results, it may be possible to obtain good characterization of terrain depressions using satellite stereo pairs, for instance, by exploiting spectral information to calculate terrain wetness indices. Such techniques, which are beyond the scope of this paper, are part of the complementary digital image processing techniques available in the remote sensing field.

Next, terrain slope, expressed here as percentage rise, was analysed. This variable is widely used for spatial modelling (for example, in rainfall–runoff models) and as a component in erosion assessment and the Universal Soil Loss Equation (e.g. Nigel and Rughooputh 2013). The longitudinal slope for each plot was calculated over the length of the plot, as shown by the arrow in Fig. 2. There was little variation (1%) between the two elevation models for the longitudinal slope, showing that the satellite model was able to correctly represent the general slope trend of an agricultural plot. However, we found that the WV-2 DEM had difficulties estimating mean and maximum slopes, with MAEs of 27% and 22%, respectively. On rolling terrain (Plot 3), slope calculation performed well, with small errors of 3% and 4% in the mean and maximum values (in contrast to flatter plots). In general, the satellite DEM underestimated slopes, confirming the spatial pattern found in the transversal profiles (Fig. 4), where a loss of local vertical details can be seen in the WV-2 model.

For runoff flow paths, the mean, maximum and total lengths of the water streams extracted were calculated. The high relative errors obtained indicate that the WV-2 DEM also underestimated this variable. The total length calculation had the best overall performance, with 18% relative error, suggesting that WV-2 DEMs may have utility for deriving composite hydrological variables such as drainage density, which is commonly used in spatial hydrological modelling. By contrast, as illustrated in Fig. 6, the WV-2 DEM was unable to accurately determine either the location or the length of flow paths. Such information is necessary to understand the internal hydrological configuration of a plot and design surface drainage structures such as inlet and drainage wells (Beaulieu et al. 2007). To evaluate quantitatively the representation of flow paths, we used two indicators, namely, (a) the drainage outlet and (b) the drainage area of the flow path. For each indicator, relative error (e) and mean absolute relative error (|e|) were calculated for each side of the plot (northwest, southeast, northeast and southwest).

Fig. 6 Elevation, runoff flow paths and drainage points for the satellite- and GPS-based DEM for each plot. Flow paths are shown as black lines, and drainage outlets are shown as: those with a total drainage area of (i) less than 10% and (ii) more than 10% of plot area. A flow accumulation threshold of 100 cells (i.e. 100 m²) was used to extract the flow paths.

For the drainage outlet indicator, the variable used was the number of drainage outlets per side of the plot (Table 6). This indicator shows a high MAE of 76%, with the WV-2 DEM mostly overestimating the number of drainage outlets per side of each plot. In line with the other indicators and variables analysed, the higher the longitudinal slope, the better the performance by the WV-2 DEM, with relative error ranging from 111% in Plot 1 to 43% in Plot 3.

Table 6 Number of drainage outlets derived per side of each plot for the flow paths extracted from the GPS- and WV-2-based DEM (sides with the highest slope are shown in italics).

Plot side	Plot 1			Plot 2			Plot 3
	GPS	WV-2	e	GPS	WV-2	e	GPS	WV-2	e
Northwest	3	5	67%	1	2	100%	12	14	17%
Southeast	7	4	−43%	3	4	33%	16	15	−6%
Northeast	3	7	133%	2	5	150%	2	5	150%
Southwest	4	12	200%	6	5	−17%	3	3	0%
|e| for each plot	111%			75%			43%
Mean |e|	76%

The drainage area indicator, which represents the cumulative area drained towards an outlet, also presented a high overall relative error (90%, Table 7). In keeping with the general trend of better performance of the WV-2 DEM in locations with steeper slopes, sample Plot 3 had the lowest error (19%).

Table 7 Total drainage surface area (m²) per side of each plot for the flow paths extracted from the GPS- and WV-2-based DEM (sides with the highest slope are shown in italics).

Plot side	Plot 1			Plot 2			Plot 3
	GPS	WV-2	e	GPS	WV-2	e	GPS	WV-2	e
Northwest	8085	12969	60%	375	2056	448%	18158	16733	−8%
Southeast	24908	4567	−82%	8197	7352	−10%	10642	11659	10%
Northeast	466	1360	192%	559	676	21%	1862	1397	−25%
Southwest	9680	24948	158%	5658	3964	−30%	795	1051	32%
|e| for each plot	123%			127%			19%
Mean |e|	90%

3.3 General discussion

The performance of the WV-2 DEM in representing plot slopes is better in steep terrain (Plot 3) than flat or gentle terrain (plots 1 and 2). The smaller relative errors obtained in the former show the potential applicability of the WV-2 DEM for environmental and hydrological modelling in such terrain. Multi-temporal VHSR WV-2 DEMs could be highly useful in those parts of the landscape that usually have the steepest slopes, such as riverine and riparian areas. Indeed, in such areas, the slopes usually vary from moderate to high, often within a short distance, and the interactions between agricultural lands and water streams must be monitored regularly to protect and improve the quality of aquatic and environmental resources. Moreover, from a single WorldView-2 stereo pair, users can have access to its eight multispectral bands, and a DEM can be extracted as well, using semi-automated methods. This fact could let organizations optimize their data acquisition costs (with prices of around Can$1.5 ha^-1 for both), the spectral bands and the DEM. Using other data sources these costs might increase to Can$1.3 ha^-1 only for the aerial photos and from Can$0.75 to 4 ha^-1 (depending upon the point density desired) for elevation information coming from LiDAR measurements.

The overall performance of the WV-2 DEM in correctly representing the terrain depressions of the three plots was 72% (mean absolute relative error |e|). Plot 2, which is characterized by gentle slopes (mean ~6%), had the lowest |e| (48%). This means the WV-2 DEM may be better at detecting terrain depressions in intermediate terrain than flatter or steeper terrain. On the whole, however, the performance of the WV-2 DEM in extracting terrain depression information was unsatisfactory. The multispectral imaging capabilities of the satellite may offer an alternative means to extract more accurate depression data, as such regions are usually characterized by higher soil water content, a parameter that can usually be detected in the spectral information of the images (Lyon et al. 2003).

The flow paths of surface runoff in vector format were used to extract three indicators, namely, mean length, maximum length and total length of flow paths. Only total length performed well in all sample plots, showing a mean absolute relative error (|e|) of 18%. This variable can be used to extract other composite variables, such as the drainage density of a plot (useful as an indicator of landscape dissection) and its runoff potential (describing the plot as a whole). The other two flow path indicators, mean length and maximum length, had |e| of 47% and 57%, respectively, confirming that the picture of the internal connectivity and distribution of runoff channels extracted from the satellite DEM is not reliable at the plot scale (i.e. inside the plot).

The WV-2 DEM also performed poorly in drainage outlet analysis, producing a mean absolute relative error (|e|) of 76% across the three plots, with Plot 3 having the lowest value (43%). The overall performance of the WV-2 DEM was lower in drainage area extraction than drainage outlet extraction. However, while the drainage area indicator shows a mean absolute error of 90% across the three plots, the lowest value recorded (Plot 3, 19%) was better than the 43% recorded for drainage outlet extraction. This suggests that drainage area information extracted from WV-2 DEMs may be useful as a parameter in environmental modelling in terrain with a slope generally >5% (and characterizing the plot as a whole).

Vertical error metrics such as RMSE rely on error propagation theory, where it is assumed that vertical errors have a normal distribution, a condition that is not always found in elevation data (Liu et al. 2012). This approach ignores the already-proved influence of land morphology in RMSE calculations (Wilson 2012) and thus reduces the usefulness of RMSE as a proxy for the quality of micro-topography and surface hydrological variables extracted.

The disparities between the drainage outlet and drainage area indicators extracted and the actual parameters of the sites can be seen in Fig. 6. Both indicators show their best performance in Plot 3, revealing that the extraction of hydrological parameters is better in steeper terrain than in flat terrain (except probably for the extraction of terrain depressions, which, as discussed above, was better in Plot 2). Figure 7 shows the global mean absolute relative errors (by plot) and the corresponding longitudinal slope of the plot. The global mean absolute relative error (MARE) of a plot was calculated using the mean of the absolute errors of the three hydrological variables extracted in that plot. We found an inverse proportional relationship between MARE and longitudinal slope, i.e. plots with a low longitudinal slope had an overall higher MARE in all hydrological variables extracted. Plot 3 had the lowest global MARE (37%). Indeed, some variables performed well in Plot 3, with relative errors of less than 20%, indicating that WV-2 DEMs may be useful in such terrain settings.

Fig. 7 Mean absolute relative errors and longitudinal slopes of sample plots.

The satellite-derived elevation model had difficulty extracting accurate surface hydrology features in flat terrain, mainly due to technical issues such as satellite geometric and spectral limitations, as well as automated DEM production algorithms (Callow et al. 2007, Nardi et al. 2008). Further research comparing hydrological variables and terrain complexity indices, as suggested by Huaxing and Tang (2008), is needed to define more meaningful accuracy metrics to assess the quality of micro-topographic hydrological variables extracted.

To better understand the influence of the particular algorithms, a number of studies have compared the performance of different software packages when calculating parameters for hydrological modelling (Zhang et al. 1999, Endreny and Wood 2001, Zhou and Liu 2004, Wilson 2012). When comparing the satellite DEM against the high-precision GPS-RTK dataset, the results show that slopes and flow path lengths were underestimated. In contrast, depressions, drainage areas and the number of outlets were overestimated. The implications of this behaviour in the hydrological modelling outcomes could be the increase in the predicted runoff volume, smaller size of predicted watersheds and sub-watersheds. Moreover, some erratic calculations of solar irradiation and evapotranspiration might occur, as well as a decrease of the mean depth of water table. When modelling erosion and soil moisture, overestimations in erosion rates, biased estimations of erosion-prone hotspots, erroneous identification of hydraulic works locations and overestimation of humid areas could occur.

The results produced by satellite-extracted DEMs, in hydrological studies vary widely, depending on the methodology and the weights of the topographic variables in the hydrological modelling, but they follow the tendencies found in several previously cited studies. For example, Kenward et al. (2000) found that even a minimal decrease in slope values can predict mean annual runoff volumes in a river basin 0.3% to 7% higher than expected. It is worth mentioning that vertical accuracy plays a key role in the hydrological variables produced, in contrast to the horizontal resolution, which is commonly associated with better, more precise results.

Nevertheless, WorldView-2 DEMs offer a practical and fast means to obtain accurate DEMs at high spatial resolution. Elevation accuracy metrics such as the root mean square error and the linear error at the 90th percentile are currently the standards upon which commercial remotely sensed elevation models are based, but further research is required to incorporate new metrics for assessing the quality of micro-topographic representation in such elevation models.

At agricultural plot scales (inside the plot), even with a vertical accuracy of 1 m in a WV-2 DEM product, not all hydrological variables are correctly extracted, proving that there is gap between vertical accuracy metrics and the quality of hydrological variables. Overall, for the three plots, the variables that were most accurately extracted were total length of runoff flow paths (18%), maximum slope (22%) and mean slope (27%). Using sample Plot 3 as a proxy for terrain with rolling and higher slopes (>5%), the variables that presented potential for use were mean slope (3%) and maximum slope (4%). The variables drainage area (19%) and total length of runoff flow paths (20%) also had comparatively good performances, making them potentially useful for representing the plot as a whole, but not for inside-plot representation. Steeper terrain is generally found in hilly and mountainous regions, but steep slopes can also be found in riverine areas of low-relief watersheds, such as the La Chevrotière River basin (the study area), where the landscape is usually dissected to accommodate streams and rivers. Consequently, a satellite DEM can serve as an important decision-making tool to help watershed organizations, government offices and landscape planners characterize, monitor and preserve the natural resources of riverine areas.

Spectral information is one of the most important features that separate a VHSR satellite like WorldView-2 from other options for the extraction of DEMs. Some of the weaknesses of the stereoscopic extraction of elevation can be overcome by developing complementary methods to extract information using both spectral and stereoscopic information. Our examination of the extraction of hydrological variables at plot scales confirmed the limitations of satellite DEMs in flat terrain, but it also showed that satellite DEMs have potential usefulness in terrain with steeper slopes, where better results were obtained. In order to exhaustively exploit the potential of VHSR stereoscopic satellite imagery for surface hydrological variable extraction and environmental modelling, further research is needed using alternative flow routing algorithms, other land uses, spectral information and DEM modification techniques.

4 CONCLUSION

Stereoscopic satellite imagery is recognized nowadays as a practical data source for extracting elevation models. In addition, this technology offers a unique source for both spectral and elevation information at very high spatial resolution (pixel size ≤1 m). Moreover, its high temporal resolution (revisit time) is very convenient for evaluating and monitoring small to medium study areas, and it is very cost-effective, when compared to aerial photographs and LiDAR. The study results show an overall vertical accuracy of 0.45 m in the WV-2 DEM when assessed against ground-truth GPS-RTK high-density and high-precision (0.1 m) point datasets over three agricultural plots with nearly level (0–0.5%), gentle (0.5–2%) and rolling (5–10%) topographic slopes. These results confirm the potential of WV-2 stereo image pairs to extract elevation information at high (1-m) spatial resolution. For the extraction of hydrological variables, the WV-2 DEM performed better on the plot with rolling slopes (5–10%), extracting quality indicators such as total length of flow paths and drainage area with relative errors lower than 20%. This suggests that satellite-derived DEMs have potential for use in environmental modelling applications in areas with rolling or steep slopes. However, the results also indicate some limitations of the WV-2 DEM in the extraction of terrain slope, drainage points and terrain depressions in areas of flatter slopes (<5%).

Acknowledgements

The authors also wish to thank Magali Wirtensohn and Jimmy Poulin for their support in the field campaign and François Riendeau for his efforts in coordinating the generation and delivery of the satellite digital elevation model.

Additional information

Funding

This research was made possible by the support of the Ministry of Agriculture of Quebec (Ministère de l’Agriculture, des Pêcheries et de l’Alimentation, MAPAQ), Project PSIA #810313.