Inverse distance weighting

ALGLIB, a free and commercial open source numerical library, provides the best implementation of the inverse distance weighting algorithm (IDW). ALGLIB is available in multiple programming languages, including C++, C#, Java, and Python.

Our implementation of IDW algorithm includes numerous improvements over both textbook and modified Shepard's methods, making it competitive with ALGLIB's thin plate splines and other scattered interpolation methods.

Using IDW models

Algorithms

Inverse distance weighting has several attractive properties:

• quick and memory-efficient model construction
• either complete absence of tunable parameters or just one parameter (search radius)
• the ability to work in N-dimensional space for N's beyond 3
• the ability to interpolate scattered data on any grid

ALGLIB includes several implementations of the inverse distance weighting algorithm:

• Shepard's method - This is the original version of IDW which uses all available points for model evaluation. It is the most popular version despite having scalability and stability problems (discussed below).
• modified Shepard's method - This is an improved version that uses only points within a search radius, giving O(logN) model evaluation complexity. It's also a popular version, and though it solves scalability problems of the original IDW, it still has issues with interpolation quality.
• modified stabilized method (MSTAB) - This is the latest improvement of the IDW algorithm; it overcomes limitations of two previous methods.

Both the original and modified Shepard's methods suffer from the same problem: on datasets with non-uniform point densities, they tend to give either too little or too much weight to distant nodes. The former results in models with big, flat spots; the latter results in strange interpolation artifacts (see "benchmarks" section below). This problem is surprisingly difficult to solve by playing with the weight formula because any solution tends to flip between two extremes.

MSTAB is our improvement of the modified Shepard's method which borrows ideas from ALGLIB's hierarchical RBF's: the algorithm fits a sequence of smoothed IDW models with a progressively decreasing search radius and smoothing amount, and it then combines these models into a single IDW interpolant. Initial IDW layers capture global trends, and subsequent layers concentrate on local features of the interpolation landscape. The overall setup ensures balanced influence of distant and nearby points.

API overview

IDW interpolation functionality is provided by the idw subpackage. The first step is model construction, which involves the following steps:

• the model builder is initialized using idwbuildercreate
• the algorithm is chosen (preferably with idwbuildersetalgomstab)
• the dataset is added using idwbuildersetpoints
• the model is built by idwfit

After the model is built, the following operations can be performed:

• one-off RBF evaluation with idwcalc and its specialized versions (convenience wrappers for 1..3-dimensional problems or thread-safe idwtscalcbuf)
• model serialization with idwserialize/idwunserialize

The ALGLIB Reference Manual entry for idw subpackage also includes several examples of using IDW models.

Programming languages supported

ALGLIB supports many programming languages, including C++, C#, Java, Python, and others:

• C++ version. Highly optimized, portable, self-contained library (x86/x64/ARM; works on Windows, Linux, and POSIX systems)
• C# version. Unified C# API for two backends: a fully managed C# implementation and a high-performance native computational core written in C
• Java, Python, and other languages. A wrapper for a high-performance native computational core written in C. Seamless integration with the target programming language

A distinctive feature of ALGLIB is that it provides the same API in all programming languages. This is achieved through our exclusive automatic code translation and wrapper generation technology.

Benchmarks

Model quality

We mentioned at the beginning of this article that both original and modified Shepard's methods have issues with interpolation quality and that our MSTAB algorithm solves the problem by using multilayered models.

The original IDW breaks down as soon as the points count increases beyond several hundred. Generally, the modified Shepard's method performs well when used on datasets of any size with uniform points distribution. Uniform geometry, even if highly irregular, is the easiest one—almost any point is surrounded by neighbors in all directions; one can achieve good results by using a search radius proportional to the mean inter-point distance. However, even the modified Shepard's method starts to perform poorly when confronted with highly nonuniform geometries:

Nonuniform geometries often arise in geospatial applications. The plot below visualizes results of three IDW algorithms (original, modified, and MSTAB) applied to a Gaussian-like bump function f=1/(1+5·x²+5·y²) densely sampled across seven parallel lines:

One can see that the original algorithm (on the left) produces strange ripple-like artifacts. As soon as we move away from interpolation nodes, the model rapidly shifts toward mean value over the entire dataset. The reason is that the total weight assigned to distant points outweighs nearby ones.

The modified Shepard's method (middle) with a small search radius goes to the other extreme—the model is mostly dominated by nearby points; it is almost flat when approaching its nearest neighbor. When the search radius is increased, the model shows somewhat better behavior until it rapidly breaks down in a way similar to that of the original IDW algorithm.

Finally, one can see that our MSTAB algorithm produces the best model, which is free from aforementioned artifacts.