Two dimensional cubic spline interpolation - does the order matter?

Hi,

I am doing two-dimensional interpolation, using cubic splines. I have an M-by-N matrix of function values: z(m,n), and a vector of x-values (of length M) and a vector of y values (of length N). The task is to calculate the function value at (x*,y*), where this point is somewhere in the interior of the two-dimensional grid. In order to have a twice-differentiable interpolating function, I am using cubic splines.

My understanding is that Matlab's interp2 procedure does something that is commonly referred to as bicubic splines. Essentially, this is what it does. First, it performs M one-dimensional splines across the rows of the table, and for each of those splines, it calculates the function values at (xm,y*), for m=1,2,...,M. Then, it does one additional one-dimensional spline down the newly created column {(xm,y*)}.

My question is this: does it matter in which order do I choose the dimensions? In other words, if the first step was N one-dimensional splines across the columns of the table, followed by a one-dimensional spline across the newly created row - would I still get the same result?

My conjecture is that the answer is yes. When I do this numerically, the order does not seem to matter, I get practically the same answer. Also, it is straightforward to prove that if the interpolation was linear (bilinear), then again, the order would not matter. For cubic interpolation, the answer is less obvious. I tried to see what the literature says, but I could not find anything definitive.

Also, if the order does not matter in two-dimensions, it should not matter in three or more dimensions either, right?

Any thoughts?

Thank you,

Zsolt

Answer 1

[Matt J]

No, the order should not matter for any separable interpolator in any dimension, except for small floating point noise. Floating point errors do depend on the order of the computations.

[Matt J]

Zsolt Commented:

Matt J,

Thanks for the quick response.

In response to your answer, I would like to point out that the second interpolation gets its inputs from the first M interpolations. In other words, the functional values over which one interpolates in the second step are calculated in the first step. To say it yet another way, the second step is dependent on the first.

But perhaps this isn't what separability of interpolation is about. I have to say that I am not familiar with the term. But aside from that, in light of what I just said above, do you think your answer still holds?

Zsolt

[Matt J]

Zsolt,

One way to see the order independence is to note that 1D interpolation is a linear operation on the column vector of data x being interpolated. Thus, it can be represented as a matrix multiplication

 result = K*x

This is true of all interpolation weighting schemes, linear, cubic, spline, or whatever.

As you pointed out, 2D interpolation of a matrix X correponds to 1D interpolation across the rows of X, hence a linear transformation K1*X of all the columns of X, followed by a 1D interpolation across the columns, hence

result = (K1*X)*K2.'

But because of the associativity of matrix multiplication, this can also be expressed

result = K1*(X*K2.')

equivalent to interpolating across the columns first, then the rows. So, the order doesn't matter.

[Zsolt Macskasi]

Aha, that is a very good point. Thanks a lot.