## RÉGRESSION LINÉAIRE
## Armel MARTIN - Septembre 2021
## distribué sous licence CC-by-nc-sa
import numpy as np
import matplotlib.pyplot as plt

## Quartet d'Anscombe - N°1
x = np.array([10., 8., 13., 9., 11., 14., 6, 4, 12, 7, 5]) 
y = np.array([8.04, 6.95, 7.58, 8.81, 8.33, 9.96, 7.24, 4.26, 10.84, 4.82, 5.68]) 

# Regression linéaire y = ax+b
a,b = np.polyfit(x,y,1)

# ou calcul direct: a = cov(x,y) / var(x) et b = ym - a*xm
xm, ym = np.mean(x), np.mean(y)
var_x = np.mean((x - xm)**2) # ou np.var(x), estimateur biaisé par défaut
cov_xy = np.mean((x-xm)*(y-ym))
a2 = cov_xy / var_x
b2 = np.mean(y) - a2 * np.mean(x)

# Coef de corrélation: r = Cov(x,y) / sqrt(var(x)*var(y))
var_y = np.mean((y - ym)**2)
r = cov_xy/np.sqrt(var_x*var_y)

# Valeurs calculées par le modèle
y_reg = b + a*x

# Analyse des composantes de variance
SST = np.mean((y - np.mean(y))**2) # variance totale
SSE = np.mean((y_reg - np.mean(y))**2) # variance expliquée (modèle)
#SSR = np.mean((y - y_reg)**2) # variance résiduelle (-> 'res' donné par polyfit)

# Coeff de détermination (R²)
R2 = SSE / SST
# On vérifie que: r**2 = R2 (valable uniquement pour un modèle affine)

# Affichage des résultats en console
print('b=',b,' ; a=',a,' ; R2=',R2,' ; r=',r)

## Fabriquer l'équation de la courbe (type string)
equation = 'y = '+format(a,'.3e')+' x + '+format(b,'.3e')
#equation = 'y = {0:.3e} x + {1:.3e}'.format(a,b)

plt.figure(figsize=[6,3])
plt.plot(x,y,'o',color='orange')
plt.plot(x,y_reg,'-k')
plt.xlabel('x')
plt.ylabel('y')
plt.title("Quartet d'Anscombe - N°1")
plt.annotate(equation,(4,10))
plt.grid()
plt.savefig('ma_figure.png')
#plt.tight_layout()
plt.show()
#plt.clf()