A comparative analysis of different DEM interpolation methods.

Arun, Pattathal Vijayakumar

Introduction

Remote sensing techniques are being effectively used as a tool for decision making in various fields because of its spatial analysis and display capabilities. The utility of decision making processes are significantly improved using 3D geographical models as they facilitate effective visualization. Digital Elevation Models (DEMs) are the generally adopted data structures for storing topographic information and are usually interpolated to establish the values for entire terrain points. DEM is an array representation of squared cells (pixels) with an elevation value associated to each pixel (Peralvo 2009). DEMs can be obtained from contour lines, topographic maps, field surveys, photogrammetry techniques, radar interferometry, and laser altimetry (Peralvo 2009). Different interpolation methods applied over the same data sources may result in different results and hence it is required to evaluate the comparative suitability of these techniques.

Interpolation techniques are based on the principles of spatial autocorrelation, which assumes that closer points are more similar compared to farther ones. Literature reveals a great deal of interpolation methods which are generally classified as local and global approaches. Local methods predict value of an unknown point based on the values of neighborhood pixels. Prominent local methods found in literature include Inverse Distance Weighting (IDW), local polynomial, Natural Neighbor (NN), and Radial Basis Functions (RBFs). On the other hand, global interpolation methods such as polynomial interpolation functions use all the available sample points to generate predictions for a particular point. These methods facilitate to evaluate and remove global variations caused by physical trends in the data.

Kriging is a geo statistical interpolation method that utilizes variogram which depend on the spatial distribution of data rather than on actual values. Kriging weights are derived using a data-driven weighting function to reduce the bias towards input values, and it provides best interpolation when good variogram models are available (Pincock, Allen & Hol 2008). IDW approach is a local deterministic interpolation technique that calculates the value as a distance-weighted average of sampled points in a defined neighborhood (Burroughs, McDonnell 1998). It considers that points closer to the query location will have more influence, and weights the sample points with inverse of their distance from the required point (Johnston et al. 2001).

Natural neighbor interpolation finds the closest subset of input samples to a query point and applies weights to them based on proportionate areas (Sibson 1981). It is a local deterministic method and interpolated heights are guaranteed to be within the range of the samples used. It does not produce peaks, pits, ridges or valleys that are not already present in the input samples and adapts locally to the structure of the input data. It does not require input from the user and works equally well for regularly as well as irregularly distributed data (Watson 1992). Spline interpolation approach uses mathematical function to minimize the surface curvature and produces a smooth surface that exactly fits the input points. Topo to Raster method uses an interpolation technique specifically designed to create a surface that more closely represents a natural drainage surface and preserves both ridgelines as well as stream networks (Hutchinson et al. 1989).

Zimmerman et al. (1999) showed that Kriging yielded better estimations of altitude than inverse distance weighting (IDW) irrespective of the landform type and sampling pattern. This result is attributed to the ability of Kriging to adjust itself to the spatial structure of the data. However, in other studies (Weber, Englund 1992; Gallichand, Marcotte 1993; Brus et al. 1996; Declercq 1996; Aguilar et al. 2005), neighborhood approaches such as IDW or RBFs were found to be as accurate as Kriging or even better. Topo to Raster interpolation method is specifically designed for the creation of hydrologically correct terrain surfaces.

In this paper, we evaluate the comparative suitability of different interpolation techniques based on their accuracy and sensitivity to terrain variations. Performance of different interpolation methods namely IDW, ordinary Kriging (KRG), Topo to Raster, NN and Spline have been evaluated with reference to the study area. Generally available DEMs for Indian terrain namely SRTM and SOI DEMs are also evaluated based on the contours generated at different intervals.

1. Data resources

Investigations have been conducted over MANIT campus and surrounding areas of Bhopal city in India; variation of the terrain, spread over more than 1000 acres made it optimal for the analyses. Satellite images of Bhopal along with SOI & SRTM DEMs have been used for comparative analysis of various methodologies. Details of the satellite data used for these investigations are summarized in Table 1. The ground truthing information has been collected using Differential Geographic Positioning System (DGPS) survey conducted over Bhopal during October 2012.

