Modélisation

Régression linéaire

La fonction linregress() du module scipy.stats retourne les paramètres de la régression linéaire d'un nuage de points.

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import linregress

x = np.array([0, 1.01, 2.02, 2.99, 3.98])    # Données en abscisse
y = np.array([10.32, 7.91, 6.13, 3.97, 2.25]) # Données en ordonnée


a, b, rho ,_ ,_ = linregress(x, y)    # Régression linéaire
print("a = ", a)                      # Affichage de coefficient directeur
print("b = ", b)                      # Affichage de l'ordonnée à l'origine
print("rho = ", rho)                  # Affichage du coefficient de corrélation

x_mod = np.linspace(0, 4, 50)         # Nouvelle abscisse pour la modélisation
y_mod = a*x_mod + b                   # Ordonnées de la fonction affine

plt.plot(x_mod, y_mod, '-')           # Tracé de la droite
plt.plot(x, y, '+')                   # Tracé des points et de la fonction affine
plt.title('Régression linéaire')      # Titre
plt.xlabel('x')                       # Etiquette en abscisse
plt.xlim(-1,5)                        # Echelle en abscisse
plt.ylabel('y')                       # Etiquette en ordonnée
plt.ylim(0, 12)                       # Echelle en ordonnée
plt.grid()
plt.show()                            # Affichage
alternate text
>>> %Run
a =  -2.0203420706406234
b =  10.156684141281247
rho =  -0.9986214342318739
  • a est le coefficient directeur.

  • b est l'ordonnée à l'origine.

  • rho est le coefficient de corrélation linéaire (pas significatif en sciences physiques).

Modélisation à partir d'un polynôme

La fonction polyfit() du paquet Numpy retourne les paramètres d'un modèle polynomial d'un nuage de points.

Cas d'une parabole

\[y = a\cdot x^2 + b \cdot x + c\]
import numpy as np
import matplotlib.pyplot as plt

# DONNEES
t = [0.0, 0.04, 0.08, 0.12, 0.16, 0.2, 0.24, 0.28, 0.32, 0.36, 0.4, 0.44, 0.48, 0.52, 0.56, 0.6, 0.64, 0.68, 0.72, 0.76, 0.8, 0.84, 0.88, 0.92]
x = [-0.953328037081172, -0.879995111151852, -0.799995555592592, -0.716662685218364, -0.636663129659105, -0.559996888914815, -0.479997333355555, -0.393331148166358, -0.313331592607099, -0.233332037047839, -0.149999166673611, -0.066666296299383, 0.013333259259877, 0.096666129634105, 0.179999000008333, 0.259998555567592, 0.343331425941821, 0.426664296316049, 0.506663851875308, 0.586663407434568, 0.663329648178858, 0.743329203738117, 0.819995444482407, 0.893328370411728]
y = [-0.046666407409568, 0.069999611114352, 0.166665740748457, 0.253331925937654, 0.326664851866975, 0.389997833351389, 0.433330925945988, 0.469997388910648, 0.486663962985494, 0.493330592615432, 0.489997277800463, 0.469997388910648, 0.433330925945988, 0.38666451853642, 0.323331537052006, 0.249998611122685, 0.156665796303549, 0.053333037039506, -0.063332981484414, -0.189998944453241, -0.333331481496913, -0.486663962985494, -0.65332970373395, -0.789995611147685]

# MODELISATION
a, b, c = np.polyfit(x, y, 2)  # Modélisation par un polynome de degré 2
print("a = ", a)
print("b = ", b)
print("c = ", c)

# CONSTRUCTION DU MODELE
x_mod = np.array(x)                # Création d'un tableau Numpy pour x
y_mod = a*x_mod**2 + b*x_mod + c   # Création d'un tableau Numpy pour y du modèle

# COURBES
plt.plot(x, y, "b+", label = "trajectoire")    # Courbe des mesures
plt.plot(x_mod, y_mod, 'r-', label = "modèle") # Courbe du modèle
plt.xlabel("x")                                # Etiquette en abscisse
plt.ylabel("y")                                # Etiquette en ordonnée
plt.title("Trajectoire et modèle associé")     # Titre
plt.legend()                                   # Affichage de la légende
plt.grid()                                     # Affichage de la grille
plt.show()                                     # Affichage de la figure
alternate text
>>> %Run
a =  -1.028916427645116
b =  -0.47659343034449314
c =  0.4413962087861902

Cas d'une fonction affine (régression linéaire)

Il est également possible de procéder à une régression linéaire avec la fonction polyfit().

import numpy as np
import matplotlib.pyplot as plt

