Is it possible to do multidimensional interpolation on a set of scattered data?

I’m completely stumped by a seemingly simple problem, and have been for a couple of years. Can any of you smart matlabbers help?

I want to do a multidimensional interpolation of a set of data. So, use interpn, right? (I tried using griddatan and I’ve been waiting for hours while it’s “busy”, trying just one interpolation (not the thousands I need).)

My data are stored in a large 2D matrix. They are not completely ordered, like in the example below. My actual problem has 9 dimensions of independent data, but in the example below, assume the first four columns are independent variables and the fifth column is the dependent variable that I seek an interpolation value for.

To use interpn, I need to supply the function with a multidimensional array – with the number of dimensions matching the number of independent variables (right?). So for this example, I’d need a 4D array. Example:

X1  X2  X3  X4  Y
24  1  8  30  0.73
24  2  9  40  0.04
24  3  9  50  0.18
24  4  10  60  0.88
24  1  11  70  0.18
24  2  12  80  0.94
24  3  8  90  0.36
24  4  9  100  0.58
24  1  10  110  0.86
24  2  10  120  0.49
26  3  10  130  0.03
26  4  11  140  0.90
26  1  12  150  0.03
26  2  8  160  0.12
26  3  9  170  0.85
26  4  10  180  0.35
26  1  11  190  0.51
26  2  12  200  0.05
26  3  12  210  0.49
26  4  13  220  0.45

How do I put them in a 4D array? Ignoring X4 for the moment, and looking at a slice of the array where X3 = 8, we would have this array:

X1\X2  1  2  3  4
24  0.73    0.36  
26    0.12

There are empty cells. Is that OK? Otherwise how could we make a multidimensional array from scattered data?

Building that array by hand is impractical with a large array. Is there some automated way to do it? Am I making this harder than it should be? Or is it so hard that I should do something else?

Thanks!!!

[Image Analyst]

Can you attach your data file?

[Canoe Commuter]

Sure. Here's an Excel file, with a header showing the variable names (X1 to X9 and Y) and a column with unique IDs for each of the 3217 points.

There's also a .m file with the data entered as a 3217 x 10 array (can't upload a .mat file to this site). Columns 1 to 9 are the independent variables (X1 to X9) and column 10 is the dependent variable (Y).

Thanks for looking at it!

Answer 1

[Canoe Commuter]

After, some consideration I eventually decided to take a different approach. Instead of trying to do a multidimensional interpolation within my large data set, I decided to create neural network models with the data set, and model any points that I need. There are advantages and disadvantages to this approach, but both approaches introduce a small amount of error, and the neural net approach is doable, whereas I never succeeded in the interpolation approach, despite the excellent guidance from Jeremy, Matt J, and Kelley Kearney.

If you are trying to figure out multidimensional interpolation, there are some good points and tools below.

[Matthias]

Hi, I am working on a similar problem. Do you perhaps have some literature suggestion on the neural network approach that you have used?

Answer 2

[Kelly Kearney]

Does this function do what approximately what you're looking for: vec2grid.m? Right now, I've only written it to deal with 2D and 3D arrays, but it could easily be expanded to include an arbitrary number of dimensions.

data = [...
24  1  8  30  0.73
24  2  9  40  0.04
24  3  9  50  0.18
24  4  10  60  0.88
24  1  11  70  0.18
24  2  12  80  0.94
24  3  8  90  0.36
24  4  9  100  0.58
24  1  10  110  0.86
24  2  10  120  0.49
26  3  10  130  0.03
26  4  11  140  0.90
26  1  12  150  0.03
26  2  8  160  0.12
26  3  9  170  0.85
26  4  10  180  0.35
26  1  11  190  0.51
26  2  12  200  0.05
26  3  12  210  0.49
26  4  13  220  0.45];

 [x1, x2, x3, y] = vec2grid(data(:,1), data(:,2), data(:,3), data(:,5))

 x1 =

    24
    26

 x2 =

     1
     2
     3
     4

 x3 =

     8
     9
    10
    11
    12
    13

 y(:,:,1) =

         0.73          NaN         0.36          NaN
          NaN         0.04         0.18         0.58
         0.86         0.49          NaN         0.88
         0.18          NaN          NaN          NaN
          NaN         0.94          NaN          NaN
          NaN          NaN          NaN          NaN

 y(:,:,2) =

          NaN         0.12          NaN          NaN
          NaN          NaN         0.85          NaN
          NaN          NaN         0.03         0.35
         0.51          NaN          NaN          0.9
         0.03         0.05         0.49          NaN
          NaN          NaN          NaN         0.45