2. Methodology

2.1. Comparative analysis of interpolation methods

Commonly used interpolation approaches have been evaluated with reference to the study area and adopted methodology is summarized in (Fig. 1). DGPS survey has been conducted over the study area to collect three-dimensional coordinates of around 1000 sample and test points in WGS-84 datum. Collected raw data has been pre-processed using GNSS software to remove various errors and to calibrate the readings at centimeter level accuracy. The processed data (GCPs) has been imported in the ArcGIS environment and plotted to a shape file. About 680 GCPs were used as sample points to generate the DEM and rest were used as test points to estimate accuracy of interpolation. Raster surface has been generated from reference DEM using different interpolation methods namely IDW, Kriging, NN, Topo to raster and Spline. Accuracies of generated surfaces have been evaluated using 320 reference GCPs as test points. Visual analyses as well as statistical parameters have been adopted for comparative evaluation of the interpolated surfaces. In the visual analysis, DEM generated heights were verified in the ground by field visit using GPS. Mathematical analysis has been done by calculating the deviations of interpolated height values from corresponding observed values in terms of root mean square error (RMSE).

[FIGURE 1 OMITTED]

2.2. Comparative analysis of SRTM and SOI DEM

Comparative suitability of SRTM and SOI DEMs has been analyzed with reference to the generation of contours. Contours of the study areas have been digitized from SOI Topo sheet no. 55E7 & 55E8 and contour heights were recorded in the attribute table. SOI DEM has been generated from corresponding contours using Kriging interpolation technique in the ARCGIS environment. Contours with interval 10 m, 5 m, 2 m and 1 m were generated from SRTM as well as SOI DEM using Arc GIS 3D analyst extension. Comparative analysis has been done with reference to the nature and number of contours generated from DEMs. Further, visual analysis has been conducted based on the 3D view generated from the two DEMs. Satellite images were draped over the DEMs using Virtual GIS viewer in ERDAS and were analyzed at different exaggeration levels.

3. Results and discussions

3.1. Comparitive analysis of interpolation methods

We have investigated the comparative performance of different interpolation techniques with reference to various terrain contexts. Visual comparisons as well as mathematical analyses have been conducted. Visual comparison of slope map generated using different interpolation techniques is presented in (Fig. 2).

DGPS survey data revealed that Kriging approach performed accurately in average cases when compared to others. Interpolated heights at different test points (points having coordinates from DGPS survey) have been also compared for the five different methods and results are summarized in Table 2.

Table 2 reveals that different approaches produce varied results over the same points. Interpolated height values for different methods at each test points have been plotted. Deviations of interpolated height values from the actual values (DGPS observed) at each test points gives a better understanding about the performance of each method and reveals a better performance of Kriging approach.

In order to investigate the sensitivity of interpolation methods to the nature of terrain, the test GCPs were divided into two zones namely mild slope and steep slope areas. Average RMSE values of the test points have been also calculated with reference to terrain variations and are summarized in Table 3. IDW and Kriging have been found to adjust themselves to the terrain variations when compared to other methods. Topo to Raster has been found to yield a better performance for ridges as well as stream areas.

The investigations have shown that interpolation results vary with variation in spatial structure and terrain nature of input data. As far as our data is concerned, we have more samples at slope areas than at plane areas. Kriging and NN were found to perform well in these contexts and can be adopted for geomorphologically smooth and small areas. In stream and ridge line areas, Topo to Raster method has shown lowest RMSE value. The NN method has shown nearly optimal values over smooth surfaces, i.e. second lowest.

[FIGURE 2 OMITTED]

This trend in RMSE values of Kriging have continued even for steep slope areas as well as for areas covering both steep and mild slopes. IDW and NN method has been found to be good for interpolation of geomorphologically smooth areas. Kriging methods take into consideration autocorrelation structures of elevations in order to define optimal weights. The method requires a skilled user with geostatistical knowledge. Spline-based methods fit a minimum-curvature surface through the input points, and ensure preservation of trend in the sample data along with rapid changes in gradient or slope.

