IDW, Kriging, Splines or TIN? Choosing an Interpolator

Problem statement

Every GIS offers four or five interpolation methods and no guidance about which one to use. The literature suggests the choice is consequential; the measurements suggest it is less consequential than almost everything else you could do instead.

On 500 sample points taken from a real 30 m elevation model, tested against 12,000 held-out cells:

kriging (fitted spherical model)   40.2 m RMSE
TIN, linear                        46.7 m
IDW, power 2                       54.0 m
nearest neighbour                  65.1 m

Kriging wins by 25% over IDW. Meanwhile, doubling the sample from 500 to 1,000 points improved IDW from 58.4 m to 43.3 m β€” a 26% gain, from the same method.

The method choice is worth roughly one doubling of the sample. Both are worth having; only one of them is under your control after the survey is done.

Quick answer

method use it when cost
Nearest neighbour the variable is categorical, or you need exact blockiness worst accuracy for continuous data
IDW you need something defensible in ten lines, samples are dense and even bullseyes; no uncertainty estimate
TIN / linear you want a surface that refuses to extrapolate no values outside the convex hull
Kriging the variogram shows real structure and you need uncertainty needs a variogram; slow; variance is not a confidence interval
Splines / RBF the surface is genuinely smooth (a physical field) can overshoot wildly near sharp changes
A decision path from the type of variable through sample density and the need for uncertainty to one of nearest neighbour, IDW, TIN or kriging.
Three questions settle it: is the variable continuous, do you need uncertainty, and can you afford a variogram?

Step-by-step solution

1. Rule out interpolation entirely if there is no spatial structure

Fit a variogram first. If semivariance is flat from the shortest lag, nearby samples are no more alike than distant ones and the best estimate everywhere is the mean. No method fixes that.

2. Categorical variable? Use nearest neighbour and stop

Soil class, land cover, geological unit β€” these have no meaningful average. Nearest neighbour (a Voronoi tessellation) is the only defensible choice, and the blockiness is honest: it says "this is the nearest observation" rather than pretending to a smooth transition that was never measured.

3. Need a number today? Use IDW, with k set properly

IDW is a weighted average with weights 1/d^p. Its parameters have a broad optimum in one and a sharp cliff in the other. Measured at 500 samples:

power  0.5     1      2      3      4      6
RMSE  72.2   64.4   54.0   50.7   50.9   53.2

k      1      3      6     12     24     48    500
RMSE  64.5  52.5   51.5   54.0   59.6   66.9   90.1

Power barely matters between 2 and 4. Neighbour count matters a great deal, and the failure is on the high side: using all 500 points cost 75% accuracy against using the nearest 6, because distant samples drag every prediction towards the global mean.

Default to power 2 and k between 6 and 12, and tune k by cross-validation before touching the power.

4. Want a surface that will not extrapolate? Use a TIN

Linear interpolation over a Delaunay triangulation produces values only inside the convex hull of the samples. Measured coverage:

samples   TIN coverage of the study area
    100        85.2%
    250        93.5%
    500        97.3%
  1,000        99.0%

Everything outside is NaN. That is often exactly what you want β€” see My interpolated surface extends far beyond the data. It also means TIN's RMSE is not comparable with methods that fill the whole area, because it declined to answer where the answer was hardest.

5. Need uncertainty, or the best accuracy? Use kriging

Kriging derives its weights from the fitted variogram rather than from a typed-in exponent, which is why it won here. It also returns a variance per prediction.

Read that variance carefully. Measured on the same 12,000 held-out cells:

mean kriging variance     4,536 mΒ²  ->  predicted sd 67.3 m
actual RMSE                                          40.2 m

The kriging variance overstated the real error by 68%. It depends only on the geometry of the samples and the fitted variogram, never on their values β€” so it is a good relative map of where the surface is weak, and a poor absolute confidence interval.

IDW error against power and against neighbour count, showing a broad optimum in power and a steep penalty for too many neighbours.
Power is forgiving between 2 and 4. Neighbour count is not: all 500 points cost 75% accuracy against the nearest six.

