np.polyfit() — Ajuste de curvas con NumPy Polyfit

La función polyfit() acepta tres valores de entrada diferentes: x, y y el grado polinómico. Los argumentos x e y corresponden a los valores de los puntos de datos que queremos ajustar, en los ejes x e y, respectivamente. El tercer parámetro especifica el grado de nuestra función polinómica. Por ejemplo, para obtener un ajuste lineal, utiliza el grado 1.

¿Qué es el ajuste de curvas?

El ajuste de curvas consiste en construir una función matemática que sea capaz de ajustarse a unos puntos de datos concretos. La mayoría de las veces, la ecuación de ajuste está sujeta a restricciones; además, también es posible hacer conjeturas iniciales para proporcionar puntos de partida útiles para la estimación de los parámetros de ajuste, este último procedimiento tiene la ventaja de reducir el trabajo computacional. En este artículo exploraremos la función de NumPy polyfit(), que permite crear funciones de ajuste polinómico de una manera muy simple e inmediata.

Regresión lineal

El tipo de ajuste más simple es el ajuste lineal (una función polinómica de primer grado), en el que los puntos de datos se ajustan utilizando una línea recta. La ecuación simplificada de una línea recta es:

y = mx + n

Donde “m” se llama pendiente y “n” ordenada en el origen o término independiente (“intercepción” del inglés literal). Cuando aplicamos un ajuste lineal, básicamente estamos buscando los valores de los parámetros “m” y “n” que producen el mejor ajuste para nuestros puntos de datos. En Numpy, la función polyfit() es una herramienta muy intuitiva y poderosa para ajustar puntos de datos; veamos cómo ajustar una serie aleatoria de puntos de datos con una línea recta.

En el siguiente ejemplo, queremos aplicar un ajuste lineal a algunos puntos de datos, descritos por los arrays x e y. La función polyfit() acepta tres valores de entrada diferentes: x, y y el grado polinómico. Mientras que x e y corresponden a los valores de los puntos de datos que queremos ajustar, en los ejes x e y, respectivamente, el tercer parámetro especifica el grado de nuestra función polinómica. Dado que queremos un ajuste lineal, vamos a especificar un grado igual a 1. La salida de la función polyfit() será una lista que contiene los parámetros de ajuste: el primero es el que en la función se multiplica por el término de mayor grado; los demás siguen ese orden decreciente.

import numpy as np
from numpy import random  #it will be useful for generating some random noise (on purpose) in the data points that we want to fit
import matplotlib.pyplot as plt  #for plotting the data

#---LINEAR FIT----

#generate the x array
x = np.linspace(0,60,60) # generate an array of 60 equally space points

#generate the y array exploiting the random.randint() function to introduce some random noise
y = np.array([random.randint(i-2, i+2) for i in x]) #each element is a random number with value between +-2 the respective x axis value

#Applying a linear fit with .polyfit()
fit = np.polyfit(x,y,1)
ang_coeff = fit[0]
intercept = fit[1]
fit_eq = ang_coeff*x + intercept  #obtaining the y axis values for the fitting function

#Plotting the data
fig = plt.figure()
ax = fig.subplots()
ax.plot(x, fit_eq,color = 'r', alpha = 0.5, label = 'Linear fit')
ax.scatter(x,y,s = 5, color = 'b', label = 'Data points') #Original data points
ax.set_title('Linear fit example')
ax.legend()
plt.show()

Como ya hemos dicho, la variable fit contendrá los parámetros de ajuste. El primero es la pendiente (coeficiente de grado 1), el último el término independiente. En este punto, para representar gráficamente nuestro ajuste, tenemos que elaborar los valores del eje y a partir de los parámetros obtenidos, utilizando los valores originales del eje x. En el ejemplo, este paso se describe mediante la definición de la variable fit_eq. Lo último que queda es representar los datos y la ecuación de ajuste. El resultado es:

Ajuste polinómico de segundo grado

En este segundo ejemplo, vamos a crear un ajuste polinómico de segundo grado. Las funciones polinómicas de este tipo describen una parábola en el plano xy; su ecuación general es:

y = ax2 + bx + c

donde a, b y c son los parámetros de la ecuación que estimamos al generar una función de ajuste. Los puntos de datos que vamos a ajustar en este ejemplo representan la trayectoria de un objeto que ha sido lanzado desde una altura desconocida. Aprovechando la función polyfit(), ajustaremos la trayectoria del objeto que cae y obtendremos también una estimación de su velocidad inicial en la dirección x, v0.

#-----POLYNOMIAL FIT----
x = np.array([1.2,2.5,3.4,4.0,5.4,6.1,7.2,8.1,9.0,10.1,11.2,12.3,13.4,14.1,15.0]) # x coordinates
y = np.array([24.8,24.5,24.0,23.3,22.4,21.3,20.0,18.5,16.8,14.9,12.8,10.5,8.0,5.3,2.4]) # y coordinates
fit = np.polyfit(x, y, 2)
a = fit[0]
b = fit[1]
c = fit[2]
fit_equation = a * np.square(x) + b * x + c

#Plotting
fig1 = plt.figure()
ax1 = fig1.subplots()
ax1.plot(x, fit_equation,color = 'r',alpha = 0.5, label = 'Polynomial fit')
ax1.scatter(x, y, s = 5, color = 'b', label = 'Data points')
ax1.set_title('Polynomial fit example')
ax1.legend()
plt.show()

Una vez inicializados los arrays x e y que definen la trayectoria del objeto, aplicamos la función polyfit(), esta vez insertando “2” como grado de la función de ajuste polinómico. Esto se debe a que la trayectoria de un objeto que cae puede describirse mediante un polinomio de segundo grado; en nuestro caso, la relación entre las coordenadas x e y viene dada por:

y = y0 – ½ (g/ v02)x2

donde y0 es la posición inicial (la altura a partir de la cual se ha lanzado el objeto), g la aceleración de la gravedad (9,81 m/s2) y v0 la velocidad inicial (m/s) en la dirección x (visita https://es.wikipedia.org/wiki/Ecuaciones_para_un_cuerpo_en_ca%C3%ADda_libre para obtener más detalles). Luego asignamos a las variables a, b y c el valor de los 3 parámetros de ajuste y definimos fit_equation, la ecuación polinómica que se trazará; el resultado es:

Si ahora imprimimos los tres parámetros de ajuste a, b y c, obtenemos los siguientes valores: a = -0.100 , b = 0.038, c = 24.92. En la ecuación que describe la trayectoria de un cuerpo que cae no existe el término b; como el ajuste es siempre una aproximación al resultado real, siempre obtendremos un valor para todos los parámetros; sin embargo, observaremos que el valor de nuestro término b es mucho menor que los demás y puede ser despreciado en cierto modo, al comparar nuestro ajuste con la ecuación que describe la física del problema. El término c representa la altura inicial (y0) mientras que el término a describe la cantidad – ½ (g/v02). Por lo tanto, la velocidad inicial v0 viene dada por:

v0=2-g2aDando como resultado el valor final de v0 = 6,979 m/s.