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

Estoy teniendo problemas para obtener scipy.interpolate.UnivariateSpline para usar cualquier suavizado al interpolar. Basándome 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:

Splines sin alisar.

He intentado esto con un rango de valores de s (0.01, 0.1, 1, 2, 5, 50), así como ponderaciones explícitas, establecidas en la misma cosa (1.0) o aleatorias. Todavía no puedo obtener suavizado, y el número de nudos es siempre el mismo que el número de puntos de datos. En particular, estoy buscando valores atípicos como ese 4to punto (7990.4664106277542, 5851.6866463790966) para ser suavizado.

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

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

La documentación para 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 

A partir de esto, se puede deducir que los valores "razonables" para el suavizado, si no se pasan ponderaciones explícitas, 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 sensibles para s dependen, por supuesto, 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.

La respuesta de @Zhenya de establecer nudos manualmente entre puntos de datos fue demasiado aproximada para ofrecer buenos resultados en datos ruidosos sin ser selectiva 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 conteo de conglomerados y parece hacer un buen trabajo de suavizado (de hecho, muy suave).

 # 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() 

introduzca la descripción de la imagen aquí

Si bien no conozco ninguna biblioteca que lo haga por usted de la mano, probaría un poco más de bricolaje: empezaría por hacer una spline con nudos entre los puntos de datos sin procesar, tanto en x como en y En su ejemplo particular, tener un solo nudo entre los puntos 4 y 5 debería hacer el truco, ya que eliminaría la enorme derivada en torno a x=8000 .

Tuve problemas para ejecutar la respuesta de BigChef, aquí hay una variación que funciona en Python 3.6:

 # Imports 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 import numpy 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(key=lambda li: li[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, int(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.show()