python - interpolate - Interpolar valores de NaN en una matriz numpy
scipy interpolate (8)
¿Existe una manera rápida de reemplazar todos los valores NaN en una matriz numpy con (por ejemplo) los valores linealmente interpolados?
Por ejemplo,
[1 1 1 nan nan 2 2 nan 0]
se convertiría en
[1 1 1 1.3 1.6 2 2 1 0]
Necesitaba un enfoque que también llene los NaN al comienzo del final de los datos, lo que parece no ser la respuesta principal.
La función que ideé utiliza una regresión lineal para completar los NaN. Esto supera mi problema:
import numpy as np
def linearly_interpolate_nans(y):
# Fit a linear regression to the non-nan y values
# Create X matrix for linreg with an intercept and an index
X = np.vstack((np.ones(len(y)), np.arange(len(y))))
# Get the non-NaN values of X and y
X_fit = X[:, ~np.isnan(y)]
y_fit = y[~np.isnan(y)].reshape(-1, 1)
# Estimate the coefficients of the linear regression
beta = np.linalg.lstsq(X_fit.T, y_fit)[0]
# Fill in all the nan values using the predicted coefficients
y.flat[np.isnan(y)] = np.dot(X[:, np.isnan(y)].T, beta)
return y
Aquí hay un caso de uso de ejemplo:
# Make an array according to some linear function
y = np.arange(12) * 1.5 + 10.
# First and last value are NaN
y[0] = np.nan
y[-1] = np.nan
# 30% of other values are NaN
for i in range(len(y)):
if np.random.rand() > 0.7:
y[i] = np.nan
# NaN''s are filled in!
print y
print linearly_interpolate_nans(y)
O basándose en la respuesta de Winston
def pad(data):
bad_indexes = np.isnan(data)
good_indexes = np.logical_not(bad_indexes)
good_data = data[good_indexes]
interpolated = np.interp(bad_indexes.nonzero()[0], good_indexes.nonzero()[0], good_data)
data[bad_indexes] = interpolated
return data
A = np.array([[1, 20, 300],
[nan, nan, nan],
[3, 40, 500]])
A = np.apply_along_axis(pad, 0, A)
print A
Resultado
[[ 1. 20. 300.]
[ 2. 30. 400.]
[ 3. 40. 500.]]
Para datos bidimensionales, la griddata
de SciPy funciona bastante bien para mí:
>>> import numpy as np
>>> from scipy.interpolate import griddata
>>>
>>> # SETUP
>>> a = np.arange(25).reshape((5, 5)).astype(float)
>>> a
array([[ 0., 1., 2., 3., 4.],
[ 5., 6., 7., 8., 9.],
[ 10., 11., 12., 13., 14.],
[ 15., 16., 17., 18., 19.],
[ 20., 21., 22., 23., 24.]])
>>> a[np.random.randint(2, size=(5, 5)).astype(bool)] = np.NaN
>>> a
array([[ nan, nan, nan, 3., 4.],
[ nan, 6., 7., nan, nan],
[ 10., nan, nan, 13., nan],
[ 15., 16., 17., nan, 19.],
[ nan, nan, 22., 23., nan]])
>>>
>>> # THE INTERPOLATION
>>> x, y = np.indices(a.shape)
>>> interp = np.array(a)
>>> interp[np.isnan(interp)] = griddata(
... (x[~np.isnan(a)], y[~np.isnan(a)]), # points we know
... a[~np.isnan(a)], # values we know
... (x[np.isnan(a)], y[np.isnan(a)])) # points to interpolate
>>> interp
array([[ nan, nan, nan, 3., 4.],
[ nan, 6., 7., 8., 9.],
[ 10., 11., 12., 13., 14.],
[ 15., 16., 17., 18., 19.],
[ nan, nan, 22., 23., nan]])
Lo estoy usando en imágenes 3D, operando en rebanadas 2D (4000 rebanadas de 350x350). Toda la operación aún demora alrededor de una hora: /
Primero definamos una función de ayuda simple para que sea más sencillo manejar índices e índices lógicos de NaNs :
import numpy as np
def nan_helper(y):
"""Helper to handle indices and logical indices of NaNs.
Input:
- y, 1d numpy array with possible NaNs
Output:
- nans, logical indices of NaNs
- index, a function, with signature indices= index(logical_indices),
to convert logical indices of NaNs to ''equivalent'' indices
Example:
>>> # linear interpolation of NaNs
>>> nans, x= nan_helper(y)
>>> y[nans]= np.interp(x(nans), x(~nans), y[~nans])
"""
return np.isnan(y), lambda z: z.nonzero()[0]
Ahora el nan_helper(.)
Ahora se puede utilizar como:
>>> y= array([1, 1, 1, NaN, NaN, 2, 2, NaN, 0])
>>>
>>> nans, x= nan_helper(y)
>>> y[nans]= np.interp(x(nans), x(~nans), y[~nans])
>>>
>>> print y.round(2)
[ 1. 1. 1. 1.33 1.67 2. 2. 1. 0. ]
---
Aunque puede parecer primero un poco exagerado especificar una función separada para hacer cosas como esta:
>>> nans, x= np.isnan(y), lambda z: z.nonzero()[0]
eventualmente pagará dividendos.
Por lo tanto, cada vez que trabaje con datos relacionados con NaN, solo debe encapsular todas las funcionalidades (nuevas relacionadas con NaN) necesarias, bajo alguna función auxiliar específica (s). Su código base será más coherente y legible, ya que sigue expresiones idiomáticas fácilmente comprensibles.
La interpolación, de hecho, es un buen contexto para ver cómo se hace el manejo de NaN, pero también se utilizan técnicas similares en otros contextos.
Puede ser más fácil cambiar la forma en que se generan los datos en primer lugar, pero si no:
bad_indexes = np.isnan(data)
Crea una matriz booleana que indique dónde están los nans
good_indexes = np.logical_not(bad_indexes)
Crear una matriz booleana que indique dónde se encuentra el área de valores válidos
good_data = data[good_indexes]
Una versión restringida de los datos originales con exclusión de los nans
interpolated = np.interp(bad_indexes.nonzero(), good_indexes.nonzero(), good_data)
Ejecute todos los índices incorrectos a través de la interpolación
data[bad_indexes] = interpolated
Reemplace los datos originales con los valores interpolados.
Se me ocurrió este código:
import numpy as np
nan = np.nan
A = np.array([1, nan, nan, 2, 2, nan, 0])
ok = -np.isnan(A)
xp = ok.ravel().nonzero()[0]
fp = A[-np.isnan(A)]
x = np.isnan(A).ravel().nonzero()[0]
A[np.isnan(A)] = np.interp(x, xp, fp)
print A
Imprime
[ 1. 1.33333333 1.66666667 2. 2. 1. 0. ]
Simplemente use numpy logical y there where statement para aplicar una interpolación 1D.
import numpy as np
from scipy import interpolate
def fill_nan(A):
''''''
interpolate to fill nan values
''''''
inds = np.arange(A.shape[0])
good = np.where(np.isfinite(A))
f = interpolate.interp1d(inds[good], A[good],bounds_error=False)
B = np.where(np.isfinite(A),A,f(inds))
return B
Sobre la base de la respuesta de , modifiqué su código para convertir también listas que constan solo de NaN
a una lista de ceros:
def fill_nan(A):
''''''
interpolate to fill nan values
''''''
inds = np.arange(A.shape[0])
good = np.where(np.isfinite(A))
if len(good[0]) == 0:
return np.nan_to_num(A)
f = interp1d(inds[good], A[good], bounds_error=False)
B = np.where(np.isfinite(A), A, f(inds))
return B
Además simple, espero que sea de utilidad para alguien.