[Cedric Wannaz]

I was using blah as well !! ;)

[Kelly Kearney]

Yes, ever since Matlab introduced all the mlint stuff, my code has been littered with orange-highlighted blah's. :-) But too many potential collaborators use older versions of Matlab for me to switch to ~.

[Canoe Commuter]

Thanks, Kelly! That's a cool function. I have modified it to deal with a 9D array, as you suggested. I do this in the hopes that if I have a 9D array, I can use interpn, and that it will execute fairly quickly. Maybe very optimistic.

Anyway, the problem I run into now is that I get an error "Error using sub2ind. Out of range subscript." I think this is caused by the fact that not all of my independent data are integers, as they were in the simple example that I provide above. Am I back to square one?

[Kelly Kearney]

No, non-integers shouldn't be a problem. I've pasted a modified version of my function below.

However , I don't suggest you test it with your data... you'll probably freeze your computer. Now that I've taken a look at the Excel file you attached above, I really don't think this is what you want to do, because your resulting matrix will be huge (~270 GB double array), and very, very, sparse (only 3217 non-zero values in a 4x5x9x73x29x12x7x21x54 array, and keep in mind that Matlab's sparse matrix data type doesn't support more than 2 dimensions).

So perhaps we should start over on this problem... what exactly is the end product that you want? Is it really the full 9-dimensional behemoth that I describe above, with all values filled in based on that tiny sampling of 3217 values? Or something else?

***

And for future reference on the n-dimensional regridding:

 function varargout = vec2grid(varargin)
 %VEC2GRID Grid ND data without interpolation
 %
 % [xg, yg] = vec2grid(x1, x2, ..., y);
 % [xg, yg] = vec2grid(xy);
 % [xg1, xg2, ... yg] = vec2grid(...);
 %
 % This function reformats a list of n-dimensional datpoints to a grid,
 % similar to the results of the griddata function but without any
 % interpolation. It reshapes the data using all unique values in the 1st,
 % 2nd, 3rd, etc. as rows, columns, pages, etc. of the final data grid.     
 %
 % Input variables:
 %
 %   x#:     coordinates of data points, vectors.  Each input corresponds to
 %           one dimension of the data 
 %
 %   y:      y coordinates of data points, vector the same size as x(s)
 %
 %   xy:     npoint x ndim+1 array, where each row represents one point.
 %           The first ndim columns are the coordinates of the points, and
 %           the final column holds the data value at that coordinate
 %
 % Output variables:
 %
 %   xg#:    npoint x 1 cell array of unique values of the #th dimension of
 %           yg 
 %
 %   yg:     ndim-dimensional array, y values rearranged into a grid
 %
 %   xg:     1 x ndim cell array holding all xg# vectors
 %
 % Example:data
 %
 % [x,y] = meshgrid(1:5,0:2:6);
 % x = x(:); y = y(:);
 % holes = rand(size(x)) < .2;
 % x(holes) = []; y(holes) = [];
 % v = 1:length(x);
 % [xg,yg,vg] = vec2grid(x,y,v)
 % xg =
 %      1
 %      2
 %      3
 %      4
 %      5
 % yg =
 %      0
 %      2
 %      4
 %      6
 % vg =
 %      1   NaN   NaN   NaN    13
 %      2     5     8    11    14
 %      3     6     9   NaN    15
 %      4     7    10    12    16

 % Copyright 2008 Kelly Kearney

 %-------------------------
 % Check input
 %-------------------------

 if nargin == 1
     x = num2cell(varargin{1}(:,1:end-1), 1);
     y = varargin{1}(:,end);
 else
     if ~all(cellfun(@isvector, varargin)) || any(diff(cellfun(@length, varargin)))
         error('Inputs must be vectors of the same length');
     end
     x = varargin(1:end-1);
     y = varargin{end};
 end

 % error(nargchk(3, 4, nargin));

 %-------------------------
 % Create grid
 %-------------------------

 ndim = length(x);

 [xunique, blah, xidx] = cellfun(@unique, x, 'uni', 0);
 nx = cellfun(@max, xidx);

 ynew = nan(nx);
 ynew(sub2ind(nx, xidx{:})) = y;

 if nargout == (ndim+1)
     varargout = cell(ndim+1,1);
     tmp = [xunique ynew];
     [varargout{:}] = deal(tmp{:});
 elseif nargout == 2
     varargout{1} = xunique;
     varargout{2} = ynew;
 end

