interpolate - spline python

scipy.interpolate.UnivariateSpline no suaviza independientemente de los parámetros

Question

Tengo problemas para obtener scipy.interpolate.UnivariateSpline para usar cualquier suavizado al interpolar. Basado en la página de la función así como en algunas publicaciones anteriores , creo que debería proporcionar suavizado con el parámetro s .

Aquí está mi código:

# Imports
import scipy
import pylab

# Set up and plot actual data
x = [0, 5024.2059124920379, 7933.1645067836089, 7990.4664106277542, 9879.9717114947653, 13738.60563208926, 15113.277958924193]
y = [0.0, 3072.5653360000988, 5477.2689107965398, 5851.6866463790966, 6056.3852496014106, 7895.2332350173638, 9154.2956175610598]
pylab.plot(x, y, "o", label="Actual")

# Plot estimates using splines with a range of degrees
for k in range(1, 4):
    mySpline = scipy.interpolate.UnivariateSpline(x=x, y=y, k=k, s=2)
    xi = range(0, 15100, 20)
    yi = mySpline(xi)
    pylab.plot(xi, yi, label="Predicted k=%d" % k)

# Show the plot
pylab.grid(True)
pylab.xticks(rotation=45)
pylab.legend( loc="lower right" )
pylab.show()

Aquí está el resultado:

Intenté esto con un rango de valores de s (0.01, 0.1, 1, 2, 5, 50), así como pesos explícitos, establecidos en la misma cosa (1.0) o aleatorizados. Todavía no puedo suavizar, y la cantidad de nudos es siempre la misma que la cantidad de puntos de datos. En particular, busco valores atípicos como el cuarto punto (7990.4664106277542, 5851.6866463790966) para suavizar.

¿Es porque no tengo suficientes datos? Si es así, ¿hay una función de spline similar o técnica de clúster que pueda aplicar para lograr el suavizado con estos pocos puntos de datos?

Answer 1

Si bien no estoy al tanto de ninguna biblioteca que lo haga por ti, intentaré un poco más de enfoque DIY: comenzaría haciendo una spline con nudos entre los puntos de datos brutos, tanto en x como en y En su ejemplo particular, tener un solo nudo entre los puntos 4 y 5 debería ser el truco, ya que eliminaría la enorme derivada en torno a x=8000 .

Answer 2

La respuesta de @Zhenya de establecer nudos manualmente entre los puntos de datos era demasiado difícil para ofrecer buenos resultados en datos ruidosos sin ser selectivo sobre cómo se aplica esta técnica. Sin embargo, inspirado por su sugerencia, he tenido éxito con la agrupación Mean-Shift del paquete scikit-learn. Realiza la autodeterminación del recuento del clúster y parece hacer un buen trabajo de suavizado (muy sencillo en realidad).

# Imports
import numpy
import pylab
import scipy
import sklearn.cluster

# Set up original data - note that it''s monotonically increasing by X value!
data = {}
data[''original''] = {}
data[''original''][''x''] = [0, 5024.2059124920379, 7933.1645067836089, 7990.4664106277542, 9879.9717114947653, 13738.60563208926, 15113.277958924193]
data[''original''][''y''] = [0.0, 3072.5653360000988, 5477.2689107965398, 5851.6866463790966, 6056.3852496014106, 7895.2332350173638, 9154.2956175610598]

# Cluster data, sort it and and save
inputNumpy = numpy.array([[data[''original''][''x''][i], data[''original''][''y''][i]] for i in range(0, len(data[''original''][''x'']))])
meanShift = sklearn.cluster.MeanShift()
meanShift.fit(inputNumpy)
clusteredData = [[pair[0], pair[1]] for pair in meanShift.cluster_centers_]
clusteredData.sort(lambda pair1, pair2: cmp(pair1[0],pair2[0]))
data[''clustered''] = {}
data[''clustered''][''x''] = [pair[0] for pair in clusteredData]
data[''clustered''][''y''] = [pair[1] for pair in clusteredData]

# Build a spline using the clustered data and predict
mySpline = scipy.interpolate.UnivariateSpline(x=data[''clustered''][''x''], y=data[''clustered''][''y''], k=1)
xi = range(0, round(max(data[''original''][''x'']), -3) + 3000, 20)
yi = mySpline(xi)

# Plot the datapoints
pylab.plot(data[''clustered''][''x''], data[''clustered''][''y''], "D", label="Datapoints (%s)" % ''clustered'')
pylab.plot(xi, yi, label="Predicted (%s)" %  ''clustered'')
pylab.plot(data[''original''][''x''], data[''original''][''y''], "o", label="Datapoints (%s)" % ''original'')

# Show the plot
pylab.grid(True)
pylab.xticks(rotation=45)
pylab.legend( loc="lower right" )
pylab.show()

Answer 3

Respuesta corta: debe elegir el valor para s más cuidado.

La documentación de UnivariateSpline establece que:

Positive smoothing factor used to choose the number of knots. Number of 
knots will be increased until the     smoothing condition is satisfied:
sum((w[i]*(y[i]-s(x[i])))**2,axis=0) <= s

De esto se puede deducir que los valores "razonables" para suavizar, si no se transfieren pesos explícitos, están alrededor de s = m * v donde m es el número de puntos de datos v la varianza de los datos. En este caso, s_good ~ 5e7 .

EDITAR : los valores razonables para s dependen, por supuesto, también del nivel de ruido en los datos. Los documentos parecen recomendar elegir s en el rango (m - sqrt(2*m)) * std**2 <= s <= (m + sqrt(2*m)) * std**2 donde std es el estándar desviación asociada con el "ruido" que desea suavizar.