MyCarta
MyCarta

Reputation: 818

Difference in x,y parameters for scipy interpolate RectBivariateSpline and interp2d

If I want to interpolate the data below:

from scipy.interpolate import RectBivariateSpline, interp2d
import numpy as np

x1 = np.linspace(0,5,10)
y1 = np.linspace(0,20,20)
xx, yy = np.meshgrid(x1, y1)
z = np.sin(xx**2+yy**2) 

with interp2d this works:

f = interp2d(x1, y1, z, kind='cubic')

however if I use RectBivariateSpline with the same x1, y1 parameters:

f = RectBivariateSpline(x1, y1, z)

I get this error:

   TypeError                                 Traceback (most recent call last)
<ipython-input-9-3da046e1ebe0> in <module>()
----> 1 f = RectBivariateSpline(x, y, z)

C:\...\Local\Continuum\Anaconda\lib\site-packages\scipy\interpolate\fitpack2.pyc in __init__(self, x, y, z, bbox, kx, ky, s)
    958             raise TypeError('y must be strictly ascending')
    959         if not x.size == z.shape[0]:
--> 960             raise TypeError('x dimension of z must have same number of '
    961                             'elements as x')
    962         if not y.size == z.shape[1]:

TypeError: x dimension of z must have same number of elements as x

I'd have to switch the sizes of x, y like this to have it work:

x2 = np.linspace(0,5,20)
y2 = np.linspace(0,20,10)
f = RectBivariateSpline(x2, y2, z)

Is there a reason for this behavior - or something I am not understanding?

Upvotes: 0

Views: 2630

Answers (1)

Jon Custer
Jon Custer

Reputation: 704

Well, the reason is that the parameters to the two functions are, as you have noted, different. Yes, this makes it really hard to just switch out one for the other, as I well know.

Why? In general it was a clear design decision to break backward compatibility with the new object-oriented spline functions, or at least not worry about it. Certainly, for large grid sizes there is significant space savings with not having to pass x and y as 2D objects. Frankly, I have found in my code that once this initial barrier is overcome, I'm much happier using the spline objects. For example, with the UnivariateSpline object, getting the derivative(s) is easy, as is the integral.

It would appear that, going forward, the SciPy folks will focus on the new objects, so you might contemplate just moving to them now. They are the same base functionality, and have additional methods that provide nice benefits.

EDIT - clarify what 'broke' between the two.

From the SciPy manual on interp2d you get the code snippet:

from scipy import interpolate
x = np.arange(-5.01, 5.01, 0.25)
y = np.arange(-5.01, 5.01, 0.25)
xx, yy = np.meshgrid(x, y)
z = np.sin(xx**2+yy**2)
f = interpolate.interp2d(x, y, z, kind=’cubic’)

This can be, unfortunately, potentially misleading since both x and y are the same length, so z will be a square matrix. So, lets play with this a bit:

x = np.linspace(0,5,11) 
y = np.linspace(0,20,21) # note different lengths
z = x[None,:].T + y*y # need broadcasting
xx,yy = np.meshgrid(x,y) # this is from the interp2d example to compare
zz = xx + yy*yy

These now have different shapes: shape(z) is (11,21) and shape(zz) is (21,11). In fact, they are the transpose of each other, z == zz.T. Once you realize this, it all becomes clearer - going from interp2d to RectBivariateSpline swapped the expected axes. Pick one instantiation of the splines (I've opted for the newer ones), and you have picked a particular set of axes to keep clear in your head. To mix them together, a simple transpose will work as well, but can get to be a headache when you go back through your code a month or more from now.

Upvotes: 1

Related Questions