MATLAB: Interp1 or its alternatives, working for matrices without loop

interpolationvectorization

Hi All,

I have matrices X and Y, with size(x) = [30,30000], and size(y) = [30,30000]. Columnar data in x are different between different columns, and so is y. xi is interpolation vector, with size(xi) = [60, 1]. For linear interpolation, I know I can do:

for i=1: 30000

output(:,i) = interp1(x(:,i), y(:,i), xi);

end

But it is slow to go through the 30000 loops. What is worse, x and y are daily data, and I need to do it annually. Therefore, using loop is a way too slow. I found that shared scripts like ScaleTime exists, which however do not accept x with different columnar values.

Question 1:

Is there a vecterize method to do the interpolation (linear)?

Question 2:

What about vecterized interpolation using cubic or pchip?

Many thanks,

John

Best Answer

So you're going over column by column, creating columns of 60 out of columns of 30. Assuming your xi are just exactly between the elements and exactly on the elements, you can just assume the "in between" interpolations are just the average of the two elements on either side. Thus you can just use a sliding mean filter. Use conv2(y(:,i), [.5; .5]) to give you the "in between" values, then just interleave it with the existing values. Exercise for the reader: interleave the convolution output with the original matrix. (It's like 3 lines of code or so.) Let me know if you can't figure it out and need help.

If your xi are spaced out at non-periodic oddball locations then you can't use conv2.

Related Solutions

MATLAB: Fast interp1 with ‘spline’

Well, the simple answer is to not use interp1. Just use spline and then ppval.

The problem is that for a spline interpolant, every time you call interp1, it must effectively call spline inside on your data. By calling spline once, this makes it more efficient.

For example...

x = linspace(-1,1,100);
y = exp(x);
xev = rand(1,1000)*2 - 1;
S = spline(x,y);
timeit(@() interp1(x,y,xev,'spline'))
ans =
    0.0003005
timeit(@() ppval(S,xev))
ans =
   0.00025177

Not as big of a difference as I thought it might be, but some.

You MIGHT also gain some time if you are willing to pre-interpolate the function to a finer interval, so that then you could do linear interpolation. Since linear interpolant will be faster to do, this should see some gain too.

x = linspace(-1,1,10);
y = exp(x);
xfine = linspace(-1,1,1000);
yfine = interp1(x,y,xfine,'spline');
xev = rand(1,1000)*2 - 1;
timeit(@() interp1(x,y,xev,'spline'))
ans =
   0.00028546
timeit(@() interp1(xfine,yfine,xev))
ans =
   0.00022313

You can also see some gain by using pchip instead of spline, as that is a faster way to build the spline, though sometimes not quite as smooth.

timeit(@() interp1(x,y,xev,'pship'))
ans =
   0.00026064

Finally I recall seeing some tools on the file exchange that tried to give a speedup for interp1, but they were mostly for linear interpolation.

Finally, while you could certainly use a mexfile to improve the time, that presumes that your skills are sufficient to write a spline interpolant, and to do so efficiently in C. I'm afraid that compiled MATLAB code would gain you nothing here.

MATLAB: Index exceeds the number of array elements (11)

h = 0.01; % time step
dxdt = @(x,y) y;
dydt = @(x,y) x-x^3;
m = 100 ;
x = zeros(m,1) ;
y = zeros(m,1) ; 
x(1)=1;
y(1)=1;
for i=1:1:m-1
    x(i+1) = x(i) + h*dxdt(x(i), y(i));
    y(i+1) = y(i) + h*dydt(x(i), y(i));
end
% Data Matrices X
 X = [x y x.^2 x.*y x.^2.*y x.*y.^2 y.^2 y.^3 y.*x.^3 x.*y.^3 y.^3 x.^4 y.*x.^4 y.^4];

Best Answer

Related Solutions

MATLAB: Fast interp1 with ‘spline’

MATLAB: Index exceeds the number of array elements (11)

Related Question