3

I have a 2D array with some NaN values. I would like to inpaint (interpolate) those values using the locations where I have data. The array looks like the one below.

If possible I would like to do the interpolation so that, as I move away from non-NaN values, I get increasingly closer to the value 0.

How can I do this?

I read about gridddata, but it seems to be designed to work with unstructured N-dim data. I also read the answers in other threads, but I think their starting point is different. 

array([[        nan,         nan,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [        nan,         nan,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [        nan,         nan,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [        nan,         nan,         nan,         nan,         nan,
                nan,         nan,         nan,         nan,         nan],
       [ 1.        ,  0.        ,  1.        ,  0.        ,  0.25      ,
                nan,  0.        ,         nan,         nan,         nan],
       [        nan,  0.        ,         nan,  0.25      ,  0.66666667,
         0.25      ,  0.66666667,  0.        ,  1.        ,         nan],
       [ 0.        ,  0.5       ,  0.66666667,  0.8       ,  0.66666667,
         0.8       ,  0.5       ,  0.83333333,         nan,         nan],
       [ 0.625     ,  0.5625    ,  0.9       ,  0.8       ,  0.8       ,
         0.83333333,  0.57142857,  0.66666667,  0.5       ,         nan],
       [        nan,  1.        ,  0.71428571,  0.85714286,  1.        ,
         1.        ,  1.        ,         nan,         nan,         nan],
       [        nan,         nan,         nan,         nan,  1.        ,
         1.        ,         nan,         nan,         nan,         nan]])
Community
  • 1
  • 1
Amelio Vazquez-Reina
  • 91,494
  • 132
  • 359
  • 564

1 Answers1

3

There are dozens of possible approaches based on what kind of interpolation technique you would like to use. In fact, as your data is rather surrounded by NaNs I would rather think about it as a function smoothing then interpolating. If you want to get closer to zero the more away you are from the not NaNs in terms of euclidean distance on your 2d map I would suggest something like:

  1. Consider each not NaN data point X[i,j] as a Gaussian centered in [i,j], with variance=1, scaled so its pdf( [i,j] ) = X[i,j], so f_ij( [a,b] ) = X[i,j] * exp( -|| [a,b] - [i,j] ||^2/2 ).
  2. For each NaN data point X[a,b] set X[a,b] = sum( f_ij( [a,b] ) ) where sumation is performed over all [i,j] indices of not NaN data points

As a result you get something like a "density estimation", and by changing the variance (which I suggested to use =1) you can modify the "speed of vanishing" the values.

So the code would be just a one loop over all NaNs, and for each of them you loop through all not NaNs and sum the gaussians values.

It would sth like this:

nans    = np.array( np.where(  np.isnan(X) ) ).T
notnans = np.array( np.where( ~np.isnan(X) ) ).T
for p in nans:
    X[p[0],p[1]] = sum( X[q[0],q[1]]*np.exp(-(sum((p-q)**2))/2) for q in notnans )
lejlot
  • 64,777
  • 8
  • 131
  • 164
  • 1
    If you do choose to go down this route, you should know that Scipy has a class for radial basis functions - see [here](http://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.Rbf.html), and [this example here](http://wiki.scipy.org/Cookbook/RadialBasisFunctions) – ali_m Oct 06 '13 at 16:27
  • IS there ready functions or standard short code for ordinary kriging interpolation method in the same way? THanks – Spider Apr 23 '14 at 22:29
  • @ali_m: in that example x, y, z need to have the same size. In this case x, y are 1d arrays while z is 2d array! – hoang tran Feb 20 '15 at 17:13
  • @hoangtran Well in this case you would just use the non-nan elements in the 2D array of z-values as your 1D vector `z`, and `x` and `y` would be their corresponding column/row positions. – ali_m Feb 20 '15 at 17:23