x = np.array([0, 1.01, 2.02, 2.99, 3.98])    # Données en abscisse
y = np.array([10.32, 7.91, 6.13, 3.97, 2.25]) # Données en ordonnée

# MODELISATION
a, b = np.polyfit(x, y, 1)  # Modélisation par un polynome de degré 2
print("a = ", a)
print("b = ", b)

x_mod = np.linspace(0, 4, 50)         # Nouvelle abscisse pour la modélisation
y_mod = a*x_mod + b                   # Ordonnées de la fonction affine

plt.plot(x_mod, y_mod, '-')           # Tracé de la droite
plt.plot(x, y, '+')                   # Tracé des points et de la fonction affine
plt.title('Régression linéaire')      # Titre
plt.xlabel('x')                       # Etiquette en abscisse
plt.xlim(-1,5)                        # Echelle en abscisse
plt.ylabel('y')                       # Etiquette en ordonnée
plt.ylim(0, 12)                       # Echelle en ordonnée
plt.grid()
plt.show()                            # Affichage
alternate text
>>> %Run
a =  -2.020342070640623
b =  10.156684141281245

Modélisation à partir d'une fonction quelconque

La fonction curve_fit() du module scipy.optimize réalise un ajustement par rapport à une fonction mathématique quelconque et retourne les paramètres de cette fonction.

Considérons la charge d'un condensateur (initialement déchargé) de capacité \(C\) à travers un résistance \(R\) sous tension constante \(E\) :

\[u(t) = E \times (1-e^{-{t}/{\tau}}) \hspace{1cm}\text{avec}\hspace{1cm} \tau = RC\]
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

# DONNEES
t = [0, 1, 2, 3, 4, 5, 6, 7, 8, 9]
u = [0., 3.935, 6.321, 7.769, 8.647, 9.179, 9.502, 9.698, 9.817, 9.889]

# DEFINITION DE LA FONCTION QUELCONQUE
def fct(t, E, tau):
   return E*(1-np.exp(-t/tau))         # Expression du modèle

# MODELISATION
(E, tau), pcov = curve_fit(fct, t, u)   # Détermination des paramètres du modèle
print("E = ", E)                        # Affichage de A
print("tau =", tau)                     # Affichage de tau

# CONSTRUCTION DU MODELE
t_mod = np.linspace(0, 10, 50)
u_mod = fct(t_mod, E, tau)

# COURBES
plt.plot(t_mod, u_mod, label="Modèle")  # Courbe du modèle
plt.plot(t, u, "r+", label="Tension")   # Nuage des points
plt.title("Charge d'un condensateur")   # Titre
plt.xlabel("t (s)")                     # Etiquette en abscisse
plt.ylabel("u (V)")                     # Etiquette en ordonnée
plt.legend()                            # Affichage de la légende
plt.grid()                              # Affichage de la grille
plt.show()                              # Affichage de la figure
alternate text
>>> %Run
E =  9.99999510282223
tau = 1.9999259182304618

Interpolation

La fonction interp1d() du module scipy.interpolate retourne une fonction Python pouvant estimer la position d'un point par interpolation.

import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import interp1d

# DONNEES
x = [0, 1, 2, 3, 4, 5, 6, 7, 8, 9]
y = [1., 0.720, 0.511, 0.371, 0.260, 0.190, 0.133, 0.099, 0.069, 0.048]

# INTERPOLATION
f = interp1d(x, y, kind='cubic')  # ou kind = 'linear'

# CONSTRUCTION DE LA FONCTION
x_mod = np.linspace(0, 9, 50)    # Nouvelle abscisse pour f(x)
y_mod = f(x_mod)                 # Calcul des ordonnées de f(x)

# COURBES
plt.plot(x_mod, y_mod, label="Interpolation")  # Tracé de la fonction
plt.plot(x, y, 'r+', label="Points")           # Tracé des points
plt.xlabel("x")                                # Etiquette en abscisse
plt.ylabel("y")                                # Etiquette en ordonnée
plt.legend()                                   # Affichage de la légende
plt.grid()                                     # Affichage de la grille
plt.show()                                     # Affichage de la figure
alternate text
>>> %Run
>>> f(1)
array(0.72)
>>> f(10)
ValueError: A value in x_new is above the interpolation range.
  • La fonction f() est capable d'estimer par interpolation la position de n'importe quel point.

  • L'interpolation n'est naturellement pas possible en dehors de la plage des valeurs de x.

>>> f.x
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
>>> f.y
array([1., 0.72, 0.511, 0.371, 0.26 , 0.19 , 0.133, 0.099, 0.069, 0.048])
  • La fonction f() est bien définie par rapport au nuage des points passé en argument à la fonction interp1d().