3.2. Comparative analysis of SRTM and SOI DEM

We have investigated the accuracy of DEMs namely SRTM and SOI with reference to contour extraction. Contours have been generated using 3D analyst extension of Arc GIS software and outcomes of these investigations are tabulated in Table 4.

From the table, it is evident that the contours generated from SOI DEM are sparse while that from SRTM are comparatively denser. Therefore we can conclude that SOI DEM is having very poor data quality compared to SRTM.

The suitability of DEMs has also been evaluated based on the comparative visualization of 3D models generated from these DEMs at different exaggeration levels as given in Figure 3.

[FIGURE 3 OMITTED]

Visual comparison also reveals that SRTM is performing better than the SOI DEMs. Reason behind the poor performance of the SOI DEM may be attributed to its construction from 1:50000 scale topographic maps. Open source SRTM data is giving more reliability and accuracy than the SOI DEM due to the usage of Radar technology.

Conclusions

The generated DEMs are found to be sensitive to height interpolation methods as well as the terrain nature. Investigations revealed that the Krigging method performs better when compared to other contemporary methods in most contexts. DEM generated from the DGPS data was found to be better than the DEM available from SOI or SRTM data. Number of contours extracted from SRTM DEM was found to be better than that from SOI DEM, which may be attributed to the better accuracy of SRTM data source. Krigging has been found to adapt itself to terrain variations while Topo to Raster is found preferable for streams and ridge lines.

doi: 10.3846/20296991.2013.859821

Caption: Fig. 1. Methodology for comparative analysis of interpolation methods

Caption: Fig. 2. Slope maps generated using different interpolation methods

Caption: Fig. 3. 3D surfaces generated from SOI & SRTM DEM at different exaggeration Levels

Received 21 July 2013; accepted 09 December 2013

References

Aguilar, F. J.; Aguera, F.; Aguilar, M. A.; Carvajal, F. 2005. Effects of terrain morphology, sampling density, and interpolation methods on grid DEM accuracy, Photogrammetric Engineering and Remote Sensing 71: 805-816.

Brus, D. J.; Gruijter, J. J.; Marsman, B. A.; Visschers, R.; Bregt, A. K.; Breeuwsma, A. 1996. The performance of spatial interpolation methods and choropleth maps to estimate properties at points: a soil survey case study, Environmetrics 7: 1-16. http://dx.doi.org/10.1002/ (SICI)1099-095X(199601)7:1<1::AID-ENV157>3.0.CO;2-Y

Burrough, P. A.; McDonnell, R. A. 1998. Principles of geographical information systems. New York: Oxford University Press, 333-335.

Declercq, F. A. N. 1996. Interpolation methods for scattered sample data: accuracy, spatial patterns, processing time, Cartography and Geographic Information Systems 23: 128-144. http://dx.doi.org/10.1559/152304096782438882

Gallichand, J.; Marcotte, D. 1993. Mapping clay content for subsurface drainage in the Nile delta, Geoderma 58: 165-179. http://dx.doi.org/10.1016/0016-7061(93)90040-R

Hutchinson, M. F. 1989. Calculation of hydrologically sound digital elevation models, in Third International Symposium on Spatial Data Handling, 23-28 January, 1989, Sydney, Australia 3(1): 120-127.

Peralvo, M. 2009. Influence of DEM interpolation methods in Drainage Analysis, CE 394K GIS in Water Resources, 301-312.

Pincock, Allen & Holt. 2008. About Kriging Pincock, Perspectives 88(1): 213-218.

Sibson, R. 1981. A brief description of Natural Neighbor Interpolation, in Interpolating multivariate data, Chapter 2. New York: John Wiley & Sons, 21-36.

Watson, D. 1992. Contouring: a guide to the analysis and display of spatial data. London: Pergamon Press, 120-123.