Code examples

Example 1 β€” pick the method by measuring, not by reading

import numpy as np
from scipy.interpolate import griddata, RBFInterpolator
from scipy.spatial import cKDTree


def bake_off(points, values, test_points, test_values, k=12):
    """Run every candidate on your own data and rank them."""
    def rmse(pred):
        ok = np.isfinite(pred)
        return float(np.sqrt(np.mean((pred[ok] - test_values[ok]) ** 2))), \
            float(ok.mean())

    tree = cKDTree(points)

    def idw(power):
        d, i = tree.query(test_points, k=min(k, len(points)))
        d = d.reshape(len(test_points), -1); i = i.reshape(len(test_points), -1)
        w = 1.0 / np.maximum(d, 1e-9) ** power
        return (w * values[i]).sum(1) / w.sum(1)

    candidates = {
        "nearest":     griddata(points, values, test_points, method="nearest"),
        "tin_linear":  griddata(points, values, test_points, method="linear"),
        "tin_cubic":   griddata(points, values, test_points, method="cubic"),
        "idw_p2":      idw(2),
        "idw_p3":      idw(3),
        "rbf_thin_plate": RBFInterpolator(
            points, values, neighbors=min(50, len(points)),
            kernel="thin_plate_spline")(test_points),
    }

    rows = [(name, *rmse(pred)) for name, pred in candidates.items()]
    for name, err, cov in sorted(rows, key=lambda r: r[1]):
        print(f"  {name:16} RMSE {err:7.2f}   coverage {cov:6.1%}")
    return candidates

Always print the coverage column beside the error. A method that answers 85% of the map with a low error has not beaten one that answered all of it.

Example 2 β€” tuning k properly, which is the parameter that matters

import numpy as np
from scipy.spatial import cKDTree


def tune_neighbours(points, values, candidates=(3, 6, 9, 12, 18, 24, 36)):
    """Leave-one-out error against neighbour count."""
    tree = cKDTree(points)
    best = None

    for k in candidates:
        # k+1 because the nearest neighbour of a sample is itself
        d, i = tree.query(points, k=min(k + 1, len(points)))
        d, i = d[:, 1:], i[:, 1:]
        w = 1.0 / np.maximum(d, 1e-9) ** 2
        pred = (w * values[i]).sum(1) / w.sum(1)
        err = float(np.sqrt(np.mean((pred - values) ** 2)))
        print(f"  k={k:3d}  LOO RMSE {err:7.2f}")
        if best is None or err < best[1]:
            best = (k, err)

    print(f"  best k = {best[0]}")
    return best[0]

Excluding the point itself is the whole trick. Forget the [:, 1:] slice and every prediction is the sample's own value, LOO RMSE is zero, and the tuning silently picks the largest k.

Example 3 β€” a hybrid that is often the right answer

import numpy as np
from scipy.spatial import cKDTree


def local_mean_fallback(points, values, targets, radius, k=12, power=2):
    """IDW inside the variogram range, the local mean beyond it, NaN far out."""
    tree = cKDTree(points)
    d, i = tree.query(targets, k=min(k, len(points)))
    d = d.reshape(len(targets), -1); i = i.reshape(len(targets), -1)

    w = 1.0 / np.maximum(d, 1e-9) ** power
    predicted = (w * values[i]).sum(1) / w.sum(1)

    nearest = d[:, 0]
    inside = nearest <= radius
    far = nearest > radius * 3

    # beyond the range, distance weighting is meaningless β€” average instead
    predicted[~inside] = values[i[~inside]].mean(axis=1)
    predicted[far] = np.nan

    print(f"  {inside.mean():.1%} interpolated, "
          f"{(~inside & ~far).mean():.1%} local mean, {far.mean():.1%} masked")
    return predicted

Beyond the variogram range, distance carries no information, so weighting by it is theatre. Switching to a plain local mean and then to NaN is more honest than letting 1/dΒ² produce a confident number four kilometres from anything.

Explanation

Why they all agree more than the literature suggests