[Canoe Commuter]

Thanks again, Kelly!

What I'm after is a way to interpolate my data sets. Each set has nine independent variables. I'd like to be able to arbitrarily choose a 9D location, and find an interpolated value at that location. It seems like interpn is intended for this purpose, but the input data need to be in a 9D array.

I'm starting to think that this is too computationally intensive. Part of the rationale for doing this was an assumption that interpolation is faster than modeling new data, but perhaps this assumption fails when the dimensionality increases. Can confirm this?

[Kelly Kearney]

The interpn function isn't going to help you here; it will incorporate the NaNs into its interpolation, so your resulting interpolated values are going to be NaN for most of the coordinate space.

Do you want a nearest-neighbor interpolation or something else (linear, etc)? If the former, there could be several shortcuts around using griddatan.

If the latter, it might help to manipulate your data a little, and remove some of the dimensions you're interpolating over. For example, it looks like your first 3 variables are pretty grid-like, while the remaining dimensions are more scattered. Perhaps you could aggregate your data into dim1 x dim2 x dim3 groups, then interpolate separately within each of those "slices."

Whether any of this will be faster than rerunning your model is entirely dependent on that model.

Answer 3

[Matt J]

So, use interpn, right?

No, not if your data is scattered. And rather than griddatan, scatteredInterpolant() is probably what would be recommended as the latest and greatest, if you have a sufficiently recent MATLAB release.

If your data can always be viewed as gridded data with missing elements, and the idea is to to fill the missing data with something, you could try this FEX file

http://www.mathworks.com/matlabcentral/fileexchange/27994-inpaint-over-missing-data-in-1-d-2-d-3-d-n-d-arrays

I haven't used it myself, but it's getting favorable reviews. However, you would probably need a lot less missing data for this to work than in your example. In the example, your data essentially forms a line in 4-dimensional space. That's not a lot to extrapolate from.

[Canoe Commuter]

Thanks, Matt. I'm using R2013a, which includes scatteredInterpolant, but it looks like that can only be used for 2D or 3D scattered data sets. I haven't yet checked out the FEX file.

Based on the wisdom of the commenters in this thread, I am currently working on reshaping my problem to be several 4D problems (dimensions 4 to 9 vary one-at-a-time, so this makes sense). As Kelly points out, dimensions 1 to 3 are not scattered. I'll update when I make progress...

Answer 4

[Jeremy]

For 2 dimensional data (surface of scattered points), I've found gridfit to be an all-star. Not sure how well that will translate to your 9 dimensional problem...

[Canoe Commuter]

Thanks, Jeremy. gridfit looks like a great tool, but as you suggested, it's only for 2D data sets.