Sommaire
Ce tutoriel est en cours de rédaction / approbation.
Régression linéaire
Lorsque l’on a un jeu de données expérimentales permettant de relier deux grandeurs physiques, il est intéressant de trouver la loi qui fait le lien entre elles. Une première loi à tester est une loi affine de type \(y = a \cdot x + b\). On appelle alors cette opération une régression linéaire.
Il existe deux méthodes pour cela : la première utilisant la bibliothèque Numpy (voir tutoriel : xxx ) , l’autre utilisant la bibliothèque Scipy que nous allons aborder dans ce tutoriel.
Set de données expérimentales
Dans l’exemple suivant, on s’intéresse à une série de données acquise lors d’un travail expérimental autour de l’interféromètre de Michelson (voir TP LEnsE – Semestre 5 ).
import numpy as np n_anneaux = np.array([1,2,3,4,5,6,7,8,9,10]) dist = np.array([25.77,25.92,26.33,26.67,27.00,27.58,28.04,28.36,28.75,29.09])*10**(-3)-25.425*10**(-3) # Les incertitudes doivent bien sûr avoir la même dimension que les données ! in_n_anneaux = [0.5]*10 # Incertitude sur n_anneaux / génération d'un vecteur de 10 valeurs in_dist = [0.005*10**(-3)]*10 # Incertitude sur dist
Méthode avec SciPy.stats
Cette régression est réalisée à l’aide de la bibliothèque scipy et son module stats. La fonction linregress renvoie alors dans un vecteur la pente et l’origine de la fonction affine a.x + b résultante de la régression.
from scipy import stats regress = stats.linregress(dist,n_anneaux)
On peut extraire la pente de la droite obtenue en prenant le premier élément retournée par la fonction linregress et l’origine en 0 avec le deuxième élément.
pente = regress[0]
origine = regress[1]
# On peut alors afficher les valeurs de ces deux paramètres
print(f'Pente = {pente} et Origine = {origine}')
Pente = 2539.8932087671783 et Origine = 0.6081656799144115
Il est également possible de tracer la fonction (pente*x + origine).
y_reg = pente * dist + origine
Cette courbe peut alors être ajoutée au graphe avec les points bruts.
fig = plt.figure(figsize=(15,10))
plt.errorbar(dist, n_anneaux, xerr=in_dist, yerr=in_n_anneaux, linestyle="none", label="Données brutes")
plt.plot(dist, y_reg, label = "Régression linéaire sur les données brutes\n La pente vaut : " + str(round(pente,2)) + " m⁻¹")
plt.legend()
plt.ylabel("Nombre d'anneaux")
plt.xlabel("Épaisseur e en m")
plt.title("Nombre d'anneaux en fonction de l'épaisseur")
plt.show()
Ajustements (fit) avec d’autres fonctions
Il est possible de trouver une courbe qui s’ajuste à des données dont l’évolution est régie par une autre fonction que linéaire : polynome, logarithme… La fonction curve_fit du module optimize de la bibliothèque Scipy permet ce type d’ajustement.
L’exemple utilisé ici correspond à des données qui se répartissent théoriquement selon la loi suivante : \(f(x) = a \cdot e^{-b \cdot x} + c\)
Set de données
x = np.array([ 0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0]) y = np.array([2.8, 1.6, 1.5, 0.7, 0.9, 0.8, 0.4, 0.3, 0.5]) incertitude = np.array([0.1]*9)
Définition de la fonction
Avant toute chose, il est nécessaire de définir la fonction qui établit un lien entre les deux grandeurs physiques considérées.
# On définit la fonction théorique représentant l'évolution des données
def fonction_fit(x, a, b, c):
return a * np.exp(-b * x) + c
Ajustement par la fonction curve_fit
Il est nécessaire de fixer la valeur des paramètres recherchés à des valeurs initiales sous forme d’un vecteur. Ce vecteur deviendra un paramètre d’entrée de la fonction curve_fit.
param0 = np.array([0.0, 0.0, 0.0])
Il est maintenant possible de lancer l’algorithme de recherche des paramètres optimaux à l’aide de la fonction curve_fit du module optimize de la bibliothèque Scipy.
from scipy.optimize import curve_fit # curve_fit(fonction à fit,données,données,paramètres initiaux,incertitude) values, cov = curve_fit(fonction_fit, x, y, param0, incertitude)
Cette fonction donne en sortie deux vecteurs :
- values : un tableau contenant les valeurs estimées des paramètres de la fonction
- cov : la matrice de covariance, les incertitudes étant alors les racines carrés des éléments de la diagonale
On peut alors tracer la courbe régressée à partir des valeurs de a, b et c obtenues avec la fonction curve_fit.
x_trace_fit = np.linspace(0,4,100) # Pour tracer une fontion régressée "lissée" a = values[0] b = values[1] c = values[2] y_reg = fonction_fit(x_trace_fit, a, b, c)
On trace ensuite la figure avec les données brutes et la courbe obtenue par régression.
plt.figure(figsize=(15,10))
plt.errorbar(x, y, xerr=incertitude, yerr=incertitude, label="Mesures", linestyle="None")
plt.plot(x_trace_fit, y_reg, label="Régression sur les mesures")
plt.xlabel("x")
plt.ylabel("y")
plt.title("Régression pour des points expérimentaux par f=a*exp(-bx)+c")
plt.legend()
plt.show()
On obtient alors le graphique suivant :