Weber, D.; Englund, E. 1992. Evaluation and comparison of spatial interpolators, Mathematical Geology 24: 381-391. http://dx.doi.org/10.1007/BF00891270

Zimmerman, D.; Pavlik, C.; Ruggles, A.; Armstrong, M. 1999. An experimental comparison of ordinary and universal kriging and inverse distance weighting, Mathematical Geology 31: 375-390. http://dx.doi.org/10.1023/A:1007586507433

Pattathal Vijayakumar Arun

College of Science & Technology, Phuentsholing, Bhutan

E-mail: arunpv2601@gmail.com

Pattathal Vijayakumar ARUN. He has completed his Masters from NIT-Bhopal, and is currently pursuing PhD. His main area of interest are artificial Intelligence, spatial mining, and image processing.

Table 1. Data resources description

S.NO    Image used       Resolution (m)            Satellite

1          PAN                 2.5            IRS-P5 (Cartosat-1)
2        LISS-IV               5.8                  IRS P6
3      Google Earth       0.15(Highest)
5        SOI DEM      As per 1:50,000 scale            -
                           topo sheet
6        SRTM DEM             3-ARC              Shuttle Radar

S.NO    Area    Date of procurement

1      Bhopal      November 2012
2      Bhopal     September 2012
3      MANIT             -
5      Bhopal      November 2012
6      Bhopal       August 2012

Table 2. Ellipsoidal heights at test GCPs from interpolated and
DGPS observed values

                      Ellipsoidal heights in meter

Control      IDW      Topo to   Natural    Spline   Kriging
points ID    value    Raster    Neighbor

FID-24      474.97    474.89     476.93    476.33   476.02
FID-81      476.64    476.78     476.72    477.54   476.68
FID-39      479.42    478.63     479.22    479.37   478.90
FID-7       476.78    477.30     475.59    475.66   475.83
FID-11      477.75    477.30     477.78    477.30   477.29
FID-14      478.06    477.48     478.59    480.76   479.63
FID-17      479.27    477.85     479.32    480.36   479.44
FID-71      477.38    476.26     477.18    477.65   477.64
FID-64      477.39    478.41     478.05    478.18   477.93
FID-61      479.30    480.06     479.59    477.79   479.03
FID-56      477.72    478.11     478.40    477.96   478.20
FID-45      477.93    478.93     477.57    477.36   477.43
FID-39      479.450   478.32     479.22    479.123  478.75
FID-91      473.11    471.95     473.94    474.50   474.18
FID-87      473.02    471.95     473.56    473.25   473.39
FID-30      474.32    474.83     473.67    471.40   472.58
FID-28      473.80    475.57     473.41    473.32   473.47
FID-89      473.14    471.95     472.31    472.33   472.22
FID-95      471.08    471.95     471.27    471.15   471.11
FID-34      477.07    478.01     477.90    477.14   477.43

            Ellipsoidal heights in meter

Control     DGPS observed
points ID       value

FID-24          476.72
FID-81          478.54
FID-39          477.60
FID-7           478.40
FID-11          478.24
FID-14          479.58
FID-17          479.37
FID-71          476.68
FID-64          477.28
FID-61          479.69
FID-56          475.72
FID-45          477.97
FID-39          477.60
FID-91          475.83
FID-87          476.42
FID-30          473.19
FID-28          471.82
FID-89          469.89
FID-95          472.14
FID-34          477.82

Table 3. RMSE values with reference to terrain variation

                                      RMSE values

Type of test GCPs      IDW     Topo to     NN     Spline   Kriging
used                           Raster

Mild slope areas      0.9367   0.8764    0.7288   0.9170   0.7067
Steep slope area      1.4579   1.8200    1.3477   1.3785   1.3137
Combined slope area   1.7329   2.0201    1.5322   1.6247   1.4918

Table 4. Comparison of contours

                            Number of contour generated

Type DEM used   1m interval  2m interval  5m interval  10m interval

SOI                3303         1637         663          342
SRTM              12182        11212         4274         2153