interpolation

[Helena Mitasova]

Dear Grass-users and programmers,

I have found, that recent discussion in this list on 
interpolation has caused some confusion and I feel that
some more clarification is needed.

----------------Simon writes:--------------------------------------------------
>However, it
>is probably worth noting that in most cases, splines are just a fancy form of
>curve-fitting (which is what they were for the wizened old draftspeople) and
>their use often cannot be rigorously justified in terms of the processes being
>modelled.

We are convinced that splines are not just `a fancy form of curve fitting'. 
We would suggest to read a paper by Franke (1982) : Scattered data 
interpolation : tests of some methods (Math. Comput.  v.38, p.181) where 
many methods have been tested and splines were among the most accurate ones 
and also already mentioned paper by Mike Hutchinson and P.Gessler: Splines - 
more than just a smooth interpolator, Geoderma, 1992. 

Splines by itself do not model any process, splines create a curve, surface
or volume function which reproduces the original phenomenon as closely as 
the functional form enables. HOWEVER they can be and they are used as a TOOL
for modeling processes and our and others experience is that they are
very successful in this matter. For example, 
thanks to the high accuracy of spline derivatives we were able to develop
advanced method for modeling erosion and deposition in complex terrain, 
where we use splines not only to interpolate the elevation data but also
to estimate the directional derivative of the surface representing
sediment transport capacity. (Splines seem to be unique among the
interpolation methods in giving high accuracy especially for the estimation
of the second derivatives.) Justification of splines (and also of kriging
which is mentioned later) is a typical A POSTERIORI process : one estimates
the surface from the data and then confronts it with reality using all
information available.

>On its own terms, I think that Kriging, particularly
>the more devloped forms (indicator, log-normal, etc) might be argued to
>give a more honest result for cases where the processes are complex or
>unknown, and it can also give some statistics about how good the estimate is!

Let me quote here Hutchinson (in J.D. Jasper (ed.) 1991, Data Assimilation
Systems", BMRC Research Report N.27, Bureau of Meteorology, Melbourne,p.104):

"The only serious competitor to thin plate spline technique is the method of
 kriging which includes the optimal interpolation method described by 
 Gandin(1963). Like thin plate splines, kriging does extend easily to
 higher dimensions and both methods attempt to achieve minimum error
 (optimum) interpolation. There is in fact a formal inclusion of thin plate 
 splines within kriging (Wahba,1990). The main limitation on kriging is
 that it depends critically on first estimating a spatial covariance function
 or variogram. The method is hampered by an ad hoc assumptions about the form
 that this variogram should take and there are computational difficulties in
 assessing the merit of different functional forms (Dubrule 1983). Nevertheless
 when the variogram is well chosen the results can be similar to those achieved
 by thin plate splines (Seaman and Hutchinson, 1985, Laslett et al. 1987)."
 
Our experience with kriging and splines is similar.
The form of spline representation which is being used in s.surf.tps
and s.surf.3d has a direct correspondence to kriging. The spline S is 
given by

 S(X) = T(X) + sum(i) Ei R(X,Xi)              X=(x,y) or X=(x,y,z)

where T(X) is a low order polynomial and R is a RADIAL BASIS FUNCTION which
is derived by minimizing smoothness functional and which includes a tension
parameter (sum is over the data points, Ei are expansion coefficients). 
Of course, R can be interpreted  as a generalized covariance. 
Thus everybody who uses s.surf.tps works, in fact, with a special type of
kriging with covariance derived by minimizing the bending of the 
surface under a given tension.  The tension  can be optimized by
statistical means like crossvalidation and also other tools of 
geostatistical machinery can be used here. Maros from Bratislava might
be willing to share some experience in relation of variogram and tension.


>Landscapes, for example, are the product of many processes, mostly non-linear,
>etc, so there is no good way of modelling elevations short of process 
>simulations!

 Sure. However, in almost any modeling or process simulation which uses 
 point data the interpolation is an unavoidable and important component.
 There is a lot of experience proving that in many cases splines are a very good choice.
 For example, for computing the hydrologically sound digital elevation
 models Mike Hutchinson uses SPLINES in combination with drainage enforcement
 (what is in fact a simplified flow simulation model).
 Our implemenation of an erosion/deposition model uses SPLINES in combination with 
 unit stream power approach proposed by I.D.Moore and one can find many other 
 similar examples.

>Some people like
>this fractal interpolation since it can be shown that many natural
>surfaces have scale dependent complexity, and the fractal interpolation
>allows this complexity to be preserved as the observation scale is changed.

The problem with  fractal digital terrain models which I have seen so far
was that due to the methods which were used they had approximately
the same number of pits and peaks. Water does not flow in such a
landscape but creates lots of ponds and lakes, so they were not very 
useful for modeling of processes. The problem is that terrain and also
many other natural phenomena are not a 'pure' self-similar fractals
in the sense of Mandelbrot definition. Many phenomena have fractal character
only within rather restricted range of scales or fractal dimension varies
as the scale changes or that the fractal character varies dramatically with
the position changes, but this goes beyond the topic discussed here.

Before finishing this, I would like to let you know that because of
anticipated future use of splines in various types of models we are in
a process of designing an GRASS interpolation library. This is ment to give more
flexibility both to us and to users and will allow to plug in various 
radial basis functions,  add the information on derivatives
(got from data or modeling of processes like flow), use s.semivariogram
and make it work like kriging, combine interpolation with process simulations
etc. This will also allow  more people to participate in the development of
 models using interpolation and approximation tools of GRASS.

I would like to conclude this with another quotation, this time from
a book: Wahba, G.: Spline models for observational data, 1990, CBMS-NSF,
Regional conf. series in applied mathematics, page x (Foreword):

"The body of spline methods available and under development provide a rich 
family of estimation and model building techniques that have found use in
many scientific disciplines. Today, it is hard to open an issue of the J. 
of American Statistical Association, the Annals of Statistics, ....,
without finding the word "spline" somewhere. It is certainly a pleasure 
to be associated with such a blossoming and important area of research."

Whoever needs some more info, explanation, examples or has some comments
please contact me directly at helena at zorro.cecer.army.mil
as I have a feeling that we have already taken too much space in this list.


Helena Mitasova                          
Spatial Analysis and Systems Team       
U.S.Army CERL                          
Champaign, IL 61826