All five methods produce sum(w_i * z_i) for weights that depend only on geometry. Where samples are dense, every reasonable weighting scheme puts almost all the weight on the same few nearby points, so they converge.

Where samples are sparse, they diverge β€” and where samples are sparse, all of them are wrong. Measured across distance bands from the nearest sample, IDW's error went from 18.9 m within 100 m to 334.3 m beyond a kilometre. A better interpolator moves the 334 m to perhaps 300 m. Another survey day moves it to 50 m.

Why kriging actually wins

Kriging solves for weights that minimise expected squared error under the fitted variogram. Three consequences:

  • Declustering. Two samples close together carry partly redundant information, and kriging downweights them jointly. IDW counts both at full weight.
  • A principled search radius. The variogram range says where correlation ends; kriging's weights fall to near zero there naturally.
  • Anisotropy. If the structure is directional, the variogram model can say so and the weights follow.

The 25% gain measured here is a fair expectation for a well-behaved continuous variable with a clean variogram. It is not free: you must fit and inspect the variogram, and a bad fit produces a worse surface than IDW.

Why splines are the risky choice

Radial basis functions and thin-plate splines fit a smooth function through the samples exactly. Where data are dense and the underlying field is genuinely smooth β€” a pressure field, a geoid β€” they are excellent.

Where the field has a step, they overshoot: to pass through both sides of a cliff smoothly, the fitted surface must swing beyond both. The result is elevations above the summit and below the valley floor, with no warning. Always check the output range against the input range; if the surface exceeds the sample extremes, a spline has overshot.

Why nearest neighbour deserves more respect than it gets

It came last here at 65.1 m, only 20% behind IDW at 54.0 m, for a fraction of the thinking. And it has two properties nothing else has: it never invents a value that was not measured, and it works for categorical data.

For a first look at an unfamiliar dataset, a Voronoi map is more honest than a smooth surface. Its blockiness shows you the sampling design, which a smoothed surface hides β€” and the sampling design is the thing that actually determines your accuracy.

Kriging implying a 67.3 m standard deviation where the actual RMSE was 40.2 m.
Useful as a relative map of where the sample geometry is weak. Not a confidence interval.

Edge cases or notes

  • Categorical data has no average. Nearest neighbour only.
  • Tune k before tuning the power. Neighbour count cost 75% here; power cost 6%.
  • A large k drags predictions to the global mean and flattens the map.
  • TIN's error is not comparable with methods that fill the whole area.
  • Kriging variance is geometry, not accuracy. It overstated real error by 68% here.
  • Check splines for overshoot against the sample minimum and maximum.
  • Duplicate points break kriging (singular matrix). Aggregate first.
  • Use a projected CRS β€” every one of these methods weights by distance.
  • Anisotropy is invisible to IDW. If structure is directional, only kriging can use that.

FAQ

Which interpolation method is most accurate?

On the data measured here, kriging at 40.2 m RMSE, then linear TIN at 46.7 m, IDW at 54.0 m and nearest neighbour at 65.1 m. The ranking depends on the variable; run the bake-off on your own data.

Is kriging worth the extra effort?

It bought 25% over IDW here β€” about the same as doubling the sample size. Worth it when you also need the uncertainty map, or when the field is anisotropic.

What IDW power should I use?

Between 2 and 4; the difference was 6% here. Spend your tuning effort on the neighbour count instead, where the difference was 75%.

How many neighbours should IDW use?

Six to twelve. Using all 500 samples raised RMSE from 51.5 m to 90.1 m, because distant points pull every prediction towards the global mean.

Can I use kriging variance as a confidence interval?

Not directly. It depends only on sample geometry and the variogram, not on the data values, and it overstated the true error by 68% here. Use it as a relative map of where the surface is weak.

When should I use a spline?

When the underlying field is genuinely smooth and the samples are dense. Check for overshoot: if the surface exceeds the sample minimum or maximum, it has invented extremes.

Why does my TIN surface have NaN around the edges?

Because linear TIN interpolation only produces values inside the convex hull of the samples β€” 85.2% of the study area at 100 points, 97.3% at 500. That is the method refusing to extrapolate.