8

I am using the griddata function in scipy to interpolate 3 and 4 dimensional data. It works like a champ, except that it returns a bunch of NaNs because some of the points I need are outside the range of the input data. Given that N-d data only works with the "linear" mode interpolation anyway, it should be a snap to have griddata do an extrapolation instead of just returning NaN. Has anyone done this or found a workaround? To clarify: I have unstructured data, so I can't use any of the functions that require a regular grid. Thanks! Alex

denfromufa
  • 5,610
  • 13
  • 81
  • 138
Alex
  • 427
  • 6
  • 12
  • Would it help to fill the points outside the range with some constant value? In that case, you could just specify the fill_value – Dhara Jun 26 '12 at 19:34
  • 3
    Also, are you sure you want to extrapolate? sometimes, getting out NaNs and knowing you are going out of range is a much better choice. I have used Univariate splines from scipy, it silently extrapolates and the results can be quite "off" – Dhara Jun 26 '12 at 19:37
  • My situation is: I measure some values at a few points, and need to then calculate values at a bunch of other points through inter/extrapolation. So a constant value, or NaN really don't help. I know how temperamental the splines can be, so I was thinking linear would be a safe bet. I would like something that works on N-d data though. – Alex Jun 27 '12 at 13:28
  • This is actually a bit more nontrivial than it seems in these triangulation-based methods. For each point outside the tesselation, you need to choose *which* simplex you extrapolate from, and also have a fast algorithm to find this simplex. Also, if you want to avoid discontinuities in the extrapolation, some care is needed. After that is sorted out, however, the linear extrapolation is a simple matter. – pv. Jun 27 '12 at 16:36
  • I am just going to do a 2nd run as nearest and use that to fill in the gaps. thanks! – Alex Jul 02 '12 at 16:23

2 Answers2

5

One possibility to interpolate & extrapolate data with 3, 4 or actually any dimensions is with scipy.interpolate.Rbf

The get_data() function and plot_3d() function are attached to the end for convenience.

Example data

The example data looks like this (fourth dimension, w, is shown with a color). The data is irregularly spaced and not gridded.

enter image description here

x, y, z, w = get_data(N=200)
plot_3d(x, y, z, w)

Interpolation & extrapolation in 3d

First, let's setup new x and y coordinates. To make this more interesting, let's extrapolate to minus x and minus y direction. This forms the new x and y range of interest.

xs = np.linspace(-10, 20) # some extrapolation to negative numbers
ys = np.linspace(-10, 20) # some extrapolation to negative numbers
xnew, ynew = np.meshgrid(xs, ys)
xnew = xnew.flatten()
ynew = ynew.flatten()

Interpolation with scipy.interpolate.Rbf. Now,

from scipy.interpolate import Rbf
rbf3 = Rbf(x, y, z, function="multiquadric", smooth=5)
znew = rbf3(xnew, ynew)

plot_3d(xnew, ynew, znew)
  • There can be as many variables/dimensions as you want. The first arguments (x,y) are treated as coordinates for the nodes. The last argument before the function argument are the "values" which are to be interpolated (now: z).
  • The function argument can be used to control how to values are interpolated. This will affect the results so play with it with your data..
  • The smooth parameter can be used to smooth some noise from the data. If smooth is zero, the result is interpolating; it will go through all your data points. If it is positive value, the data is smoother. This will affect the results so play with it with your data..
  • Below are the results and the extrapolation is of course badly off. This is just to demonstrate that extrapolation is possible. You might want to fine tune function and smooth for desired results. Usually, data should not be extrapolated "much" (like in this example)

enter image description here

Adding fourth dimension

It is also possible to interpolate & extrapolate to the fourth dimension. Here is how:

rbf4 = Rbf(x, y, z, w, function="thin_plate", smooth=5)
wnew = rbf4(xnew, ynew, znew)
plot_3d(xnew, ynew, znew, wnew)
  • I created another Rbf instance for the fourth dimension and I used the znew calculated with the rbf3 (interpolation in 3d).
  • I changed the function to "thin_plate" since it seemed visually to perform better with this dataset.
  • Here is how the results looks like: enter image description here

Appendix: get_data and plot_3d

For testing purposes:

import numpy as np

def get_data():
    np.random.seed(100)
    N = 200
    maxval = 20
    x = np.random.random(N) * maxval
    y = np.random.random(N) * maxval
    z = x ** 2 + np.sqrt(y) * y - y ** 3 + np.random.random(N) + 18 * y ** 2 * 2
    w = x ** 2 - np.log(y + (x * y) ** 2)
    return x, y, z, w


def plot_3d(x, y, z, w=None, show=True):
    import matplotlib.pyplot as plt
    from mpl_toolkits.mplot3d import axes3d

    fig = plt.figure(figsize=(10, 6))
    ax = axes3d.Axes3D(fig)
    ax.scatter3D(x, y, z, c=w if not w is None else "b")
    plt.show()
Niko Föhr
  • 28,336
  • 10
  • 93
  • 96
0

Not quite sure this will work for you and it is not available yet, but in the development version of numpy there is a 'pad' array function...

https://github.com/numpy/numpy/blob/master/numpy/lib/arraypad.py

One of the options is 'linear_ramp' which extrapolates (pads) outward starting at the edge value and linearly increase/decreases to a specified end value.

It is a pure python function so you could just copy it into your path and import (untested by me though)

TimCera
  • 574
  • 3
  • 8
  • looks totally unrelated to me. here are docs, nothing related to interpolation/extrapolation https://docs.scipy.org/doc/numpy-1.10.0/reference/generated/numpy.pad.html – denfromufa Jun 21 '17 at 19:51