Statistiques et Probabilités
[!bug] Avant-propos Ce document est découpé en deux parties :
- Un cours de statistiques appliquées à la data science et à l’intelligence artificielle. Chaque notion est expliquée en détail avec définitions rigoureuses, formules complètes, valeurs théoriques, et exemples concrets tirés des domaines de la chimie organique, cinétique, photochimie, photocatalyse et photophysique.
- Un cours exhaustif sur l’analyse de régression, outil fondamental de la data science et de l’intelligence artificielle. Les régressions permettent de modéliser les relations entre variables et constituent la base de l’inférence causale. Chaque concept est développé avec rigueur mathématique et illustré par des exemples issus de la chimie quantitative, de la cinétique et de la spectroscopie.
Statistiques Descriptives
Les statistiques descriptives constituent le premier niveau d’analyse de données expérimentales. Elles permettent de résumer, visualiser et comprendre les caractéristiques essentielles d’un ensemble de mesures.
Typologie des données
[!NOTE] Définition Données catégorielles :
- Nominales : Pas d’ordre (types de réaction, solvants)
- Ordinales : Ordre naturel (faible/moyen/fort, grades)
Données numériques :
- Discrètes : Valeurs entières (nombre de photons, molécules)
- Continues : Toute valeur réelle (température, concentration)
Proportions : Agrégation de données nominales binaires en valeur numérique [0,1]
[!example]- Classification de mesures spectroscopiques Dans une étude photophysique :
- Nominal : Type de transition (π→π, n→π, CT)
- Ordinal : Intensité de fluorescence (faible < moyenne < forte)
- Discret : Nombre de bandes d’absorption
- Continu : Longueur d’onde du maximum (λmax = 452.7 nm)
- Proportion : Fraction de molécules fluorescentes (65%)
Variables statistiques
Moyenne Arithmétique
[!NOTE] Définition La moyenne arithmétique (ou espérance empirique) est la somme de toutes les valeurs divisée par le nombre total d’observations. C’est la mesure de tendance centrale la plus couramment utilisée.
\[\bar{x} = \frac{1}{n} \sum_{i=1}^{n} x_i = \frac{x_1 + x_2 + ... + x_n}{n}\]$\mu$ pour une population.
où :
- $\bar{x}$ : moyenne de l’échantillon
- $n$ : nombre d’observations
- $x_i$ : i-ème observation
Propriétés :
- Sensible aux valeurs extrêmes (outliers)
- Linéaire : $\overline{aX + b} = a\bar{x} + b$
- Somme des écarts à la moyenne = 0 : $\sum_{i=1}^{n} (x_i - \bar{x}) = 0$
- Minimise la somme des carrés des écarts
[!example]- Cinétique chimique Lors d’une étude de la dégradation photocatalytique du bleu de méthylène sous irradiation UV, on mesure la constante de vitesse apparente $k_{app}$ (en min⁻¹) sur 8 expériences répétées :
Valeurs : 0.045, 0.052, 0.048, 0.051, 0.047, 0.049, 0.053, 0.050 min⁻¹
\[\bar{k}_{app} = \frac{0.045 + 0.052 + 0.048 + 0.051 + 0.047 + 0.049 + 0.053 + 0.050}{8} = \frac{0.395}{8} = 0.049375 \text{ min}^{-1}\]Cette moyenne permet d’estimer la constante de vitesse “vraie” de la réaction.
[!example]- Spectroscopie d’absorption Pour caractériser un nouveau chromophore organique, on mesure son coefficient d’extinction molaire $\varepsilon$ à 450 nm sur 5 préparations indépendantes :
Valeurs : 12 500, 12 800, 12 300, 12 700, 12 600 L·mol⁻¹·cm⁻¹
\[\bar{\varepsilon} = \frac{12500 + 12800 + 12300 + 12700 + 12600}{5} = \frac{63900}{5} = 12780 \text{ L·mol}^{-1}\text{·cm}^{-1}\]
Médiane
[!NOTE] Définition La médiane est la valeur centrale qui sépare l’ensemble ordonné en deux parties égales : 50% des valeurs sont inférieures et 50% sont supérieures. Elle est robuste aux valeurs extrêmes.
Calcul :
- Si $n$ est impair : médiane = valeur en position $\frac{n+1}{2}$
- Si $n$ est pair : médiane = $\frac{x_{n/2} + x_{(n/2)+1}}{2}$
Propriétés :
- Robuste aux outliers
- Minimise la somme des valeurs absolues des écarts
- Médiane = Q₂ (deuxième quartile)
- Pour une distribution symétrique : médiane ≈ moyenne
[!example]- Rendement quantique de fluorescence On mesure le rendement quantique $\Phi_F$ de 9 dérivés de fluorescéine :
Valeurs ordonnées : 0.12, 0.18, 0.22, 0.25, 0.31, 0.35, 0.38, 0.42, 0.88
Position centrale : $\frac{9+1}{2} = 5$
Médiane = 0.31
La valeur 0.88 est un outlier (probablement un dérivé avec structure particulière), mais la médiane reste représentative du groupe principal.
Mode
[!NOTE] Définition Le mode est la valeur qui apparaît avec la plus grande fréquence dans l’ensemble de données. Un ensemble peut être :
- Unimodal : un seul pic de fréquence
- Bimodal : deux pics de fréquence
- Multimodal : plusieurs pics de fréquence
- Amodal : aucune valeur ne se répète
[!example]- États de spin en photochimie Dans une étude de la durée de vie de l’état triplet T₁ (en microsecondes) pour 20 molécules photoactives, on observe une distribution bimodale :
- Groupe 1 (états triplets courts) : τ ≈ 5 μs (10 molécules)
- Groupe 2 (états triplets longs) : τ ≈ 120 μs (10 molécules)
Cette bimodalité révèle deux mécanismes distincts de désactivation : ISC rapide pour le groupe 1, phosphorescence lente pour le groupe 2.
Variance et Écart-Type
[!NOTE] Définition La variance mesure la dispersion des données autour de la moyenne. L’écart-type en est la racine carrée et s’exprime dans la même unité que les données.
Variance de population : \(\sigma^2 = \frac{1}{N} \sum_{i=1}^{N} (x_i - \mu)^2\)
Variance d’échantillon (estimateur non biaisé) : \(s^2 = \frac{1}{n-1} \sum_{i=1}^{n} (x_i - \bar{x})^2\)
Écart-type : \(\sigma = \sqrt{\sigma^2} \quad \text{ou} \quad s = \sqrt{s^2}\)
Pourquoi n-1 ? Le diviseur $(n-1)$ est appelé correction de Bessel. Il compense le biais introduit par l’utilisation de $\bar{x}$ au lieu de $\mu$. C’est la notion de [[Statistiques et probabilités#7. Degrés de Liberté|degrés de libertés]] : une fois la moyenne calculée, seules $(n-1)$ valeurs sont libres.
Pourquoi élever au carré ?
- Direction nécessaire : Les écarts simples s’annulent (certains positifs, certains négatifs)
- Pondération des outliers : Les valeurs extrêmes sont amplifiées, ce qui reflète leur impact sur la variabilité
- Propriétés mathématiques : Facilite les calculs de propagation d’erreur et les tests statistiques
[!example]- Énergie de liaison détaillée On mesure l’énergie de dissociation de liaison C-H dans 6 alcanes différents par photolyse laser (en kJ/mol) :
Valeurs : 410, 415, 408, 412, 411, 414 kJ/mol
Étape 1 : Calcul de la moyenne \(\bar{x} = \frac{410 + 415 + 408 + 412 + 411 + 414}{6} = \frac{2470}{6} = 411.67 \text{ kJ/mol}\)
Étape 2 : Calcul des écarts au carré
- $(410 - 411.67)^2 = (-1.67)^2 = 2.79$
- $(415 - 411.67)^2 = (3.33)^2 = 11.09$
- $(408 - 411.67)^2 = (-3.67)^2 = 13.47$
- $(412 - 411.67)^2 = (0.33)^2 = 0.11$
- $(411 - 411.67)^2 = (-0.67)^2 = 0.45$
- $(414 - 411.67)^2 = (2.33)^2 = 5.43$
Étape 3 : Variance \(s^2 = \frac{2.79 + 11.09 + 13.47 + 0.11 + 0.45 + 5.43}{6-1} = \frac{33.34}{5} = 6.67 \text{ (kJ/mol)}^2\)
Étape 4 : Écart-type \(s = \sqrt{6.67} = 2.58 \text{ kJ/mol}\)
Interprétation : La dispersion de ±2.58 kJ/mol représente la variabilité naturelle entre différentes liaisons C-H et l’incertitude expérimentale.
Coefficient de Variation
[!NOTE] Définition Le coefficient de variation exprime l’écart-type en pourcentage de la moyenne. Il permet de comparer la variabilité relative de mesures d’unités différentes.
\[CV = \frac{s}{\bar{x}} \times 100\% = \frac{\sigma}{\mu} \times 100\%\]Interprétation :
- CV < 10% : faible variabilité (reproductibilité excellente)
- 10% < CV < 20% : variabilité modérée (reproductibilité acceptable)
- CV > 20% : forte variabilité (nécessite investigation)
[!example]- Comparaison expérimentale Expérience A - Mesure de concentration (μM) :
- $\bar{x}_A = 150$ μM, $s_A = 15$ μM
- $CV_A = \frac{15}{150} \times 100\% = 10\%$
Expérience B - Mesure d’absorbance (sans unité) :
- $\bar{x}_B = 0.85$, $s_B = 0.085$
- $CV_B = \frac{0.085}{0.85} \times 100\% = 10\%$
Bien que les unités soient différentes, les deux expériences ont la même variabilité relative.
| Concept | Notation | Formule clé | Interprétation |
|---|---|---|---|
| Moyenne | $\bar{x}$ | $\frac{\sum x_i}{n}$ | Tendance centrale |
| Variance | $s^2$ | $\frac{\sum (x_i - \bar{x})^2}{n-1}$ | Dispersion |
| Écart-type | $s$ | $\sqrt{s^2}$ | Dispersion (même unité que X) |
| Covariance | $\text{Cov}(X,Y)$ | $\frac{\sum (x_i - \bar{x})(y_i - \bar{y})}{n-1}$ | Relation directionnelle |
| Corrélation | $r_{XY}$ | $\frac{\text{Cov}(X,Y)}{s_X s_Y}$ | Relation normalisée |
Covariance
[!NOTE] Définition La covariance mesure la tendance de deux variables à varier ensemble. Une covariance positive indique que les variables augmentent ensemble, une covariance négative indique qu’elles varient en sens opposé.
Formule pour échantillon : \(\text{Cov}(X,Y) = \frac{1}{n-1} \sum_{i=1}^{n} (x_i - \bar{x})(y_i - \bar{y})\)
Formule pour population : \(\text{Cov}(X,Y) = \frac{1}{N} \sum_{i=1}^{N} (x_i - \mu_X)(y_i - \mu_Y)\)
Propriétés :
- $\text{Cov}(X,X) = \text{Var}(X)$
- $\text{Cov}(X,Y) = \text{Cov}(Y,X)$ (symétrie)
- $\text{Cov}(aX, bY) = ab \cdot \text{Cov}(X,Y)$
[!example]- Photocatalyse TiO₂ On étudie la relation entre la concentration en catalyseur TiO₂ (mg/L) et le taux de dégradation d’un polluant (% en 30 min) :
TiO₂ (X) Dégradation (Y) $(x_i - \bar{x})$ $(y_i - \bar{y})$ Produit 50 35 -100 -28.8 2880 100 52 -50 -11.8 590 150 68 0 4.2 0 200 79 50 15.2 760 250 85 100 21.2 2120 Moyennes :
- $\bar{x} = 150$ mg/L
- $\bar{y} = 63.8$ %
Covariance : \(\text{Cov}(X,Y) = \frac{2880 + 590 + 0 + 760 + 2120}{5-1} = \frac{6350}{4} = 1587.5\)
La covariance positive confirme que plus on augmente la concentration en TiO₂, plus la dégradation est efficace.
Corrélation
[!NOTE] Définition Le coefficient de corrélation de Pearson normalise la covariance par les écarts-types, donnant une mesure sans dimension comprise entre -1 et +1.
\[r = \frac{\text{Cov}(X,Y)}{s_X \cdot s_Y} = \frac{\sum_{i=1}^{n} (x_i - \bar{x})(y_i - \bar{y})}{\sqrt{\sum_{i=1}^{n} (x_i - \bar{x})^2} \cdot \sqrt{\sum_{i=1}^{n} (y_i - \bar{y})^2}}\]Interprétation :
$ r < 0.3$ : corrélation faible
$0.3 ≤ r < 0.7$ : corrélation modérée
$ r ≥ 0.7$ : corrélation forte - $r = 1$ : corrélation positive parfaite
- $r = -1$ : corrélation négative parfaite
==Attention : Corrélation n’implique pas causalité !==
![[Pasted_image_20251027095105.png]] ![[Pasted_image_20251027095119.png]]
Correlation vs. Independence?
- (X,Y) independent ⇒ Corr(X,Y) = 0
- Corr(X,Y) = 0 ⇏ (X,Y) independent (as
ronly captures linear dependence - cf Wikipedia)
Summary statistics and correlation are important, but… ⚠️ Beware of the Datasaurus 🦕 Summary statistics and correlation are not enough, you need to conduct Exploratory data analysis too!
![[datasaurus.gif]]
[!example]- Photophysique : loi de Stokes-Einstein On étudie la relation entre la viscosité du solvant (cP) et le temps de rotation d’une sonde fluorescente (ns) :
Viscosité (X) Temps rotation (Y) 0.55 0.12 0.89 0.18 1.31 0.28 2.24 0.45 8.95 1.82 Après calculs :
\[r = \frac{2.28}{3.45 \times 0.67} = 0.987\]
- $s_X = 3.45$ cP
- $s_Y = 0.67$ ns
- $\text{Cov}(X,Y) = 2.28$
La corrélation très forte (r = 0.987) confirme la validité de la relation de Stokes-Einstein pour cette sonde.
Skewness
Mesure la symétrie de la courbe, indicateur calculable facilement avec python.
![[Pasted_image_20251027092853.png]]
| Permet de mesure un écart par exemple à la [[Statistiques et probabilités#Loi Normale (Gaussienne) | loi normale]]. |
Probabilités
Expériences aléatoires
- A random experiment is a process by which we observe something uncertain
- After the experiment, the result of the random experiment is known: it is the (outcome)
- The set of all possible outcomes is called the sample space S,Ω or U (Univers 🇫🇷)
Examples:
- Toss a coin. S={H,T}
- Roll a die. S={1,2,3,4,5,6}
- Observe the number of goals in a soccer match. S={0,1,2,3,…—
When we repeat a random experiment several times, we call each one of them a trial (épreuve 🇫🇷).
Sample space S is defined based on how you define the random experiment…
❓ If the experiment is Toss a coin three times, what’s the sample space?
![[Pasted_image_20251027095943.png]]
The goal is to assign a probability to events, defined as subsets of a sample space S.
Un évènement comme E = {(H, H, H), (H, H, T)} correspond à “j’ai obtenu deux fois faces au premier lancer”, de probabilité 2/8 (car deux sous ensemble sur 8 totaux).
Variables aléatoires (PMF, E, Bernouilli process)
To analyze random experiments, we focus on some numerical aspects of the experiment.
For example, in a soccer game we may be interested in the number of goals, shots, etc…
[!NOTE] Définition A random variable X is a function from the sample space to the real numbers: X : S → R
The range of a random variable X, shown by Range(X), is the set of possible values of X.
In our example (X as the number of heads): Range(X)={0,1,2} ![[Pasted_image_20251027102954.png]]
[!note] Probability Mass Function The PMF for X is defined as $Vxi in Range(X), pmfX(xi)=P(X=xi)$ ![[Pasted_image_20251027102937.png]]
[!NOTE] Expected value Intuitively, a random variable’s expected value E[X]
represents the average of a large number of independent realizations of the random variable X
.
X
being a discrete random value:
E[X]=n∑i=1xip(xi).
E[“Number of heads in 2 tosses”]=0∗0.25+1∗0.5+2∗0.25=1
[!NOTE] Bernoulli process Take a random experiment with exactly two possible outcomes and repeat this experiment multiple times
✏️ e.g let’s define our Bernoulli experiment as: “toss a coin once with a proba
p = 0.3of getting (e.g. Heads) at each toss”and repeat this experiment a
size = 10000times (large to smooth-out noise)size = 10000 p = 0.3 np.random.binomial(n=1, p=p, size=size)
array([0, 0, 1, ..., 1, 0, 0])# Plot results after 10000 repetitions sns.histplot(np.random.binomial(n=1, p=p, size=size), kde=False);![[Pasted_image_20251027103601.png]]
| Voir les autres distributions dans [[Statistiques et probabilités#Distributions de Probabilité | la partie associée]]. |
Types de Probabilité
[!NOTE] Définition Probabilité conjointe : Probabilité que deux événements se produisent simultanément \(P(A \cap B) = P(A \text{ et } B)\)
Probabilité marginale : Probabilité d’un événement sans considérer les autres variables \(P(A) = \sum_i P(A \cap B_i)\)
Probabilité conditionnelle : Probabilité d’un événement sachant qu’un autre s’est produit \(P(A|B) = \frac{P(A \cap B)}{P(B)}\)
Union de probabilités : \(P(A \cup B) = P(A) + P(B) - P(A \cap B)\)
[!example]- Sélectivité de photoréaction Une photoréaction peut suivre deux voies :
- Voie A (isomérisation) : P(A) = 0.6
- Voie B (cyclisation) : P(B) = 0.3
- Les deux voies simultanément : P(A ∩ B) = 0.1
Probabilité qu’au moins une réaction se produise : \(P(A \cup B) = 0.6 + 0.3 - 0.1 = 0.8\)
Probabilité de cyclisation sachant qu’il y a eu isomérisation : \(P(B|A) = \frac{0.1}{0.6} = 0.167\)
Indépendance Statistique
[!NOTE] Définition Deux événements A et B sont indépendants si et seulement si :
\[P(A \cap B) = P(A) \times P(B)\]Ou de manière équivalente : \(P(A|B) = P(A)\)
L’occurrence de B ne modifie pas la probabilité de A.
[!example]- Excitation indépendante de chromophores Dans un système bi-chromophore, la probabilité d’excitation de chaque chromophore est :
- Chromophore 1 : P(C₁) = 0.3
- Chromophore 2 : P(C₂) = 0.4
Si les excitations sont indépendantes :
- P(les deux excités) = 0.3 × 0.4 = 0.12
- P(au moins un excité) = P(C₁) + P(C₂) - P(C₁ ∩ C₂) = 0.3 + 0.4 - 0.12 = 0.58
Visualisation
Quartiles et Boîtes à Moustaches
[!NOTE] Définition Les quartiles divisent un ensemble ordonné en quatre parties égales :
- Q₁ (premier quartile) : 25% des données sont inférieures
- Q₂ (deuxième quartile) : médiane, 50% des données sont inférieures
- Q₃ (troisième quartile) : 75% des données sont inférieures
Intervalle Interquartile (IQR) : \(IQR = Q_3 - Q_1\)
![[Pasted_image_20251027093745.png]]
Détection des valeurs aberrantes :
- Outlier inférieur : $x < Q_1 - 1.5 \times IQR$
- Outlier supérieur : $x > Q_3 + 1.5 \times IQR$
- Outlier extrême : utiliser 3.0 × IQR au lieu de 1.5
[!example]- Temps de vie de fluorescence On mesure le temps de vie de fluorescence (en ns) de 15 nanoparticules de CdSe/ZnS :
Données ordonnées : 8.2, 9.1, 9.5, 10.2, 10.8, 11.2, 11.5, 12.1, 12.4, 12.9, 13.3, 13.8, 14.2, 14.7, 22.5 ns
Calcul des quartiles :
- Position Q₁ : $\frac{15+1}{4} = 4$ → Q₁ = 10.2 ns
- Position Q₂ : $\frac{15+1}{2} = 8$ → Q₂ = 12.1 ns (médiane)
- Position Q₃ : $3 \times \frac{15+1}{4} = 12$ → Q₃ = 13.8 ns
IQR : \(IQR = 13.8 - 10.2 = 3.6 \text{ ns}\)
Détection des outliers :
- Limite supérieure : $13.8 + 1.5 \times 3.6 = 19.2$ ns
- La valeur 22.5 ns est un outlier (probablement une nanoparticule avec défauts de surface)
Histogrammes
A histogram is a representation of the distribution of numerical data.
It is an estimate of the probability distribution of a continuous variable.
⚠️ Histogram (continuous variable) ≠ Bar chart (categorical or discrete variable)
ax = sns.histplot(weights_df["weight"], bins=weights_df["weight"].max() - weights_df["weight"].min())
ax.set_xlabel("Weight (in lbs)")
plt.show()
![[Pasted_image_20251027091747.png]]
Histogram Bins
Instead of drawing one bar per integer in [95,215], we can create 12 bins and count weights falling into these intervals.
f, (ax1, ax2, ax3) = plt.subplots(1, 3, figsize=(15, 4))
ax1.hist(weights_df["weight"], bins=[95, 105, 115, 125, 135, 145, 155, 165, 175, 185, 195, 205, 215])
ax2.hist(weights_df["weight"], bins=12)
sns.histplot(weights_df["weight"], bins=12, ax=ax3, kde=True)
plt.show()
![[Pasted_image_20251027091836.png]]
Cumulative plots
Alternatively, we can plot the count of weights inferior to a certain value.
Instead of the counts, we can also plot the density (sums up to 1) or percentage (sums up to 100).
f, (ax1, ax2) = plt.subplots(1, 2, figsize=(15, 4))
ax1.hist(weights_df["weight"], bins=25, cumulative=True, density=True)
sns.histplot(weights_df["weight"], binwidth=5, ax=ax2, kde=True,
cumulative=True, stat="percent");
![[Pasted_image_20251027091947.png]]
La moyenne est sensible aux outliers, pas la médiane !! Par exemple pour un salaire en entretien d’embauche le salaire médian est plus révélateur que le salaire moyen.
Distributions
Loi Uniforme Continue
[!NOTE] Définition Toutes les valeurs dans l’intervalle [a, b] ont la même densité de probabilité.
Fonction de densité : \(f(x) = \begin{cases} \frac{1}{b-a} & \text{si } a \leq x \leq b \\ 0 & \text{sinon} \end{cases}\)
Espérance et variance : \(E(X) = \frac{a+b}{2} \quad,\quad \text{Var}(X) = \frac{(b-a)^2}{12}\)
[!example]- Orientation moléculaire aléatoire L’angle d’orientation θ d’une molécule adsorbée sur une surface suit une loi uniforme sur [0, 2π].
- E(θ) = π radians
- Var(θ) = π²/3 ≈ 3.29 rad²
- P(0 < θ < π/2) = (π/2)/(2π) = 1/4
Donc 25% des molécules sont orientées dans le premier quadrant.
Loi Binomiale
[!NOTE] Définition Modélise le nombre de succès dans n essais indépendants avec probabilité p constante.
Fonction de masse : \(P(X = k) = \binom{n}{k} p^k (1-p)^{n-k} = \frac{n!}{k!(n-k)!} p^k (1-p)^{n-k}\)
Paramètres :
- E(X) = np
- Var(X) = np(1-p)
- Mode = ⌊(n+1)p⌋
Conditions :
- Essais indépendants
- Probabilité constante
- Nombre fixe d’essais
[!example]- Excitation multiphotonique Dans une expérience d’absorption à deux photons, chaque molécule a p = 0.1 de probabilité d’être excitée. Sur n = 20 molécules, la probabilité qu’exactement 3 soient excitées :
\(P(X=3) = \binom{20}{3} (0.1)^3 (0.9)^{17}\) \(= \frac{20!}{3!17!} \times 0.001 \times 0.1668\) \(= 1140 \times 0.001 \times 0.1668 = 0.190\)
Environ 19% des expériences donneraient exactement trois excitations.
Vocabulary:
- An experiment with only two outcomes (success/failure) has a Bernoulli (p) distribution.
- When we repeat this Bernoulli process n times and count the numbers success $X$, we talk about a Binomial (n,p) distribution.
🧑🏻🏫 Counting the number of “successes” among n repeated experiments: Counting the number of ways to get $k$ heads (the successes) among $n$ flips (the repeated experiment) … …is equivalent to counting the number of ways to select $k$ items from a set that has $n$ distinct elements, such that the order of selection does not matter
-
If the order mattered, picking $k$ elements one-by-one among $n$ could be done in $n(n−1)…(n−k+1)$ ways ($n$ choice for the first element, $n-1$ for the second, …, $n−k+1$ for the $k$-th )
-
However, in this ordered count, any unordered set of $k$ elements have been counted $k(k-1)(k-2)…$ times ($k$ choice for the first, $k-1$ for the second, etc…)
-
Therefore, if we want the unordered count, we have to compensate for (divide by) them. Hence, the number of ways to get k successes out of n experiments is given by: \({\frac {n(n-1)\dotsb (n-k+1)}{k(k-1)\dotsb 1}}\) This is mathematically equivalent to: \(\frac{n!}{k! (n - k)!} \text{ , where } n! = 1\times 2 \times ... \times n\) and is written
- $\binom{n}{k}$ reads as
"n choose k", or"binomial coefficient for k among n" - $n!$ reads as
"n factorial"
- When the order doesn’t matter, it is a Combination.
- When the order does matter it is a Permutation.
Permutation and Combination
In combinatorics, the key idea is to determine how many ways we can choose or arrange items — and everything depends on whether order matters and whether repetition is allowed.
Permutation — when order matters
A permutation is an arrangement of items where order matters. We can think of it as the number of possible “ordered lists” we can make from a set.
There are two main cases:
-
With repetition: Each item can be chosen more than once. Example: digits in a PIN code (0000 to 9999). If there are n choices and r positions: \(P = n^r\)
-
Without repetition: Once an item is chosen, it can’t be chosen again. So our choices get reduced each time. Example: first place, second place, third place in a race. \(P(n, r) = n \times (n - 1) \times (n - 2) \times \dotsb = \frac{n!}{(n - r)!}\)
[!EXAMPLE] Permutation example There are 10 people, and we want to choose a president, a vice-president, and a secretary: \(P(10,3) = 10 \times 9 \times 8 = 720\) Because changing the order changes the result — order matters.
Combination — when order does not matter
A combination is a selection of items where order doesn’t matter — we only care about which items are chosen, not in what sequence.
Again, there are two main cases:
-
With repetition (a.k.a. combinations with replacement): Items can be chosen more than once. Example: choosing 3 scoops of ice cream out of 5 flavors — we can repeat flavors. Formula: \(C(n + r - 1, r) = \frac{(n + r - 1)!}{r!(n - 1)!}\)
-
Without repetition: Once an item is chosen, it cannot be chosen again. Example: choosing 3 students from 10 for a team. Formula: \(C(n, r) = \binom{n}{r} = \frac{n!}{r!(n - r)!}\)
[!EXAMPLE] Combination example Choosing 3 students from a group of 10: \(\binom{10}{3} = 120\) The team {Alice, Bob, Carol} is the same as {Carol, Bob, Alice} — order does not matter.
Summary Table
| Type | Order matters? | Repetition allowed? | Formula | Example |
|---|---|---|---|---|
| Permutation (with repetition) | ✅ Yes | ✅ Yes | $n^r$ | 4-digit PIN (0000–9999) |
| Permutation (without repetition) | ✅ Yes | ❌ No | $\dfrac{n!}{(n - r)!}$ | Race ranking (gold, silver, bronze) |
| Combination (with repetition) | ❌ No | ✅ Yes | $\dfrac{(n + r - 1)!}{r!(n - 1)!}$ | Ice cream flavors (3 scoops, 5 choices) |
| Combination (without repetition) | ❌ No | ❌ No | $\dfrac{n!}{r!(n - r)!}$ | Team selection (3 from 10) |
[!success] Programmation
n = 5 np.random.binomial(n=n, p=p, size=size) array([1, 0, 4, ..., 2, 2, 0]) # Plot results after counting the number of heads ("5 trials") n = 5 sns.histplot(np.random.binomial(n=n, p=p, size=size), kde=False);![[Pasted_image_20251027103823.png]]
| The sum of n Bernoulli Distributions B(p) is called a Binomial Distribution B(n,p). As n increases, Binomial Distribution B(n,p) approximates a [[Statistiques et probabilités#Loi Normale (Gaussienne) | normal distribution]] N centered on n∗p. |
Loi de Poisson
[!NOTE] Définition Modélise le nombre d’occurrences d’un événement rare pendant un intervalle fixe.
Fonction de masse : \(P(X = k) = \frac{\lambda^k e^{-\lambda}}{k!}\)
Paramètres :
- E(X) = Var(X) = λ
- Approximation de la binomiale quand n → ∞, p → 0, np = λ constant
Conditions :
- Taux constant d’occurrence
- Événements indépendants
- Probabilité faible sur petit intervalle
[!example]- Comptage de photons uniques Une photodiode détecte en moyenne λ = 4 photons par microseconde.
Probabilité de détecter exactement 2 photons : \(P(X=2) = \frac{4^2 e^{-4}}{2!} = \frac{16 \times 0.0183}{2} = 0.146\)
Probabilité de détecter au moins 1 photon : \(P(X \geq 1) = 1 - P(X=0) = 1 - e^{-4} = 1 - 0.0183 = 0.982\)
Le détecteur a 98.2% de chances de capter au moins un photon par μs.
Loi Exponentielle
[!NOTE] Définition Modélise le temps entre événements dans un processus de Poisson.
Fonction de densité : \(f(t) = \lambda e^{-\lambda t} \quad \text{pour } t \geq 0\)
Paramètres :
- E(T) = 1/λ
- Var(T) = 1/λ²
Propriété sans mémoire : P(T > s+t T > s) = P(T > t)
[!example]- Décroissance de fluorescence La durée de vie de fluorescence suit une loi exponentielle avec τ = 1/λ = 5 ns.
Probabilité qu’une molécule fluoresce avant 3 ns : \(P(T < 3) = 1 - e^{-3/5} = 1 - e^{-0.6} = 0.451\)
Temps médian de fluorescence : \(P(T < t_{1/2}) = 0.5 \Rightarrow t_{1/2} = -\tau \ln(0.5) = 5 \times 0.693 = 3.47 \text{ ns}\)
Loi Normale (Gaussienne)
[!NOTE] Définition La loi normale est centrale en statistiques grâce au théorème central limite.
Fonction de densité : \(f(x) = \frac{1}{\sigma \sqrt{2\pi}} e^{-\frac{(x-\mu)^2}{2\sigma^2}}\)
Propriétés :
- Symétrique par rapport à μ
- Points d’inflexion en μ ± σ
- 68.3% des valeurs dans [μ-σ, μ+σ]
- 95.4% dans [μ-2σ, μ+2σ]
- 99.7% dans [μ-3σ, μ+3σ]
![[Pasted_image_20251027104125.png 500]] Importance :
- Théorème central limite
- Nombreux phénomènes naturels
- Base de nombreux tests statistiques
[!example]- Largeur spectrale d’émission La longueur d’onde d’émission d’un laser suit N(μ=532 nm, σ=0.5 nm).
Probabilité d’émission entre 531 et 533 nm : \(Z_1 = \frac{531-532}{0.5} = -2, \quad Z_2 = \frac{533-532}{0.5} = 2\) \(P(531 < X < 533) = \Phi(2) - \Phi(-2) = 0.977 - 0.023 = 0.954\)
Intervalle contenant 99% des émissions : \(\mu \pm 2.58\sigma = 532 \pm 1.29 = [530.71, 533.29] \text{ nm}\)
[!NOTE] Programmation - Probability Density Function
def plot_normal_distribution(mu, variance): sigma = math.sqrt(variance) x = np.linspace(-10, 10, 100) plt.plot(x, stats.norm.pdf(x, mu, sigma), label=f"μ={mu}, σ²={variance}") plot_normal_distribution(0, 1) plot_normal_distribution(1, 2) plot_normal_distribution(-3, 5) plt.legend() plt.show()![[Pasted_image_20251027104105.png]]
Probability Density Function or Cumulative Distribution Function ?
PDF is the probability distribution of a continuous random variable. PDFs tells us the probability of an infinitely small region. We can use integration to find the probability.
CDF tells us the probability a random variable returns a value less than some specified value. It is the accumulation of the probability value up to the same value. $f(x) = Pr[X <= y]$
[!NOTE] Programmation - Cumulative Distribution Function
def plot_cumulative_normal_distribution(mu, variance): sigma = math.sqrt(variance) x = np.linspace(-10, 10, 100) plt.plot(x, stats.norm.cdf(x, mu, sigma), label=f"μ={mu}, σ²={variance}") plot_cumulative_normal_distribution(0, 1) plot_cumulative_normal_distribution(1, 2) plot_cumulative_normal_distribution(-3, 5) plt.legend() plt.show()![[Pasted_image_20251027104347.png]]
Distribution Normale et Théorème Central Limite
Standardisation (Z-Score)
[!NOTE] Définition Le Z-score transforme toute variable normale en loi normale standard N(0,1). The standard or z-score is the number of standard deviations by which the value of an observation is above or below the mean value of what is expected from its sampling distribution.
\[Z = \frac{X - \mu}{\sigma}\]Interprétation :
Z < 1 : valeur proche de la moyenne
1 < Z < 2 : valeur modérément éloignée
2 < Z < 3 : valeur très éloignée
Z > 3 : valeur exceptionnelle (outlier potentiel)
[!example]- Contrôle qualité spectroscopique Un colorant doit avoir ε = 50 000 ± 2000 L·mol⁻¹·cm⁻¹ (μ ± σ).
Un lot mesure ε = 54 500. Est-ce acceptable ?
\[Z = \frac{54500 - 50000}{2000} = 2.25\]P(Z > 2.25) = 0.0122, soit 1.22% de probabilité.
Ce lot est dans le top 1.2% mais reste dans les spécifications ( Z < 3).
Théorème Central Limite (TCL)
[!NOTE] Définition Pour n observations indépendantes de moyenne μ et variance σ², la moyenne échantillonnale suit approximativement :
\[\bar{X} \sim \mathcal{N}\left(\mu, \frac{\sigma^2}{n}\right)\]Conditions :
- Échantillons indépendants
- n ≥ 30 pour distributions non-normales
- n ≥ 10 pour distributions symétriques
Erreur standard : \(SE = \frac{\sigma}{\sqrt{n}}\)
Z-Score pour une moyenne d’échantillons : \(Z = \frac{\bar{X} - \mu}{\sigma / \sqrt{n}}\)
Implications :
- La distribution d’échantillonnage devient normale même si la population ne l’est pas
- La variance diminue avec la taille de l’échantillon
- Permet l’inférence statistique
[!example]- Moyennage de signal en spectroscopie Un spectromètre mesure un signal bruité avec σ = 10 mV. En moyennant n = 100 acquisitions :
\[SE = \frac{10}{\sqrt{100}} = 1 \text{ mV}\]Le rapport signal/bruit s’améliore d’un facteur √100 = 10.
Pour atteindre SE = 0.5 mV : \(n = \left(\frac{\sigma}{SE}\right)^2 = \left(\frac{10}{0.5}\right)^2 = 400 \text{ acquisitions}\)
Statistiques inférentielles
Dans le cadre du Théorème Central Limite (TCL) et du calcul des intervalles de confiance (IC), on cherche à estimer la moyenne réelle de la population à partir d’un échantillon. Le TCL nous assure que, pour un grand nombre d’observations indépendantes, la distribution des moyennes d’échantillons suit une loi normale, même si la population d’origine ne l’est pas.
Variables et Notations
| Symbole | Nom / Signification | Détermination |
|---|---|---|
| $\mu$ | Moyenne de la population (valeur réelle, inconnue) | Théorique — on cherche à l’estimer |
| $\overline{X}$ | Moyenne de l’échantillon | Calculée à partir des données : $\displaystyle \overline{X} = \frac{1}{n}\sum_{i=1}^{n} X_i$ |
| $\sigma^2$ | Variance de la population | Souvent inconnue — représente la dispersion réelle |
| $\sigma$ | Écart-type de la population | $\sqrt{\sigma^2}$, rarement connu en pratique |
| $s^2$ | Variance de l’échantillon | Estimation de $\sigma^2$ : $\displaystyle s^2 = \frac{1}{n-1}\sum (X_i - \overline{X})^2$ |
| $s$ | Écart-type de l’échantillon | $\sqrt{s^2}$, estimation de $\sigma$ |
| $n$ | Taille de l’échantillon | Nombre total d’observations |
| $z$ | Score normalisé (z-score) | Mesure combien d’écarts-types séparent $\overline{X}$ de $\mu$ |
| $t$ | Statistique de Student | Variante du z-score quand $\sigma$ est inconnu |
| $\sigma / \sqrt{n}$ ou $s / \sqrt{n}$ | Erreur standard de la moyenne | Mesure la dispersion attendue des moyennes d’échantillons |
Formules clés
Si $\sigma$ est connu
\(z = \frac{\overline{X} - \mu}{\sigma / \sqrt{n}}\)
Si $\sigma$ est inconnu
\(t = \frac{\overline{X} - \mu}{s / \sqrt{n}}\)
Ces scores servent à construire un intervalle de confiance : \(IC = \overline{X} \pm z_{\alpha/2} \cdot \frac{\sigma}{\sqrt{n}} \quad \text{ou} \quad IC = \overline{X} \pm t_{\alpha/2,\,n-1} \cdot \frac{s}{\sqrt{n}}\)
On peut se demander si à mesure que $n$ augmente, la distribution des moyennes d’échantillons ne tendrait pas toujours vers une forme normale, rendant ainsi $\overline{X}$ une approximation de plus en plus fiable de $\mu$ ?
Loi des grands nombres
When you consider independent random variables $X_1, \ldots, X_n$ with a common underlying probability distribution $p(\mu, \sigma)$, their average $\overline{X_n}$ becomes a strong approximation of $ \mu $ as the sample size $ n $ increases:
\[\overline{X_n} = \frac{X_1 + \ldots + X_n}{n} \xrightarrow[n \to \infty]{} \mu\]
Échantillonnage et Estimation
[!NOTE] Définition La distribution d’échantillonnage décrit la variabilité d’une statistique calculée sur différents échantillons de même taille tirés de la même population.
Pour la moyenne :
- Espérance : E($\bar{X}$) = μ
- Variance : Var($\bar{X}$) = σ²/n
- Forme : Normale si population normale ou n grand (TCL)
Pour la proportion :
- Espérance : E($\hat{p}$) = p
- Variance : Var($\hat{p}$) = p(1-p)/n
- Approximation normale si np ≥ 5 et n(1-p) ≥ 5
![[Pasted_image_20251029093346.png]]
[!example]- Taux de conversion photochimique Sur 200 molécules irradiées, 140 subissent la photoréaction (p̂ = 0.70).
Erreur standard : \(SE = \sqrt{\frac{0.70 \times 0.30}{200}} = 0.032\)
IC à 95% pour le vrai taux : \(0.70 \pm 1.96 \times 0.032 = [0.637, 0.763]\)
Le vrai taux de conversion est entre 63.7% et 76.3% avec 95% de confiance.
Variance and Standard Deviation
Variance of a Population
The variance of a population of $ N $ elements is:
\[\sigma^2 = \frac{1}{N} \sum_{i=1}^{N} (x_i - \mu)^2\]Standard Deviation of a Population
The standard deviation of a population of $ N $ elements is the square root of the variance:
\[\sigma = \sqrt{ \frac{1}{N} \sum_{i=1}^{N} (x_i - \mu)^2 }\]Sample Standard Deviation
Based on a sample of a population, a commonly used estimator of $ \sigma $ is the sample standard deviation:
\[s = \sqrt{ \frac{1}{n - 1} \sum_{i=1}^{n} (x_i - \overline{x})^2 }\]⚠️ Using $\frac{1}{n}$ instead of $\frac{1}{n-1}$ would underestimate the true population variance — this correction is known as Bessel’s correction.
Intervalle de Confiance
[!NOTE] Définition Un intervalle de confiance au niveau (1-α)% contient le paramètre vrai avec probabilité (1-α).
Pour la moyenne (σ connu) : \(IC_{1-\alpha} = \bar{x} \pm z_{\alpha/2} \frac{\sigma}{\sqrt{n}}\)
Pour la moyenne (σ inconnu) : \(IC_{1-\alpha} = \bar{x} \pm t_{\alpha/2,n-1} \frac{s}{\sqrt{n}}\)
Valeurs critiques usuelles :
- 90% : z = 1.645
- 95% : z = 1.96
- 99% : z = 2.576
[!example]- Rendement quantique moyen Sur 12 mesures de Φ_F, on obtient $\bar{x}$ = 0.65, s = 0.08.
IC à 95% :
- t₀.₀₂₅,₁₁ = 2.201
- Marge d’erreur = 2.201 × 0.08/√12 = 0.051
- IC₉₅% = [0.65 - 0.051, 0.65 + 0.051] = [0.599, 0.701]
On peut affirmer avec 95% de confiance que le vrai rendement quantique est entre 0.599 et 0.701.
Interprétation de l’IC
If we were to repeat this process and construct many other samples, 95% of the intervals produced will actually contain the true US mean pop height We’re 95% confident that [168.7 - 171.2] captures the true average height.
![[Pasted_image_20251029094444.png]]
==Don’t say “there is a 95% probability that μ__is between …” because the real μ isn’t random!==
Programmation
# We can check these figure using a Cumulative Density Function `cdf`
from scipy import stats
mu_estim = stats.norm(170, 0.63)
# use the cdf to find the probabilities associated with height values
print('% confidence interval = ', round(mu_estim.cdf(171.2) - mu_estim.cdf(168.7),2))
% confidence interval = 0.95
💡 Actually, there is a formula to find the lower bound and the upper bound of any confidence interval
(ex: 99%) ```python confidence_interval = 0.99
sup_proba = (1 + confidence_interval)/2 # 99.5% inf_proba = (1 - confidence_interval)/2 # 0.5%
mu_upper_bound = mu_estim.ppf(sup_proba) mu_lower_bound = mu_estim.ppf(inf_proba)
use the inverse of the cdf to find the heights associated with probabilities
print(‘mu_upper_bound: ‘, mu_upper_bound) print(‘mu_lower_bound: ‘, mu_lower_bound)
print(‘% confidence interval = ‘, round(mu_estim.cdf(mu_upper_bound) - mu_estim.cdf(mu_lower_bound),2))
>mu_upper_bound: 171.6227724612358
>
>mu_lower_bound: 168.3772275387642
>
>% confidence interval = 0.99
---
### Human perception
Here are some usual English names for confidence intervals
_(according to the Intergovernmental Panel for Climate Change: [IPCC](https://archive.ipcc.ch/publications_and_data/ar4/wg1/en/ch1s1-6.html) for instance)_
- 1-sigma (68%) `"likely"`
- 90% `"very likely"`
- 2-sigma (95%) `"extremely likely"`
- 3-sigma (99.7%) `"virtually certain"`
- 5-sigma: "proof" `threshold in theoretical physics`
![[Pasted_image_20251029095347.png]]
---
### Sample size `n` considerations
❓ When is n considered **large enough** for the CLT
Three cases:
- If n>30 ⇒ CLT applies **and** the sample std s can be used to approximate the true pop σ- in z-statistics
- If n>10 _and_ observations are "non-skewed" and without outliers ⇒- CLT still applies
- If the global population is known to be normally distributed ⇒ CLT always applies even with an arbitrary small n, in a sense that we can use the Gaussian distribution
❓ When is n is **small enough** to consider each draw independent, even without replacement
- n<10%×N
## Test d'Hypothèse
### Concepts Fondamentaux
> [!NOTE] Définition
> Un test d'hypothèse est une procédure statistique permettant de prendre une décision sur une propriété de la population à partir d'un échantillon.
>
> **Éléments clés :**
> - **Hypothèse nulle (H₀)** : affirmation de base à tester (pas d'effet, pas de différence)
> - **Hypothèse alternative (H₁)** : ce qu'on cherche à démontrer
> - **Statistique de test** : valeur calculée à partir des données
> - **Seuil de signification (α)** : probabilité de rejeter H₀ à tort (généralement 0.05)
> - **p-value** : probabilité d'obtenir un résultat au moins aussi extrême sous H₀
> - **Puissance (1-β)** : probabilité de détecter un effet s'il existe vraiment
>
>**Types d'erreurs :**
> - **Erreur de type I (α)** : Rejeter H₀ alors qu'elle est vraie (faux positif)
> - **Erreur de type II (β)** : Ne pas rejeter H₀ alors qu'elle est fausse (faux négatif)
>
> **Décision :**
> - Si p-value < α : on rejette H₀
> - Si p-value ≥ α : on ne rejette pas H₀ (on ne l'accepte pas pour autant)
> [!example]- Meilleure rétention avec l’ajout d’une feature
> ![[Pasted_image_20251029100539.png|500]]
>
> 1. Create Null Hypothesis $H_0$
> $$H_0: \mu = 300 \quad \text{(unchanged) in dark mode}$$
>
> 2. Create Alternative Hypothesis $H_a$
> $$H_a: \mu > 300 \quad \text{(increased) in dark mode}$$
>
> 3. Choose a Significance Level $alpha$
> $$\text{Example: } \alpha = 5\%$$
>
> 4. Compute p-value
> Suppose that $H_0$ is true, and compute the probability of observing a sample mean $overline{X}_n \ge 310$:
> $$\text{p-value} = P(\overline{X}_n \ge 310 \mid H_0)$$
> - If $\text{p-value} < \alpha$, reject $H_0$ in favor of $H_a$.
> - If $\text{p-value} > \alpha$, fail to reject $H_0$ (this does **not** mean we accept $H_0$).
>
> 5. Compute the p-value using CLT
> Since $n = 100$ is large enough, the Central Limit Theorem (CLT) applies, telling us that the distribution of sample means $overline{X}_n$ should follow a normal distribution:
> $$\overline{X}_n \approx N\left(\mu, \frac{\sigma}{\sqrt{n}}\right)$$
>
> Supposing $H_0$ is true ($mu = 300$, $sigma = 50$):
> $$\overline{X}_{100} \approx N\left(300, \frac{50}{\sqrt{100}}\right)$$
>
> ![[Pasted_image_20251029101453.png|500]]
>
> ```python
> from scipy.stats import norm
> X = norm(300, 50/(100**0.5))
> p_value = (1 - X.cdf(310));
> round(p_value,2)
> ```
> > 0.02
>
> `p-value` < 0.05 → We can safely reject the Null Hypothesis H0 ⇒ The dark mode is a real plus!
> [!example]- Validation d'un nouveau photocatalyseur
> Un nouveau photocatalyseur dopé au platine prétend améliorer le rendement de dégradation.
>
> - H₀ : μ_nouveau = μ_standard = 75% (pas d'amélioration)
> - H₁ : μ_nouveau > μ_standard (amélioration)
> - α = 0.05
>
> Après 10 tests, on obtient $\bar{x} = 82$%, s = 5%, p-value = 0.002
>
> Comme p-value (0.002) < α (0.05), on rejette H₀. Le nouveau catalyseur est significativement plus efficace.
You will also often encounter the **power** of a statistical test (the larger the better):
- Power is the probability that we will **correctly reject the null hypothesis** (if it was correct to reject it, so the probability of _not missing out_ a great feature in A/B testing or _not missing out_ an effective new drug in clinical trial)
- It correspond to P(not making a type II error)
### P-hacking
- Tester **toutes les variables possibles** puis ne publier que les résultats “significatifs”.
- Risque : obtenir une p-value < 0.05 par hasard.
- Bonne pratique : tester une hypothèse spécifique **avant de collecter les données**.
> [!example]-
> - Tester 20 solvants différents pour un produit.
> - 5% des tests pourraient montrer un rendement “significatif” par hasard → **p-hacking**.
## Test t de Student
> [!NOTE] Définition
> Le test t compare des moyennes quand la variance de la population est inconnue et/ou l'échantillon est petit (n < 30).
>
> **Test t pour un échantillon :**
> $$t = \frac{\bar{x} - \mu_0}{s/\sqrt{n}}$$
>
> **Test t pour deux échantillons indépendants :**
> $$t = \frac{\bar{x}_1 - \bar{x}_2}{\sqrt{\frac{s_1^2}{n_1} + \frac{s_2^2}{n_2}}}$$
>
> **Test t apparié :**
> $$t = \frac{\bar{d}}{s_d/\sqrt{n}}$$
>
> où $\bar{d}$ est la moyenne des différences appariées.
>
> **Degrés de liberté :**
> - Un échantillon : df = n - 1
> - Deux échantillons : df ≈ n₁ + n₂ - 2 (ou formule de Welch si variances inégales)
> - Échantillons appariés : df = n - 1
>
> ![[Pasted_image_20251029102521.png|300]]
> [!example]- Comparaison de méthodes de synthèse
> On compare deux méthodes de synthèse d'un fluorophore organique (rendements en %) :
>
> **Méthode A (thermique)** : n₁ = 8, $\bar{x}_1$ = 65%, s₁ = 4%
> **Méthode B (photochimique)** : n₂ = 8, $\bar{x}_2$ = 72%, s₂ = 5%
>
> $$t = \frac{72 - 65}{\sqrt{\frac{25}{8} + \frac{16}{8}}} = \frac{7}{\sqrt{5.125}} = 3.09$$
>
> Avec df = 14 et α = 0.05 (bilatéral), t_critique = 2.145
>
> Comme |t| = 3.09 > 2.145, la différence est significative. La méthode photochimique donne un meilleur rendement.
---
## Test du Chi-deux (χ²)
> [!NOTE] Définition
> Le test du chi-deux évalue l'indépendance entre variables catégorielles ou la conformité à une distribution théorique.
>
> **Statistique de test :**
> $$\chi^2 = \sum_{i} \frac{(O_i - E_i)^2}{E_i}$$
>
> où :
> - $O_i$ : fréquence observée
> - $E_i$ : fréquence attendue sous H₀
>
> **Degrés de liberté :**
> - Test d'ajustement : df = k - 1 - p (k catégories, p paramètres estimés)
> - Test d'indépendance : df = (r-1)(c-1) (r lignes, c colonnes)
> [!example]- Sélectivité de photoréaction
> Une photoréaction peut donner 3 produits. Théoriquement, les proportions attendues sont 50%, 30%, 20%.
> Sur 200 réactions :
>
> | Produit | Observé | Attendu | $(O-E)^2/E$ |
> |---------|---------|---------|-------------|
> | A | 110 | 100 | 1.00 |
> | B | 55 | 60 | 0.42 |
> | C | 35 | 40 | 0.63 |
>
> $$\chi^2 = 1.00 + 0.42 + 0.63 = 2.05$$
>
> Avec df = 3-1 = 2 et α = 0.05, χ²_critique = 5.99
>
> Comme χ² = 2.05 < 5.99, on ne rejette pas H₀. La distribution observée est conforme à la théorie.
---
## ANOVA
> [!NOTE] Définition
> L'ANOVA,analyse la variance, compare les moyennes de plus de deux groupes en décomposant la variance totale en variance inter-groupes et intra-groupes.
>
> **Statistique F :**
> $$F = \frac{\text{MS}_{\text{between}}}{\text{MS}_{\text{within}}} = \frac{\text{Variance inter-groupes}}{\text{Variance intra-groupes}}$$
>
> **Conditions d'application :**
> - Normalité dans chaque groupe
> - Homoscédasticité (variances égales)
> - Indépendance des observations
>
> **Degrés de liberté :**
> - Inter-groupes : df₁ = k - 1
> - Intra-groupes : df₂ = n - k
> [!example]- Optimisation de conditions photocatalytiques
> On teste l'effet du pH (4, 7, 10) sur l'efficacité de dégradation (%) :
>
> | pH 4 | pH 7 | pH 10 |
> |------|------|-------|
> | 45 | 75 | 60 |
> | 48 | 78 | 62 |
> | 43 | 73 | 58 |
> | 46 | 76 | 61 |
> | 47 | 77 | 59 |
>
> Moyennes : $\bar{x}_4$ = 45.8%, $\bar{x}_7$ = 75.8%, $\bar{x}_{10}$ = 60%
>
> Après calculs :
> - SS_between = 2346.4
> - SS_within = 54.8
> - F = (2346.4/2)/(54.8/12) = 257.1
>
> Avec F_critique(2,12,0.05) = 3.89
>
> F = 257.1 >> 3.89, donc l'effet du pH est hautement significatif.
---
## Degrés de Liberté
> [!NOTE] Définition
> Les degrés de liberté (df) représentent le nombre de valeurs indépendantes pour estimer un paramètre.
>
> **Formules courantes :**
> - Variance d'échantillon : df = n - 1
> - Test t deux échantillons : df ≈ n₁ + n₂ - 2
> - Régression linéaire simple : df = n - 2
> - Régression multiple (k variables) : df = n - k - 1
> - ANOVA à k groupes : df_between = k - 1, df_within = n - k
> - Test χ² d'ajustement : df = catégories - 1 - paramètres estimés
> - Test χ² d'indépendance : df = (r-1)(c-1)
>
> **Intuition :**
> Chaque paramètre estimé "consomme" un degré de liberté car il impose une contrainte sur les données.
> [!example]- Ajustement cinétique complexe
> On ajuste une cinétique bi-exponentielle A(t) = A₁e^(-k₁t) + A₂e^(-k₂t) sur 20 points.
>
> Paramètres à estimer : A₁, A₂, k₁, k₂ (4 paramètres)
>
> df = 20 - 4 = 16
>
> Ces 16 degrés de liberté déterminent la distribution du χ² pour évaluer la qualité de l'ajustement.
> [!example]- Correction de Welch-Satterthwaite
> Pour deux échantillons de variances inégales :
>
> $$df = \frac{\left(\frac{s_1^2}{n_1} + \frac{s_2^2}{n_2}\right)^2}{\frac{(s_1^2/n_1)^2}{n_1-1} + \frac{(s_2^2/n_2)^2}{n_2-1}}$$
>
> Avec n₁=10, s₁²=4, n₂=15, s₂²=9 :
>
> df ≈ 22.4 → arrondi à 22
---
## [[Théorème de Bayes et Vraisemblance]]
### Théorème de Bayes
> [!NOTE] Définition
> Le théorème de Bayes met à jour la probabilité d'une hypothèse à la lumière de nouvelles données.
> Le but est de voir l’indépendance des événements, que AsachantB vaut B sachant A.
>
> **Forme simple :**
> $$P(H|E) = \frac{P(E|H) \cdot P(H)}{P(E)}$$
>
> **Forme développée :**
> $$P(H|E) = \frac{P(E|H) \cdot P(H)}{P(E|H) \cdot P(H) + P(E|\neg H) \cdot P(\neg H)}$$
>
> où :
> - P(H) : probabilité a priori
> - P(E|H) : vraisemblance
> - P(H|E) : probabilité a posteriori
> - P(E) : évidence (normalisation)
L'inférence statistique repose sur une inversion conceptuelle essentielle appelée le "flip logique".
En **probabilité classique**, on fixe un paramètre θ et on calcule la probabilité d'observer différents résultats Y.
En **vraisemblance**, on inverse cette logique : les données Y sont fixes, et on évalue la plausibilité de différentes valeurs de θ.
### Vraisemblance (Likelihood)
> [!NOTE] Définition
> La vraisemblance L(θ|x) mesure la plausibilité des paramètres θ étant donné les données x.
>
> **Pour n observations indépendantes :**
> $$L(\theta|x_1,...,x_n) = \prod_{i=1}^n f(x_i|\theta)$$
>
> **Log-vraisemblance (plus pratique) :**
> $$\ell(\theta) = \ln L(\theta) = \sum_{i=1}^n \ln f(x_i|\theta)$$
Mathématiquement, L(θ|Y) = P(Y|θ), mais l'interprétation change radicalement.
→ Si la prévalence d'une maladie est θ = 7% et qu'on observe 6 cas sur 100 personnes, L(θ=0.07|Y=6) = 0.153.
Cette valeur absolue n'a aucun sens isolément, mais le ratio de vraisemblance devient informatif.
→ Si L(0.07)/L(0.08) = 1.24, alors θ = 7% est 1.24 fois plus plausible.
Pour n observations indépendantes, la vraisemblance est le produit des probabilités individuelles.
On prend systématiquement le logarithme pour transformer ces produits en sommes, évitant les problèmes numériques et simplifiant la dérivation. Le logarithme étant monotone croissant, maximiser L(θ) ou ln L(θ) revient au même.
Le numérateur du théorème de Bayes est exactement la vraisemblance. La différence : Bayes produit une probabilité normalisée en divisant par l'évidence P(E), tandis que la vraisemblance reste non normalisée puisque seuls les ratios comptent.
### Estimateur du Maximum de Vraisemblance (MLE)
L'estimateur du maximum de vraisemblance (MLE) trouve les paramètres rendant les données les plus plausibles. La procédure : écrire la log-vraisemblance, dériver, égaler à zéro, résoudre. Pour une normale avec variance connue, le MLE de la moyenne est la moyenne empirique. Le MLE possède des propriétés asymptotiques remarquables : consistance (converge vers la vraie valeur), normalité asymptotique, et efficacité (atteint la borne de Cramér-Rao, donc variance minimale). Cette efficacité fait du MLE un gold standard.
> [!NOTE] Définition
> Le MLE trouve les paramètres qui maximisent la vraisemblance des données observées.
>
> **Procédure :**
> 1. Écrire la vraisemblance L(θ)
> 2. Prendre le logarithme : ℓ(θ)
> 3. Dériver : ∂ℓ/∂θ
> 4. Égaler à zéro et résoudre
> 5. Vérifier que c'est un maximum
>
> **Propriétés asymptotiques :**
> - Consistance : θ̂ₙ → θ₀
> - Normalité : √n(θ̂ₙ - θ₀) → N(0, I⁻¹(θ₀))
> - Efficacité : variance minimale
### Statistiques Suffisantes et Paramètres Multiples
Un concept puissant émerge avec les statistiques suffisantes : on n'a pas toujours besoin de toutes les observations. Une statistique T(Y) est suffisante pour θ si elle contient toute l'information pertinente. Le critère de factorisation stipule que T est suffisante si la vraisemblance s'écrit comme produit d'une fonction de T et θ, et d'une fonction des données seules. Pour une normale avec variance connue, la somme des observations suffit : au lieu de n valeurs, un seul nombre ! Avec deux paramètres inconnus, il faut deux statistiques : somme des observations et somme des carrés. En photochimie, des milliers de spectres se résument ainsi à quelques statistiques.
Quand on estime plusieurs paramètres, la vraisemblance devient multivariée. Pour un paramètre, on a une courbe 2D ; pour deux, une surface 3D dont le sommet représente les MLE. Les courbes de niveau définissent des régions de confiance. Au-delà de trois paramètres, impossible de visualiser, mais le concept reste : chercher le maximum d'une hypersurface. La complexité croît exponentiellement, nécessitant des méthodes numériques comme Newton-Raphson ou l'algorithme EM.
### Profile likelihood et Paramètres Nuisibles
> [!NOTE] Définition
> Pour un modèle à plusieurs paramètres θ = (ψ, λ) où ψ est d'intérêt et λ est nuisible :
>
> **Profile likelihood :**
> $$L_p(\psi) = \max_\lambda L(\psi, \lambda) = L(\psi, \hat{\lambda}(\psi))$$
>
> **Intervalle de confiance par profile likelihood :**
> $$\{ψ : 2[\ell(\hat{\psi}) - \ell_p(\psi)] < \chi^2_{1,\alpha}\}$$
Dans de nombreuses situations, un paramètre ψ nous intéresse (constante de vitesse), tandis que d'autres λ (amplitude, bruit) sont des nuisances.
Le profile likelihood résout cela élégamment : pour chaque ψ, on maximise par rapport à λ, obtenant λ̂(ψ).
Ainsi Lₚ(ψ) = L(ψ, λ̂(ψ)) concentre toute l'information sur ψ.
C'est comme prendre une coupe optimale de la surface 3D.
Pour estimer μ avec σ² inconnu, on calcule σ²(μ) = (1/n)Σ(yᵢ-μ)² pour chaque μ, puis Lₚ(μ) devient fonction uniquement de μ. Cette méthode est robuste et sans perturbation des paramètres secondaires.
### Information de Fisher
L'information de Fisher I(θ) quantifie l'information contenue sur le paramètre. Information élevée = log-vraisemblance variant fortement = estimation précise. La borne de Cramér-Rao : variance minimale = 1/[nI(θ)], atteinte asymptotiquement par le MLE. L'information s'additionne : Iₙ(θ) = n·I(θ). L'erreur standard SE(θ̂) ≈ 1/√[nI(θ)] permet de planifier les expériences : pour mesurer un temps de fluorescence avec 1% de précision, il faut au moins 10,000 photons, illustrant le trade-off précision-temps-photoblanchiment.
> [!NOTE] Définition
> L'information de Fisher quantifie l'information sur θ contenue dans les données.
>
> **Définition :**
> $$I(\theta) = E\left[\left(\frac{\partial \ln f(X|\theta)}{\partial \theta}\right)^2\right] = -E\left[\frac{\partial^2 \ln f(X|\theta)}{\partial \theta^2}\right]$$
>
> **Borne de Cramér-Rao :**
> $$\text{Var}(\hat{\theta}) \geq \frac{1}{nI(\theta)}$$
>
> Un estimateur atteignant cette borne est dit efficace.
Données Y observées ↓ [Le flip logique] ↓ Fonction de vraisemblance L(θ|Y) ↓ [Log-transformation] ↓ Log-vraisemblance ℓ(θ) ↓ [Statistiques suffisantes] ← Réduction dimensionnelle ↓ Maximisation ↓ MLE θ̂ ← Estimateur optimal ↓ Information de Fisher I(θ) ← Précision ↓ SE(θ̂) = 1/√[nI(θ)] ← Incertitude ↓ Intervalles de confiance
---
# Régressions
La régression constitue l'outil statistique appliqué le plus important en science des données. Son principal avantage réside dans son **interprétabilité** : contrairement aux modèles de type "boîte noire", les coefficients de régression possèdent une signification physique directe. En chimie quantitative, on peut se demander si la régression n'est pas l'équivalent moderne de la loi de Beer-Lambert : une relation linéaire fondamentale entre cause et effet.
La régression linéaire trouve ses racines dans les travaux de Legendre (1805) et Gauss (1809) sur la méthode des moindres carrés. Initialement développée pour l'astronomie, elle s'applique aujourd'hui à tous les domaines où des relations quantitatives doivent être établies. En photochimie, par exemple, elle permet de relier l'intensité lumineuse au rendement quantique, ou la concentration en catalyseur à la vitesse de réaction.
> [!bug] **Applications en chimie**
> - Calibration spectroscopique (loi de Beer-Lambert)
> - Relations structure-activité (QSAR)
> - Cinétique chimique (équations de vitesse)
> - Thermochimie (corrélations linéaires d'énergie libre)
> - Analyse quantitative (courbes d'étalonnage)
---
## Régression Linéaire Simple
La régression linéaire simple modélise la relation entre une variable explicative $X$ (prédicteur) et une variable réponse $Y$ (variable dépendante) par une fonction affine.
### Fondements mathématiques
> [!NOTE] Définition
> La régression linéaire simple cherche à modéliser la relation entre deux variables par une droite :
>
> $$Y = \beta_0 + \beta_1 X + \varepsilon$$
>
> où :
> - $Y$ : variable réponse (dépendante)
> - $X$ : variable explicative (indépendante, prédicteur)
> - $\beta_0$ : ordonnée à l'origine (intercept)
> - $\beta_1$ : pente (slope), coefficient directeur
> - $\varepsilon$ : terme d'erreur aléatoire, $\varepsilon \sim \mathcal{N}(0, \sigma^2)$
>
> **Hypothèses du modèle :**
> 1. **Linéarité** : La relation entre $X$ et $Y$ est linéaire
> 2. **Indépendance** : Les observations sont indépendantes
> 3. **Homoscédasticité** : La variance des erreurs est constante
> 4. **Normalité** : Les erreurs suivent une loi normale
> 5. **Absence de multicolinéarité** : (pour la régression multiple)
**Interprétation physique :**
- $\beta_0$ : valeur de $Y$ lorsque $X = 0$ (pas toujours physiquement significatif)
- $\beta_1$ : variation de $Y$ pour une augmentation unitaire de $X$
- $\sigma^2$ : variance inexpliquée, représente le bruit expérimental et les facteurs non modélisés
> [!example]- Loi de Beer-Lambert
> La concentration d'un chromophore en solution est déterminée par spectroscopie UV-visible. L'absorbance $A$ est reliée à la concentration $c$ par :
>
> $$A = \varepsilon \ell c$$
>
> où $\varepsilon$ est le coefficient d'extinction molaire et $\ell$ la longueur du trajet optique.
>
> Mesures expérimentales pour un dérivé de porphyrine à 420 nm ($\ell = 1$ cm) :
>
> | Concentration (μM) | Absorbance |
> |-------------------|------------|
> | 1.0 | 0.052 |
> | 2.0 | 0.098 |
> | 3.0 | 0.151 |
> | 4.0 | 0.195 |
> | 5.0 | 0.248 |
>
> Régression linéaire : $A = 0.003 + 0.049 \times c$
>
> - $\beta_0 = 0.003$ : absorbance résiduelle (impuretés, ligne de base)
> - $\beta_1 = 0.049$ : correspond à $\varepsilon \ell = 49000$ L·mol⁻¹·cm⁻¹
---
### Méthode des Moindres Carrés Ordinaires (OLS)
> [!NOTE] Définition
> La méthode des Moindres Carrés Ordinaires (Ordinary Least Squares, OLS) minimise la somme des carrés des résidus (Sum of Squared Residuals, SSR) :
>
> $$\text{SSR} = \sum_{i=1}^{n} (y_i - \hat{y}_i)^2 = \sum_{i=1}^{n} (y_i - \beta_0 - \beta_1 x_i)^2$$
>
> où $\hat{y}_i = \beta_0 + \beta_1 x_i$ est la valeur prédite.
>
> **Solutions analytiques :**
>
> $$\beta_1 = \frac{\sum_{i=1}^{n} (x_i - \bar{x})(y_i - \bar{y})}{\sum_{i=1}^{n} (x_i - \bar{x})^2} = \frac{\text{Cov}(X,Y)}{\text{Var}(X)}$$
>
> $$\beta_0 = \bar{y} - \beta_1 \bar{x}$$
>
> **Propriétés :**
> - Estimateurs non biaisés : $\mathbb{E}[\hat{\beta}_0] = \beta_0$ et $\mathbb{E}[\hat{\beta}_1] = \beta_1$
> - Variance minimale parmi tous les estimateurs linéaires non biaisés (théorème de Gauss-Markov)
> - La droite de régression passe par le point $(\bar{x}, \bar{y})$
> - La somme des résidus est nulle : $\sum_{i=1}^{n} (y_i - \hat{y}_i) = 0$
**Pourquoi élever au carré ?**
On peut se demander si d'autres métriques de distance seraient préférables. Trois raisons justifient l'utilisation du carré :
1. **Symétrie des erreurs** : Les écarts positifs et négatifs sont traités de manière équivalente
2. **Amplification des outliers** : Les valeurs aberrantes sont fortement pénalisées, ce qui est souhaitable pour détecter les anomalies expérimentales
3. **Propriétés mathématiques** : Le carré est dérivable partout, facilitant l'optimisation analytique
La distance euclidienne (norme $L_2$) émerge naturellement de ces considérations. On pourrait utiliser d'autres normes ($L_1$ pour la régression robuste, $L_\infty$ pour minimiser l'erreur maximale), mais la norme $L_2$ possède les meilleures propriétés statistiques sous hypothèse de normalité.
> [!example]- Cinétique d'ordre un
> On étudie la photolyse d'un colorant organique sous irradiation UV. La concentration résiduelle suit une cinétique d'ordre un :
>
> $$\ln[C]_t = \ln[C]_0 - k_{app} t$$
>
> Mesures expérimentales :
>
> | Temps (min) | [C] (μM) | ln[C] |
> |------------|----------|-------|
> | 0 | 100.0 | 4.605 |
> | 5 | 77.8 | 4.354 |
> | 10 | 60.7 | 4.106 |
> | 15 | 47.2 | 3.855 |
> | 20 | 36.8 | 3.606 |
> | 25 | 28.7 | 3.357 |
>
> Régression linéaire de $\ln[C]$ vs $t$ :
>
> $$\ln[C]_t = 4.608 - 0.0503 \times t$$
>
> - $\beta_0 = 4.608 \approx \ln(100.4)$ : concentration initiale ajustée
> - $\beta_1 = -0.0503$ min⁻¹ : constante de vitesse apparente $k_{app}$
> - $t_{1/2} = \frac{\ln 2}{k_{app}} = 13.8$ min
---
### Coefficient de Détermination (R²)
> [!NOTE] Définition
> Le coefficient de détermination $R^2$ mesure la proportion de variance de $Y$ expliquée par le modèle :
>
> $$R^2 = 1 - \frac{\text{SSR}}{\text{SST}} = 1 - \frac{\sum_{i=1}^{n}(y_i - \hat{y}_i)^2}{\sum_{i=1}^{n}(y_i - \bar{y})^2}$$
>
> où :
> - $\text{SSR}$ : somme des carrés des résidus (unexplained variance)
> - $\text{SST}$ : somme totale des carrés (total variance)
> - $\text{SSE} = \text{SST} - \text{SSR}$ : somme expliquée des carrés
>
> **Forme alternative :**
>
> $$R^2 = \frac{\text{SSE}}{\text{SST}} = \left(\frac{\text{Cov}(X,Y)}{\sigma_X \sigma_Y}\right)^2 = \rho^2$$
>
> où $\rho$ est le coefficient de corrélation de Pearson.
>
> **Interprétation :**
> - $R^2 = 0$ : le modèle n'explique rien (prédiction = moyenne)
> - $R^2 = 1$ : ajustement parfait (tous les points sur la droite)
> - $0 < R^2 < 1$ : qualité intermédiaire
> - En régression simple : $R^2 = r^2$ (carré du coefficient de corrélation)
**Limites du R² :**
On peut se demander si un $R^2$ élevé garantit un bon modèle. La réponse est non ! Le $R^2$ peut être trompeur dans plusieurs situations :
- **Données non linéaires** : un $R^2$ faible peut indiquer une relation non linéaire forte
- **Extrapolation** : un bon ajustement local ne garantit pas la validité hors domaine
- **Corrélation vs causalité** : $R^2$ élevé ≠ relation causale
- **Overfitting** : en régression multiple, on peut toujours augmenter $R^2$ en ajoutant des variables
> [!example]- Relation structure-activité
> On étudie le rendement quantique de fluorescence $\Phi_F$ de 15 dérivés de coumarine en fonction de leur moment dipolaire $\mu$ (Debye) :
>
> Régression : $\Phi_F = 0.82 - 0.068 \times \mu$
>
> - $R^2 = 0.73$ : 73% de la variance du rendement quantique est expliquée par le moment dipolaire
> - 27% de variance résiduelle provient d'autres facteurs : effets stériques, couplage vibronique, interactions solvant-soluté
>
> **Interprétation physique :** L'augmentation du moment dipolaire favorise le transfert de charge intramoléculaire (ICT), ouvrant un canal de désactivation non radiatif en compétition avec la fluorescence.
---
### Coefficient de Corrélation de Pearson
> [!NOTE] Définition
> Le coefficient de corrélation de Pearson $\rho$ (ou $r$ pour un échantillon) mesure l'intensité de la relation linéaire entre deux variables :
>
> $$r = \frac{\sum_{i=1}^{n}(x_i - \bar{x})(y_i - \bar{y})}{\sqrt{\sum_{i=1}^{n}(x_i - \bar{x})^2} \sqrt{\sum_{i=1}^{n}(y_i - \bar{y})^2}} = \frac{\text{Cov}(X,Y)}{\sigma_X \sigma_Y}$$
>
> **Propriétés :**
> - $r \in [-1, 1]$
> - $r = +1$ : corrélation linéaire positive parfaite
> - $r = -1$ : corrélation linéaire négative parfaite
> - $r = 0$ : absence de corrélation linéaire (mais peut exister une relation non linéaire !)
> - Invariant par transformation affine
> - Symétrique : $r_{X,Y} = r_{Y,X}$
>
> **Relation avec la régression :**
> $$r = \text{signe}(\beta_1) \times \sqrt{R^2}$$
**Interprétation prudente :**
On peut se demander si une corrélation forte implique nécessairement une causalité. C'est l'un des pièges les plus dangereux en statistique ! La corrélation détecte une covariation, pas une relation de cause à effet. Trois scénarios sont possibles :
1. $X$ cause $Y$ (causalité directe)
2. $Y$ cause $X$ (causalité inverse)
3. $Z$ cause $X$ et $Y$ (variable confondante)
L'exemple classique : la vente de glaces et le nombre de noyades sont fortement corrélés, mais ni l'un ne cause l'autre. La température (variable confondante) explique les deux phénomènes.
> [!example]- Matrice de corrélation spectroscopique
> Pour un ensemble de chromophores organiques, on mesure plusieurs propriétés photophysiques :
>
> | | $\lambda_{abs}$ | $\varepsilon$ | $\Phi_F$ | $\tau_F$ |
> |--|----------|-------|-------|-------|
> | $\lambda_{abs}$ | 1.00 | -0.12 | -0.85 | 0.78 |
> | $\varepsilon$ | -0.12 | 1.00 | 0.23 | -0.18 |
> | $\Phi_F$ | -0.85 | 0.23 | 1.00 | -0.91 |
> | $\tau_F$ | 0.78 | -0.18 | -0.91 | 1.00 |
>
> **Observations :**
> - Forte anticorrélation $\lambda_{abs}$ ↔ $\Phi_F$ ($r = -0.85$) : les chromophores absorbant dans le rouge ont généralement des rendements quantiques plus faibles (gap énergétique réduit, couplage vibronique accru)
> - Forte anticorrélation $\Phi_F$ ↔ $\tau_F$ ($r = -0.91$) : attendue par la relation $\Phi_F = \frac{k_r}{k_r + k_{nr}}$ et $\tau_F = \frac{1}{k_r + k_{nr}}$
> - Faible corrélation $\varepsilon$ avec les autres : le coefficient d'extinction dépend principalement de facteurs électroniques locaux (moment de transition) plutôt que de la dynamique de l'état excité
---
### Inférence Statistique sur les Coefficients
> [!NOTE] Définition
> Les coefficients de régression sont des estimateurs qui possèdent leurs propres distributions d'échantillonnage.
>
> **Écarts-types des coefficients :**
>
> $$SE(\beta_1) = \frac{s}{\sqrt{\sum_{i=1}^{n}(x_i - \bar{x})^2}} = \frac{s}{\sigma_X \sqrt{n-1}}$$
>
> $$SE(\beta_0) = s\sqrt{\frac{1}{n} + \frac{\bar{x}^2}{\sum_{i=1}^{n}(x_i - \bar{x})^2}}$$
>
> où $s = \sqrt{\frac{\text{SSR}}{n-2}}$ est l'écart-type résiduel.
>
> **Test de Student pour $\beta_1$ :**
>
> Hypothèses : $H_0: \beta_1 = 0$ vs $H_1: \beta_1 \neq 0$
>
> $$t = \frac{\hat{\beta}_1}{SE(\beta_1)} \sim t_{n-2}$$
>
> **Intervalle de confiance à (1-α)% :**
>
> $$\beta_1 \in \left[\hat{\beta}_1 - t_{n-2,\alpha/2} \times SE(\beta_1), \hat{\beta}_1 + t_{n-2,\alpha/2} \times SE(\beta_1)\right]$$
**Interprétation des p-values :**
La p-value associée à $\beta_1$ répond à la question : "Quelle est la probabilité d'observer une pente au moins aussi extrême si la vraie relation était nulle ?"
- $p < 0.001$ : évidence très forte d'une relation
- $p < 0.05$ : évidence significative (seuil conventionnel)
- $p > 0.05$ : absence d'évidence (≠ preuve d'absence !)
On peut se demander si le seuil de 0.05 est toujours pertinent. En physique des particules, on exige $p < 3 \times 10^{-7}$ (5σ) pour une découverte ! Le seuil dépend du contexte et des conséquences d'une erreur.
> [!example]- Test de significativité d'une constante de vitesse
> On étudie l'effet de la concentration en photocatalyseur [Cat] sur la constante de vitesse apparente d'une réaction de dégradation.
>
> Régression : $k_{app} = 0.012 + 0.087 \times [\text{Cat}]$
>
> Résultats du test :
> - $SE(\beta_1) = 0.0043$ min⁻¹·mM⁻¹
> - $t = \frac{0.087}{0.0043} = 20.2$
> - Avec $df = n-2 = 18$, la valeur critique $t_{18, 0.025} = 2.10$
> - $p < 0.0001$
>
> **Conclusion :** La concentration en catalyseur a un effet hautement significatif sur la vitesse de réaction. L'intervalle de confiance à 95% pour $\beta_1$ est [0.078, 0.096] min⁻¹·mM⁻¹.
---
### Analyse des Résidus
> [!NOTE] Définition
> Les résidus sont les écarts entre valeurs observées et prédites :
>
> $$e_i = y_i - \hat{y}_i = y_i - (\beta_0 + \beta_1 x_i)$$
>
> **Résidus standardisés :**
>
> $$r_i = \frac{e_i}{s\sqrt{1 - h_i}}$$
>
> où $h_i$ est le levier (leverage) de l'observation $i$.
>
> **Propriétés des résidus :**
> - $\sum_{i=1}^{n} e_i = 0$
> - $\sum_{i=1}^{n} e_i x_i = 0$
> - Si le modèle est correct : $e_i \sim \mathcal{N}(0, \sigma^2)$ indépendants
**Graphiques diagnostiques :**
1. **Résidus vs valeurs ajustées** : détecte l'hétéroscédasticité et la non-linéarité
2. **QQ-plot des résidus** : vérifie la normalité
3. **Résidus vs ordre temporel** : détecte l'autocorrélation
4. **Distance de Cook** : identifie les observations influentes
> [!example]- Diagnostic de régression pour une courbe d'étalonnage
> Calibration d'un spectromètre UV-visible pour doser le Fe(II) par la méthode de la phénanthroline.
>
> **Observation des résidus :**
>
> Résidus vs [Fe(II)] : les résidus montrent une structure en forme de parabole, suggérant une déviation de la loi de Beer-Lambert aux fortes concentrations (phénomène d'agrégation ou d'écrantage).
>
> **Action corrective :**
> - Limiter le domaine de linéarité : [Fe(II)] < 5 ppm
> - Ou utiliser un modèle polynomial : $A = \beta_0 + \beta_1 c + \beta_2 c^2$
>
> QQ-plot : Les résidus suivent approximativement une droite, validant l'hypothèse de normalité.
---
### Python
#### Régression simple
**Tracer une régression simple :**
```python
import seaborn as sns
import matplotlib.pyplot as plt
# Régression avec intervalle de confiance à 95%
sns.regplot(x='concentration', y='absorbance', data=df, ci=95)
plt.xlabel('Concentration (μM)')
plt.ylabel('Absorbance')
plt.show()
Matrice de corrélation :
# Calculer et visualiser les corrélations
corr_matrix = df.corr()
sns.heatmap(corr_matrix, annot=True, cmap='coolwarm',
center=0, vmin=-1, vmax=1, square=True)
plt.title('Matrice de Corrélation')
plt.show()
# Extraire une corrélation spécifique
r = corr_matrix.loc['lambda_abs', 'phi_F']
print(f"r = {r:.3f}, R² = {r**2:.3f}")
Régression avec Statsmodels
[!NOTE] Statsmodels
statsmodelsest la bibliothèque Python de référence pour la modélisation statistique et l’inférence.API Formula (syntaxe R, intuitive) :
import statsmodels.formula.api as smf model = smf.ols(formula='Y ~ X1 + X2', data=df).fit()API Standard (matricielle) :
import statsmodels.api as sm X = sm.add_constant(df[['X1', 'X2']]) # Ajoute l'intercept model = sm.OLS(df['Y'], X).fit()
Exemple d’utilisation :
import statsmodels.formula.api as smf
# Ajuster le modèle
model = smf.ols(formula='ln_C ~ temps', data=df).fit()
# Paramètres
print(model.params)
# > Intercept 4.608
# > temps -0.050
# R² et p-values
print(f"R² = {model.rsquared:.4f}")
print(f"p-value (pente) = {model.pvalues['temps']:.2e}")
# Résumé complet
print(model.summary())
Vérification de la formule de l’erreur standard :
n = len(df)
residuals = model.resid
s = residuals.std()
sigma_X = df['temps'].std()
SE_beta1_manual = s / (sigma_X * np.sqrt(n - 1))
SE_beta1_statsmodels = model.bse['temps']
# Les deux valeurs doivent être identiques
Analyse des Résidus
Graphiques diagnostiques essentiels :
import matplotlib.pyplot as plt
from scipy import stats
y_pred = model.predict()
residuals = model.resid
fig, axes = plt.subplots(2, 2, figsize=(12, 10))
# 1. Résidus vs Fitted (homoscédasticité)
axes[0, 0].scatter(y_pred, residuals, alpha=0.6)
axes[0, 0].axhline(y=0, color='r', linestyle='--')
axes[0, 0].set_xlabel('Valeurs ajustées')
axes[0, 0].set_ylabel('Résidus')
axes[0, 0].set_title('Résidus vs Valeurs Ajustées')
# 2. QQ-plot (normalité)
stats.probplot(residuals, dist="norm", plot=axes[0, 1])
axes[0, 1].set_title('QQ-Plot')
# 3. Histogramme
axes[1, 0].hist(residuals, bins=15, edgecolor='black', alpha=0.7)
axes[1, 0].set_xlabel('Résidus')
axes[1, 0].set_title('Distribution des Résidus')
# 4. Résidus standardisés
std_residuals = residuals / residuals.std()
axes[1, 1].scatter(range(len(std_residuals)), std_residuals)
axes[1, 1].axhline(y=0, color='r', linestyle='--')
axes[1, 1].axhline(y=2, color='orange', linestyle=':')
axes[1, 1].axhline(y=-2, color='orange', linestyle=':')
axes[1, 1].set_title('Résidus Standardisés')
plt.tight_layout()
plt.show()
Test de normalité :
from scipy.stats import shapiro
stat, p_value = shapiro(residuals)
print(f"Shapiro-Wilk: p-value = {p_value:.4f}")
# Si p > 0.05 : résidus normaux
Régression Linéaire Multiple
La régression multiple étend le cadre de la régression simple à plusieurs variables explicatives. Elle permet de modéliser l’effet simultané de multiples facteurs et de contrôler les variables confondantes.
Formulation Matricielle
[!NOTE] Définition La régression linéaire multiple modélise $Y$ comme combinaison linéaire de $p$ prédicteurs :
\[Y = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + \cdots + \beta_p X_p + \varepsilon\]Forme matricielle :
\[\mathbf{Y} = \mathbf{X}\boldsymbol{\beta} + \boldsymbol{\varepsilon}\]où :
- $\mathbf{Y}$ : vecteur $(n \times 1)$ des réponses
- $\mathbf{X}$ : matrice $(n \times (p+1))$ des prédicteurs (avec colonne de 1 pour l’intercept)
- $\boldsymbol{\beta}$ : vecteur $((p+1) \times 1)$ des coefficients
- $\boldsymbol{\varepsilon}$ : vecteur $(n \times 1)$ des erreurs
Solution OLS :
\[\hat{\boldsymbol{\beta}} = (\mathbf{X}^T\mathbf{X})^{-1}\mathbf{X}^T\mathbf{Y}\]Cette solution nécessite que $\mathbf{X}^T\mathbf{X}$ soit inversible, donc que $\mathbf{X}$ soit de rang plein (pas de multicolinéarité stricte).
Interprétation des coefficients partiels :
Dans un modèle multiple, $\beta_j$ représente l’effet de $X_j$ sur $Y$ en maintenant tous les autres prédicteurs constants. C’est un coefficient partiel, qui diffère généralement du coefficient simple obtenu en régressant $Y$ sur $X_j$ seul.
On peut se demander si cette interprétation “toutes choses égales par ailleurs” est toujours réalisable expérimentalement. En général, non ! Les variables sont souvent corrélées entre elles, rendant difficile l’isolation de l’effet pur d’un seul facteur. C’est l’essence du problème de la multicolinéarité.
[!example]- Modélisation QSAR d’activité antioxydante On cherche à prédire l’activité antioxydante (IC₅₀, μM) de composés phénoliques en fonction de descripteurs moléculaires :
- $X_1$ : nombre de groupes hydroxyle
- $X_2$ : coefficient de partition octanol-eau (log P)
- $X_3$ : potentiel d’ionisation (eV)
- $X_4$ : énergie de dissociation O-H (kJ/mol)
Modèle ajusté (n = 45 composés) :
\[\log(\text{IC}_{50}) = 12.3 - 0.87 X_1 + 0.31 X_2 - 0.46 X_3 - 0.053 X_4\]
- $R^2 = 0.82$ : le modèle explique 82% de la variance
- $\beta_1 = -0.87$ : chaque OH supplémentaire diminue IC₅₀ (augmente l’activité)
- $\beta_3 = -0.46$ : un PI faible (facilité d’oxydation) favorise l’activité
Ces coefficients partiels permettent de guider la conception de nouveaux antioxydants.
Multicolinéarité
La multicolinéarité survient lorsque deux ou plusieurs prédicteurs sont fortement corrélés entre eux. C’est l’un des problèmes les plus insidieux en régression multiple.
[!NOTE] Définition Multicolinéarité stricte : Relation linéaire exacte entre prédicteurs
\[X_j = \sum_{k \neq j} \alpha_k X_k\]La matrice $\mathbf{X}^T\mathbf{X}$ devient singulière (non inversible), rendant impossible le calcul des coefficients.
Multicolinéarité forte (non stricte) : Corrélation élevée mais non exacte entre prédicteurs.
Conséquences :
- Instabilité des coefficients : Les $\beta_j$ deviennent très sensibles aux variations de l’échantillon
- Erreurs standard gonflées : Les $SE(\beta_j)$ augmentent drastiquement
- p-values non fiables : Tests individuels perdent leur puissance
- Interprétation ambiguë : Impossible de distinguer l’effet de $X_j$ de celui de $X_k$
Ce qui n’est PAS affecté :
- $R^2$ et qualité prédictive globale (tant qu’on prédit dans le même espace)
- Coefficients des variables non colinéaires
- Test F global
Pourquoi est-ce problématique ?
On peut se demander si la multicolinéarité est vraiment un problème si la qualité prédictive reste bonne. La réponse dépend de l’objectif. Pour la prédiction pure, elle n’est pas gênante tant qu’on reste dans le domaine des données. Mais pour l’interprétation et l’inférence causale, elle est désastreuse : les coefficients partiels n’ont plus de sens physique clair.
Imaginons deux catalyseurs photocatalytiques toujours utilisés dans les mêmes proportions. Le modèle ne peut pas distinguer l’effet de l’un de l’effet de l’autre. Les coefficients oscillent de manière erratique entre échantillons, rendant toute conclusion impossible.
Détection : Facteur d’Inflation de la Variance (VIF)
[!NOTE] Définition Le VIF quantifie l’inflation de la variance d’un coefficient due à la multicolinéarité :
\[\text{VIF}(X_j) = \frac{1}{1 - R_j^2}\]où $R_j^2$ est le $R^2$ de la régression de $X_j$ sur tous les autres prédicteurs.
Interprétation :
- $\text{VIF} = 1$ : aucune corrélation avec les autres prédicteurs
- $\text{VIF} > 5$ : multicolinéarité modérée, surveiller
- $\text{VIF} > 10$ : multicolinéarité problématique, action nécessaire
- $\text{VIF} > 100$ : multicolinéarité sévère
Intuition : Si $X_j$ est parfaitement prédictible par les autres ($R_j^2 \to 1$), alors $\text{VIF} \to \infty$ : l’incertitude sur $\beta_j$ explose.
Solutions à la multicolinéarité :
- Retirer une variable redondante : Si deux prédicteurs mesurent essentiellement la même chose
- Combiner les variables : Créer un indice composite (ACP)
- Augmenter la taille de l’échantillon : Réduit l’erreur standard (mais ne résout pas le problème conceptuel)
- Régularisation : Ridge ou LASSO (pénalise les coefficients)
- Centrer les variables : Réduit la corrélation induite par l’intercept
Calcul des VIF en Python :
from statsmodels.stats.outliers_influence import variance_inflation_factor as vif
# Standardiser les données (IMPORTANT pour VIF)
df_scaled = df.copy()
for col in df.columns:
mu = df[col].mean()
sigma = df[col].std()
df_scaled[col] = df_scaled[col].apply(lambda x: (x-mu)/sigma)
# Calculer les VIF
vif_df = pd.DataFrame()
vif_df["features"] = df_scaled.columns
vif_df["vif_index"] = [vif(df_scaled.values, i) for i in range(df_scaled.shape[1])]
# Trier par VIF décroissant
round(vif_df.sort_values(by="vif_index", ascending=False), 2)
Règle empirique : VIF ≥ 10 indique une multicolinéarité problématique.
[!example]- Multicolinéarité dans un modèle QSAR On étudie la solubilité aqueuse de composés organiques en fonction de descripteurs moléculaires :
- $X_1$ : nombre d’atomes de carbone
- $X_2$ : masse molaire
- $X_3$ : log P (coefficient de partition)
- $X_4$ : surface moléculaire accessible
Matrice de corrélation :
$X_1$ $X_2$ $X_3$ $X_4$ $X_1$ 1.00 0.98 0.89 0.95 $X_2$ 0.98 1.00 0.91 0.97 $X_3$ 0.89 0.91 1.00 0.88 $X_4$ 0.95 0.97 0.88 1.00 Calcul des VIF :
Variable VIF $X_1$ 47.2 $X_2$ 83.1 $X_3$ 12.4 $X_4$ 29.6 Diagnostic : Multicolinéarité sévère ! $X_1$, $X_2$, et $X_4$ sont redondants (tous mesurent la “taille” moléculaire).
Solution retenue :
- Éliminer $X_1$ et $X_2$
- Garder $X_3$ (lipophilie, concept distinct) et $X_4$ (taille, VIF acceptable après élimination de $X_1$ et $X_2$)
- Modèle final : $\log(S) = \beta_0 + \beta_3 X_3 + \beta_4 X_4$
- VIF finaux : $\text{VIF}(X_3) = 3.1$, $\text{VIF}(X_4) = 3.1$ → acceptable
R² Ajusté
[!NOTE] Définition Le $R^2$ ajusté pénalise l’ajout de variables non informatives :
\[R^2_{adj} = 1 - \frac{(1-R^2)(n-1)}{n-p-1} = 1 - \frac{\text{SSR}/(n-p-1)}{\text{SST}/(n-1)}\]où :
- $n$ : nombre d’observations
- $p$ : nombre de prédicteurs (hors intercept)
Propriétés :
- $R^2_{adj} \leq R^2$
- $R^2_{adj}$ peut diminuer si on ajoute une variable peu informative
- Plus adapté pour comparer des modèles avec différents nombres de prédicteurs
- Peut être négatif si le modèle est très mauvais
Pourquoi ajuster ?
Le $R^2$ classique augmente toujours (ou reste constant) quand on ajoute une variable, même si celle-ci n’apporte aucune information réelle. C’est un phénomène mécanique : plus de paramètres libres = meilleur ajustement sur les données d’entraînement, au risque de sur-apprentissage (overfitting).
Le $R^2$ ajusté corrige ce biais en intégrant une pénalité proportionnelle au nombre de variables. On peut se demander si cette pénalité est suffisante. Pour la sélection de modèles, des critères plus sophistiqués existent : AIC (Akaike Information Criterion) et BIC (Bayesian Information Criterion), qui pénalisent encore plus fortement la complexité.
[!example]- Comparaison de modèles cinétiques On cherche à modéliser la constante de vitesse d’une réaction photocatalytique en fonction de plusieurs facteurs expérimentaux.
Modèle Prédicteurs $R^2$ $R^2_{adj}$ AIC M1 [Cat] 0.68 0.67 -45.2 M2 [Cat], pH 0.81 0.80 -58.7 M3 [Cat], pH, T 0.83 0.81 -57.1 M4 [Cat], pH, T, [O₂] 0.84 0.81 -55.4 Analyse :
- M1 → M2 : amélioration substantielle ($\Delta R^2_{adj} = +0.13$)
- M2 → M3 : amélioration marginale ($\Delta R^2_{adj} = +0.01$)
- M3 → M4 : $R^2_{adj}$ stagne, AIC se dégrade (overfitting)
Conclusion : Le modèle M2 offre le meilleur compromis entre qualité d’ajustement et parcimonie (principe du rasoir d’Occam).
Test F Global
[!NOTE] Définition Le test F global teste si le modèle dans son ensemble apporte une amélioration significative par rapport au modèle nul (prédiction par la moyenne) :
Hypothèses : $H_0: \beta_1 = \beta_2 = \cdots = \beta_p = 0$ vs $H_1:$ au moins un $\beta_j \neq 0$
Statistique F :
\[F = \frac{\text{MSE}_{\text{model}}}{\text{MSE}_{\text{residual}}} = \frac{\text{SSE}/p}{\text{SSR}/(n-p-1)} = \frac{R^2/p}{(1-R^2)/(n-p-1)}\]sous $H_0$ : $F \sim F_{p, n-p-1}$
où :
- $\text{MSE}_{\text{model}}$ : variance expliquée par variable
- $\text{MSE}_{\text{residual}}$ : variance résiduelle
Décision :
- Si $F > F_{crit}$ ou $p < \alpha$ : rejet de $H_0$, le modèle est globalement significatif
- Attention : un F significatif ne garantit pas que tous les coefficients individuels le sont !
[!example]- Test F pour un modèle de photodégradation Modèle de prédiction du temps de demi-vie $t_{1/2}$ (min) d’un polluant organique en fonction de trois facteurs : [TiO₂], intensité lumineuse $I$, et pH.
Résultats ANOVA :
Source df SS MS F p-value Régression 3 1245.7 415.2 87.4 <0.0001 Résidus 36 171.2 4.8 Total 39 1416.9 Calcul :
- $F = \frac{415.2}{4.8} = 87.4$
- $F_{crit}(3, 36, 0.05) = 2.87$
- $F » F_{crit}$, donc $p < 0.0001$
Conclusion : Le modèle à trois prédicteurs explique de manière hautement significative la variance de $t_{1/2}$. Au moins un (et probablement tous) des facteurs a un effet réel.
Python
Modèle à plusieurs prédicteurs :
import statsmodels.formula.api as smf
# Formule : séparer les prédicteurs par +
model = smf.ols(formula='k_app ~ Cat + pH + T', data=df).fit()
print(model.params)
print(f"R² = {model.rsquared:.4f}")
print(f"R² ajusté = {model.rsquared_adj:.4f}")
Graphiques de régression partielle :
import statsmodels.api as sm
import matplotlib.pyplot as plt
fig = plt.figure(figsize=(15, 5))
sm.graphics.plot_partregress_grid(model, fig=fig)
plt.show()
Chaque graphique montre l’effet d’une variable en maintenant les autres constantes.
Variables catégorielles (encodage automatique) :
# Encoder automatiquement avec C()
model_cat = smf.ols(formula='phi_F ~ C(solvent)', data=df).fit()
print(model_cat.params)
# > Intercept 0.820 (référence : premier solvant alphabétique)
# > C(solvent)[T.toluene] -0.060 (différence vs référence)
# > C(solvent)[T.acetone] -0.280
# > C(solvent)[T.methanol] -0.510
# Sans intercept pour voir toutes les moyennes
model_no_int = smf.ols(formula='phi_F ~ C(solvent) - 1', data=df).fit()
print(model_no_int.params)
# > C(solvent)[hexane] 0.820
# > C(solvent)[toluene] 0.760
# > C(solvent)[acetone] 0.540
# > C(solvent)[methanol] 0.310
Régression Polynomiale
La régression polynomiale permet de modéliser des relations non linéaires tout en conservant le cadre mathématique de la régression linéaire (linéarité par rapport aux paramètres).
[!NOTE] Définition La régression polynomiale de degré $d$ modélise $Y$ par un polynôme de $X$ :
\[Y = \beta_0 + \beta_1 X + \beta_2 X^2 + \cdots + \beta_d X^d + \varepsilon\]Transformation : On crée de nouvelles variables $X_2 = X^2, X_3 = X^3$, etc., et on applique la régression multiple classique.
\[Y = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + \cdots + \beta_d X_d + \varepsilon\]où $X_k = X^k$.
Le modèle reste linéaire en $\boldsymbol{\beta}$, donc toutes les méthodes OLS s’appliquent !
Choix du degré :
- Sous-ajustement (underfitting) : $d$ trop faible, structure systématique dans les résidus
- Sur-ajustement (overfitting) : $d$ trop élevé, bruit modélisé comme signal
- Compromis biais-variance
Attention à la multicolinéarité intrinsèque !
On peut se demander si les puissances successives de $X$ ne vont pas créer de la multicolinéarité. C’est effectivement le cas ! $X, X^2, X^3$ sont naturellement très corrélés. La solution : centrer et réduire les données avant calcul des puissances, ou utiliser des polynômes orthogonaux (Legendre, Chebyshev).
[!example]- Courbe de Stern-Volmer non linéaire En spectroscopie de fluorescence, l’ajout d’un quencher (Q) diminue l’intensité de fluorescence. La relation de Stern-Volmer classique prédit une décroissance linéaire :
\[\frac{I_0}{I} = 1 + K_{SV}[Q]\]Mais expérimentalement, on observe souvent une courbure, indiquant deux populations de fluorophores (accessibles et inaccessibles au quencher).
Modèle polynomial :
\[\frac{I_0}{I} = 1 + K_1[Q] + K_2[Q]^2\]Données expérimentales (quencher : iodure de potassium, fluorophore : tryptophane en protéine) :
[KI] (mM) $I_0/I$ 0 1.00 5 1.82 10 2.95 15 4.38 20 6.15 25 8.28 Régression linéaire : $I_0/I = 1.03 + 0.289[Q]$, $R^2 = 0.981$
Mais résidus structurés (forme parabolique) → sous-ajustement !
Régression quadratique : $I_0/I = 1.01 + 0.201[Q] + 0.0035[Q]^2$, $R^2 = 0.9998$
Résidus aléatoires, ajustement excellent. Le terme quadratique capture l’hétérogénéité des sites de liaison.
Python
Créer et ajuster un modèle polynomial :
import statsmodels.formula.api as smf
# Créer les termes polynomiaux
df['X_squared'] = df['X']**2
df['X_cubed'] = df['X']**3
# Modèle linéaire
model_lin = smf.ols('Y ~ X', data=df).fit()
# Modèle quadratique
model_quad = smf.ols('Y ~ X + X_squared', data=df).fit()
# Modèle cubique
model_cubic = smf.ols('Y ~ X + X_squared + X_cubed', data=df).fit()
# Comparer
print(f"Linéaire: R²={model_lin.rsquared:.4f}, AIC={model_lin.aic:.2f}")
print(f"Quadratique: R²={model_quad.rsquared:.4f}, AIC={model_quad.aic:.2f}")
print(f"Cubique: R²={model_cubic.rsquared:.4f}, AIC={model_cubic.aic:.2f}")
Visualisation rapide :
import numpy as np
X_smooth = np.linspace(df['X'].min(), df['X'].max(), 100)
df_pred = pd.DataFrame({'X': X_smooth, 'X_squared': X_smooth**2})
plt.scatter(df['X'], df['Y'], label='Données')
plt.plot(X_smooth, model_quad.predict(df_pred), 'r-', label='Quadratique')
plt.legend()
plt.show()
Régression Logistique
La régression logistique étend le cadre de la régression aux variables réponses binaires (0/1, succès/échec).
Fonction Logistique
[!NOTE] Définition La régression logistique modélise la probabilité qu’un événement se produise :
\[P(Y=1|X) = \frac{1}{1 + e^{-(\beta_0 + \beta_1 X_1 + \cdots + \beta_p X_p)}} = \frac{e^{\eta}}{1 + e^{\eta}}\]où $\eta = \beta_0 + \beta_1 X_1 + \cdots + \beta_p X_p$ est le prédicteur linéaire.
Transformation logit (log-odds) :
\[\text{logit}(p) = \ln\left(\frac{p}{1-p}\right) = \eta = \beta_0 + \beta_1 X_1 + \cdots + \beta_p X_p\]Le logit transforme une probabilité $p \in [0,1]$ en un nombre réel $\eta \in (-\infty, +\infty)$.
Interprétation des coefficients :
\[\beta_j = \frac{\partial \text{logit}(p)}{\partial X_j}\]$\beta_j > 0$ : $X_j$ augmente la probabilité de $Y=1$
$\beta_j < 0$ : $X_j$ diminue la probabilité de $Y=1$
Odds Ratio (rapport de cotes) :
L’exponentielle des coefficients donne les odds ratios :
\[\text{OR} = e^{\beta_j}\]Si $X_j$ augmente de 1 unité, le ratio $\frac{P(Y=1)}{P(Y=0)}$ est multiplié par $e^{\beta_j}$.
On peut se demander si cette métrique est intuitive. L’odds ratio est plus naturel dans certains contextes (études cas-témoins, paris sportifs) que la probabilité directe. Un OR = 2 signifie que l’odds double, un OR = 0.5 signifie qu’il est divisé par 2.
[!example]- Prédiction du transfert de charge intramoléculaire On étudie 50 chromophores organiques et on cherche à prédire si une transition de type “transfert de charge” (CT) est observée (Y=1) ou non (Y=0) en fonction de :
- $X_1$ : différence de potentiel d’ionisation donneur-accepteur (eV)
- $X_2$ : constante diélectrique du solvant
Modèle logistique ajusté :
\[\ln\left(\frac{p}{1-p}\right) = -4.2 + 1.87 X_1 + 0.35 X_2\]Interprétation :
$\beta_1 = 1.87$ : OR = $e^{1.87} = 6.5$
Une augmentation de 1 eV du gradient électronique multiplie par 6.5 la probabilité relative de transition CT (favorise le transfert d’électron)
$\beta_2 = 0.35$ : OR = $e^{0.35} = 1.42$
Un solvant plus polaire (+1 unité de $\varepsilon$) augmente de 42% la probabilité de CT (stabilise l’état à séparation de charge)
Prédiction pour un nouveau composé :
$X_1 = 3.0$ eV, $X_2 = 10$ (acétone)
\[\eta = -4.2 + 1.87(3.0) + 0.35(10) = 4.91\] \[P(\text{CT}) = \frac{e^{4.91}}{1 + e^{4.91}} = \frac{135.6}{136.6} = 0.993\]Probabilité de 99.3% d’observer une transition CT → fortement probable.
Variables Catégorielles
Les variables catégorielles ne peuvent pas être directement insérées dans un modèle de régression. On doit les encoder numériquement.
Encodage One-Hot (Dummy Variables)
[!NOTE] Définition L’encodage one-hot transforme une variable catégorielle à $k$ modalités en $(k-1)$ variables binaires (dummy variables).
Principe : Pour une variable à $k$ niveaux, créer $(k-1)$ indicatrices :
\[D_j = \begin{cases} 1 & \text{si observation dans catégorie } j \\ 0 & \text{sinon} \end{cases}\]La catégorie omise (niveau de référence) est codée par $(0, 0, …, 0)$.
Interprétation des coefficients :
\[Y = \beta_0 + \beta_1 D_1 + \beta_2 D_2 + \cdots + \beta_{k-1} D_{k-1} + \varepsilon\]
- $\beta_0$ : moyenne du groupe de référence
- $\beta_j$ : différence entre le groupe $j$ et le groupe de référence
Pourquoi $k-1$ et non $k$ variables ?
Si on créait $k$ indicatrices, on aurait $D_1 + D_2 + \cdots + D_k = 1$ (identité), créant une multicolinéarité stricte avec l’intercept. La matrice $\mathbf{X}^T\mathbf{X}$ serait singulière. On sacrifie donc une catégorie comme référence.
On peut se demander si le choix de la référence affecte les résultats. Non pour les prédictions et le $R^2$, mais oui pour l’interprétation des $\beta_j$ ! Choisir la catégorie la plus courante ou le “contrôle” comme référence facilite l’interprétation.
[!example]- Effet du solvant sur le rendement quantique On mesure le rendement quantique de fluorescence $\Phi_F$ d’un chromophore dans 4 solvants : hexane (apolaire), toluène (aromatique apolaire), acétone (polaire aprotique), méthanol (polaire protique).
Solvant $\Phi_F$ moyen n Hexane (réf.) 0.82 8 Toluène 0.76 8 Acétone 0.54 8 Méthanol 0.31 8 Encodage :
- $D_1 = 1$ si toluène, 0 sinon
- $D_2 = 1$ si acétone, 0 sinon
- $D_3 = 1$ si méthanol, 0 sinon
Modèle :
\[\Phi_F = 0.82 - 0.06 D_1 - 0.28 D_2 - 0.51 D_3\]
- $\beta_0 = 0.82$ : rendement dans l’hexane (référence)
- $\beta_1 = -0.06$ : le toluène diminue légèrement $\Phi_F$ (quenching par interactions π-π)
- $\beta_2 = -0.28$ : l’acétone diminue fortement $\Phi_F$ (polarité favorise les états CT non radiatifs)
- $\beta_3 = -0.51$ : le méthanol diminue très fortement $\Phi_F$ (liaisons H + polarité)
Test ANOVA associé : $F(3, 28) = 127.4$, $p < 0.0001$
Le type de solvant a un effet hautement significatif sur le rendement quantique.
Interactions entre Variables
[!NOTE] Définition Une interaction survient lorsque l’effet d’une variable sur $Y$ dépend de la valeur d’une autre variable.
Modèle avec interaction :
\[Y = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + \beta_3 (X_1 \times X_2) + \varepsilon\]Interprétation :
- $\beta_1$ : effet de $X_1$ quand $X_2 = 0$ (effet principal)
- $\beta_2$ : effet de $X_2$ quand $X_1 = 0$ (effet principal)
- $\beta_3$ : modification de l’effet de $X_1$ par unité de $X_2$ (effet d’interaction)
L’effet total de $X_1$ dépend de $X_2$ :
\[\frac{\partial Y}{\partial X_1} = \beta_1 + \beta_3 X_2\]
Quand tester les interactions ?
On peut se demander si on doit systématiquement inclure tous les termes d’interaction possibles. La réponse est non ! Le nombre de termes explose rapidement (avec $p$ prédicteurs, il y a $\binom{p}{2} = \frac{p(p-1)}{2}$ interactions d’ordre 2). Il faut :
- Justification théorique : L’interaction a-t-elle un sens physico-chimique ?
- Exploration visuelle : Les effets semblent-ils synergiques ou antagonistes ?
- Tests statistiques : L’ajout de l’interaction améliore-t-il significativement le modèle ?
[!example]- Interaction pH × [Catalyseur] en photocatalyse On étudie la vitesse de dégradation d’un polluant en fonction du pH et de la concentration en TiO₂.
Modèle sans interaction :
\[k = \beta_0 + \beta_1[\text{TiO}_2] + \beta_2 \text{pH}\]Résultats : $R^2 = 0.71$, résidus structurés.
Modèle avec interaction :
\[k = \beta_0 + \beta_1[\text{TiO}_2] + \beta_2 \text{pH} + \beta_3([\text{TiO}_2] \times \text{pH})\]Résultats : $R^2 = 0.89$, $\beta_3 = -0.034$ (significatif, $p = 0.002$)
Interprétation physique :
L’effet du catalyseur dépend du pH ! À pH élevé, la surface du TiO₂ est chargée négativement (point de charge nulle ~6.5), réduisant l’adsorption du polluant. L’efficacité du catalyseur décroît donc avec le pH.
- À pH 5 : $\frac{\partial k}{\partial [\text{TiO}_2]} = \beta_1 + 5\beta_3 = 0.24 - 0.17 = 0.07$ (effet fort)
- À pH 9 : $\frac{\partial k}{\partial [\text{TiO}_2]} = \beta_1 + 9\beta_3 = 0.24 - 0.31 = -0.07$ (effet négatif !)
Cette interaction explique les résultats contradictoires de la littérature sur l’effet du pH.
Estimation par Maximum de Vraisemblance
[!NOTE] Définition Contrairement à la régression linéaire, la régression logistique ne possède pas de solution analytique fermée. On utilise le Maximum de Vraisemblance (MLE).
Fonction de vraisemblance :
\[L(\boldsymbol{\beta}) = \prod_{i=1}^{n} p_i^{y_i}(1-p_i)^{1-y_i}\]
où $p_i = P(Y_i = 1 X_i) = \frac{1}{1 + e^{-\mathbf{X}_i\boldsymbol{\beta}}}$ Log-vraisemblance :
\[\ell(\boldsymbol{\beta}) = \sum_{i=1}^{n} \left[y_i \ln(p_i) + (1-y_i)\ln(1-p_i)\right]\]Optimisation : On maximise $\ell$ par des méthodes numériques (Newton-Raphson, gradient descent, IRLS).
Les propriétés asymptotiques du MLE s’appliquent :
- Consistance
- Normalité asymptotique : $\hat{\boldsymbol{\beta}} \sim \mathcal{N}(\boldsymbol{\beta}, I^{-1})$
- Efficacité
Tests de Wald :
Pour tester $H_0: \beta_j = 0$ :
\[z = \frac{\hat{\beta}_j}{SE(\hat{\beta}_j)} \sim \mathcal{N}(0,1)\]sous $H_0$ (approximation asymptotique).
Mesures de Qualité d’Ajustement
[!NOTE] Définition Le $R^2$ classique n’est pas défini pour la régression logistique. On utilise des pseudo-R² :
McFadden’s R² :
\[R^2_{McF} = 1 - \frac{\ell(\hat{\boldsymbol{\beta}})}{\ell(\beta_0)}\]où $\ell(\beta_0)$ est la log-vraisemblance du modèle nul (intercept seul).
Interprétation : proportion de “déviance” expliquée, analogue au $R^2$ linéaire.
Valeurs typiques : 0.2-0.4 considéré comme bon ajustement.
Critères d’information :
- AIC (Akaike) : $AIC = -2\ell + 2p$
- BIC (Bayesian) : $BIC = -2\ell + p\ln(n)$
Plus ces critères sont faibles, meilleur est le modèle.
Matrice de confusion et métriques de classification :
Pour évaluer les performances prédictives, on utilise une matrice de confusion (seuil par défaut : 0.5) :
| Prédit 0 | Prédit 1 | |
|---|---|---|
| Vrai 0 | TN | FP |
| Vrai 1 | FN | TP |
- Accuracy : $\frac{TP + TN}{n}$
- Précision : $\frac{TP}{TP + FP}$
- Rappel (Sensibilité) : $\frac{TP}{TP + FN}$
- Spécificité : $\frac{TN}{TN + FP}$
- F1-score : $\frac{2 \times \text{Précision} \times \text{Rappel}}{\text{Précision} + \text{Rappel}}$
Courbe ROC et AUC :
La courbe ROC trace Sensibilité vs (1-Spécificité) pour tous les seuils possibles. L’aire sous la courbe (AUC) mesure la capacité discriminante globale :
- AUC = 0.5 : modèle aléatoire
- AUC = 1.0 : discrimination parfaite
- AUC > 0.8 : excellent modèle
[!example]- Classification de spectres IR On développe un modèle logistique pour classifier automatiquement des spectres IR de composés organiques en deux classes : “aromatiques” vs “aliphatiques”.
Variables : intensités relatives à 3000, 1600, et 1500 cm⁻¹.
Résultats (n=120, 60 par classe) :
- McFadden’s R² = 0.42
- AIC = 87.3
Matrice de confusion (seuil 0.5) :
Prédit Aliph. Prédit Arom. Vrai Aliph. 54 6 Vrai Arom. 8 52
- Accuracy = $\frac{54+52}{120} = 0.883$ (88.3%)
- Précision (Arom.) = $\frac{52}{52+6} = 0.897$
- Rappel (Arom.) = $\frac{52}{52+8} = 0.867$
- AUC = 0.93 : excellente capacité discriminante
Le modèle identifie correctement les cycles aromatiques dans 93% des cas en moyenne.
Régression Logistique en Python
Ajuster un modèle logistique :
import statsmodels.formula.api as smf
# Ajuster le modèle (Y binaire : 0 ou 1)
model_logit = smf.logit('CT ~ delta_IP + epsilon', data=df).fit()
print(model_logit.summary())
# Coefficients et Odds Ratios
coeffs = model_logit.params
odds_ratios = np.exp(coeffs)
print(f"β (delta_IP) = {coeffs['delta_IP']:.3f}")
print(f"OR (delta_IP) = {odds_ratios['delta_IP']:.3f}")
Prédictions et probabilités :
# Prédire les probabilités
df['prob_CT'] = model_logit.predict(df)
# Prédire les classes (seuil = 0.5)
df['pred_CT'] = (df['prob_CT'] > 0.5).astype(int)
Matrice de confusion :
from sklearn.metrics import confusion_matrix, classification_report
cm = confusion_matrix(df['CT'], df['pred_CT'])
print(classification_report(df['CT'], df['pred_CT']))
Courbe ROC et AUC :
from sklearn.metrics import roc_curve, roc_auc_score
fpr, tpr, thresholds = roc_curve(df['CT'], df['prob_CT'])
auc = roc_auc_score(df['CT'], df['prob_CT'])
plt.plot(fpr, tpr, label=f'AUC = {auc:.3f}')
plt.plot([0, 1], [0, 1], 'k--', label='Chance')
plt.xlabel('Taux de Faux Positifs')
plt.ylabel('Taux de Vrais Positifs')
plt.legend()
plt.show()
print(f"AUC = {auc:.3f}")
Prédiction pour un nouveau composé :
# Nouveau composé
nouveau = pd.DataFrame({'delta_IP': [3.0], 'epsilon': [10.0]})
# Probabilité de CT
prob = model_logit.predict(nouveau)[0]
print(f"Probabilité de CT = {prob:.3f} ({prob*100:.1f}%)")
Régression Non Linéaire
Pour des modèles intrinsèquement non linéaires (non linéarisables par transformation), on utilise la régression non linéaire.
[!NOTE] Définition Modèle général :
\[Y = f(X, \boldsymbol{\theta}) + \varepsilon\]où $f$ est une fonction non linéaire des paramètres $\boldsymbol{\theta}$.
Exemples de modèles :
- Michaelis-Menten : $v = \frac{V_{max}[S]}{K_M + [S]}$
- Hill : $Y = \frac{Y_{max}X^n}{K^n + X^n}$
- Décroissance bi-exponentielle : $A(t) = A_1 e^{-k_1 t} + A_2 e^{-k_2 t}$
- Équation de Stern-Volmer modifiée : $\frac{I_0}{I} = (1-f)(1 + K[Q]) + f$
Estimation : Minimisation numérique de SSR par algorithmes itératifs (Levenberg-Marquardt, Gauss-Newton).
[!example]- Cinétique bi-exponentielle de fluorescence On mesure la décroissance temporelle de fluorescence d’un système à deux fluorophores ($\tau_1$ court et $\tau_2$ long).
Modèle :
\[I(t) = A_1 e^{-t/\tau_1} + A_2 e^{-t/\tau_2}\]Données : 200 points de 0 à 50 ns.
Ajustement non linéaire :
Paramètres estimés :
- $A_1 = 4850 \pm 120$ (amplitude composante rapide)
- $\tau_1 = 2.3 \pm 0.1$ ns (temps de vie court)
- $A_2 = 1240 \pm 80$ (amplitude composante lente)
- $\tau_2 = 12.7 \pm 0.4$ ns (temps de vie long)
Qualité de l’ajustement :
- $R^2 = 0.9992$
- Résidus aléatoires, distribution normale
- $\chi^2$ réduit = 1.08 (proche de 1, bon ajustement)
Interprétation : La fraction de fluorescence rapide est $f_1 = \frac{A_1}{A_1+A_2} = 0.80$ (80%), suggérant une population majoritaire dans un environnement polaire (quenching rapide).
Validation et étude des modèles
Validation Croisée
[!NOTE] Définition La validation croisée évalue la capacité de généralisation d’un modèle à de nouvelles données.
K-Fold Cross-Validation :
- Diviser les données en $K$ sous-ensembles (folds) de taille égale
- Pour $k = 1, …, K$ :
- Entraîner le modèle sur $K-1$ folds
- Tester sur le fold $k$
- Calculer une métrique d’erreur (MSE, MAE, etc.)
- Moyenne des $K$ erreurs : performance estimée
Leave-One-Out (LOO) : Cas particulier avec $K = n$
- Avantage : utilise quasi toutes les données pour l’entraînement
- Inconvénient : coûteux en calcul ($n$ itérations)
Erreur de prédiction moyenne :
\[\text{CV}_{(K)} = \frac{1}{K}\sum_{k=1}^{K} \text{MSE}_k = \frac{1}{K}\sum_{k=1}^{K} \frac{1}{n_k}\sum_{i \in \text{fold}_k} (y_i - \hat{y}_i)^2\]
Pourquoi valider ?
On peut se demander si un $R^2$ élevé garantit de bonnes prédictions futures. Non ! Un modèle peut parfaitement ajuster les données d’entraînement (overfitting) tout en échouant sur de nouvelles données. La validation croisée simule la généralisation en testant sur des données “non vues” pendant l’entraînement.
[!example]- Validation d’un modèle QSAR Modèle de prédiction de l’activité cytotoxique (IC₅₀) de 80 dérivés de platine en fonction de 5 descripteurs moléculaires.
Régression linéaire multiple :
- $R^2$ (entraînement complet) = 0.87
- 5-Fold CV : $\text{RMSE}_{CV} = 0.42$ log(μM)
- LOO : $\text{RMSE}_{LOO} = 0.39$ log(μM)
Ajout d’un terme polynomial (degré 2) :
- $R^2$ = 0.94 (amélioration !)
- 5-Fold CV : $\text{RMSE}_{CV} = 0.68$ log(μM) (dégradation !)
Diagnostic : Le modèle polynomial sur-ajuste. Bien qu’il explique mieux les données d’entraînement, sa capacité prédictive se dégrade. Le modèle linéaire est préférable.
from sklearn.model_selection import cross_val_score, KFold
from sklearn.linear_model import LinearRegression
model = LinearRegression()
kf = KFold(n_splits=5, shuffle=True, random_state=42)
# Calculer les scores pour chaque fold
cv_scores_r2 = cross_val_score(model, X, y, cv=kf, scoring='r2')
cv_scores_rmse = np.sqrt(-cross_val_score(model, X, y, cv=kf,
scoring='neg_mean_squared_error'))
print(f"R² moyen (CV): {cv_scores_r2.mean():.4f} ± {cv_scores_r2.std():.4f}")
print(f"RMSE moyen (CV): {cv_scores_rmse.mean():.4f} ± {cv_scores_rmse.std():.4f}")
# Comparer avec l'ajustement complet
model.fit(X, y)
print(f"R² (full data): {model.score(X, y):.4f}")
from sklearn.model_selection import LeaveOneOut
loo = LeaveOneOut()
cv_scores_loo = cross_val_score(model, X, y, cv=loo, scoring='r2')
print(f"R² moyen (LOO): {cv_scores_loo.mean():.4f}")
Courbe d’Apprentissage
[!NOTE] Définition La courbe d’apprentissage trace l’erreur du modèle en fonction de la taille de l’échantillon d’entraînement.
Interprétation :
High bias (sous-ajustement) : Les deux courbes convergent vers une erreur élevée → Modèle trop simple, besoin de plus de complexité
High variance (sur-ajustement) : Large écart entre erreur d’entraînement (faible) et erreur de validation (élevée) → Modèle trop complexe, besoin de plus de données ou de régularisation
Bon équilibre : Les courbes convergent vers une erreur faible acceptable
On peut se demander si ajouter plus de données résout toujours le problème. Non ! Si le modèle est intrinsèquement trop simple (high bias), aucune quantité de données ne l’améliorera. À l’inverse, pour un modèle complexe (high variance), plus de données réduisent l’overfitting.
from sklearn.model_selection import learning_curve
train_sizes, train_scores, val_scores = learning_curve(
model, X, y,
train_sizes=np.linspace(0.1, 1.0, 10),
cv=5, scoring='neg_mean_squared_error'
)
# Convertir en RMSE et calculer moyennes
train_rmse = np.sqrt(-train_scores).mean(axis=1)
val_rmse = np.sqrt(-val_scores).mean(axis=1)
# Visualiser
plt.plot(train_sizes, train_rmse, 'o-', label='Erreur entraînement')
plt.plot(train_sizes, val_rmse, 'o-', label='Erreur validation')
plt.xlabel('Taille échantillon')
plt.ylabel('RMSE')
plt.legend()
plt.show()
Régularisation : Ridge et LASSO
La régularisation pénalise la complexité du modèle pour prévenir l’overfitting.
[!NOTE] Définition Ridge Regression (régularisation L2) :
\[\min_{\boldsymbol{\beta}} \left\{\text{SSR} + \lambda \sum_{j=1}^{p} \beta_j^2\right\} = \min_{\boldsymbol{\beta}} \left\{\|\mathbf{Y} - \mathbf{X}\boldsymbol{\beta}\|^2 + \lambda\|\boldsymbol{\beta}\|^2\right\}\]LASSO (Least Absolute Shrinkage and Selection Operator, régularisation L1) :
\[\min_{\boldsymbol{\beta}} \left\{\text{SSR} + \lambda \sum_{j=1}^{p} |\beta_j|\right\} = \min_{\boldsymbol{\beta}} \left\{\|\mathbf{Y} - \mathbf{X}\boldsymbol{\beta}\|^2 + \lambda\|\boldsymbol{\beta}\|_1\right\}\]où $\lambda \geq 0$ est le paramètre de régularisation.
Effets :
- Ridge : Réduit (shrink) tous les coefficients vers 0, mais ne les annule jamais. Garde toutes les variables.
- LASSO : Peut annuler des coefficients (sélection de variables automatique). Produit des modèles parcimonieux.
Choix de $\lambda$ : Par validation croisée (minimiser l’erreur de prédiction).
Géométrie de la régularisation :
On peut se demander si la différence entre L1 et L2 est profonde. Oui ! Géométriquement, Ridge impose une contrainte sphérique $|\boldsymbol{\beta}|^2 \leq t$, tandis que LASSO impose un diamant $|\boldsymbol{\beta}|_1 \leq t$. Les “coins” du diamant L1 favorisent l’annulation de certains coefficients, d’où la sélection de variables.
[!example]- Régularisation d’un modèle spectroscopique Prédiction de la concentration d’un analyte à partir d’un spectre IR complet (1000 longueurs d’onde = 1000 prédicteurs !).
Avec $n = 50$ échantillons et $p = 1000$ variables, la régression OLS classique est impossible ($p » n$, multicolinéarité extrême).
Ridge Regression ($\lambda$ optimisé par CV) :
- $\lambda_{opt} = 2.3$
- Tous les 1000 coefficients sont non nuls (mais très petits)
- $\text{RMSE}_{CV} = 0.18$ mg/L
LASSO ($\lambda$ optimisé par CV) :
- $\lambda_{opt} = 0.08$
- Seulement 47 coefficients non nuls (sélection de 47 longueurs d’onde informatives)
- $\text{RMSE}_{CV} = 0.16$ mg/L (légèrement meilleur)
Avantage du LASSO : Le modèle est interprétable (47 longueurs d’onde clés identifiées) et plus rapide à déployer.
[!example]- Comparaison Ridge vs LASSO avec scikit-learn Régularisation sur un problème haute dimension :
import numpy as np import pandas as pd import matplotlib.pyplot as plt from sklearn.linear_model import Ridge, Lasso, RidgeCV, LassoCV from sklearn.preprocessing import StandardScaler from sklearn.model_selection import train_test_split from sklearn.metrics import mean_squared_error, r2_score # Simuler un problème haute dimension (spectroscopie IR) np.random.seed(42) n_samples = 50 n_features = 100 # 100 longueurs d'onde # Seules 10 features sont vraiment informatives X = np.random.randn(n_samples, n_features) true_coef = np.zeros(n_features) true_coef[:10] = np.random.randn(10) * 5 # Coefficients non nuls y = X @ true_coef + np.random.randn(n_samples) * 0.5 # Standardiser scaler = StandardScaler() X_scaled = scaler.fit_transform(X) # Split train/test X_train, X_test, y_train, y_test = train_test_split( X_scaled, y, test_size=0.3, random_state=42 ) print("="*70) print("RÉGULARISATION : RIDGE VS LASSO") print("="*70) print(f"Nombre d'observations: {n_samples}") print(f"Nombre de features: {n_features}") print(f"Features informatives: 10") print() # Ridge avec validation croisée pour λ optimal alphas_ridge = np.logspace(-3, 3, 50) ridge_cv = RidgeCV(alphas=alphas_ridge, cv=5) ridge_cv.fit(X_train, y_train) print(f"Ridge - λ optimal: {ridge_cv.alpha_:.4f}") # LASSO avec validation croisée lasso_cv = LassoCV(alphas=None, cv=5, random_state=42, max_iter=10000) lasso_cv.fit(X_train, y_train) print(f"LASSO - λ optimal: {lasso_cv.alpha_:.4f}") print() # Prédictions y_pred_ridge = ridge_cv.predict(X_test) y_pred_lasso = lasso_cv.predict(X_test) # Métriques print("="*70) print("PERFORMANCES SUR TEST SET") print("="*70) models_compare = { 'Ridge': (ridge_cv, y_pred_ridge), 'LASSO': (lasso_cv, y_pred_lasso) } for name, (model, y_pred) in models_compare.items(): rmse = np.sqrt(mean_squared_error(y_test, y_pred)) r2 = r2_score(y_test, y_pred) n_nonzero = np.sum(model.coef_ != 0) print(f"\n{name}:") print(f" RMSE: {rmse:.4f}") print(f" R²: {r2:.4f}") print(f" Coefficients non nuls: {n_nonzero}/{n_features}") # Visualisation des coefficients fig, axes = plt.subplots(1, 3, figsize=(18, 5)) # Vrais coefficients axes[0].stem(range(n_features), true_coef, basefmt=' ') axes[0].set_xlabel('Index de la feature', fontsize=11) axes[0].set_ylabel('Coefficient', fontsize=11) axes[0].set_title('Vrais Coefficients', fontsize=12, fontweight='bold') axes[0].axhline(y=0, color='red', linestyle='--', alpha=0.5) axes[0].grid(alpha=0.3) # Ridge axes[1].stem(range(n_features), ridge_cv.coef_, basefmt=' ') axes[1].set_xlabel('Index de la feature', fontsize=11) axes[1].set_ylabel('Coefficient', fontsize=11) axes[1].set_title(f'Ridge (λ={ridge_cv.alpha_:.3f})\n{np.sum(ridge_cv.coef_ != 0)} coef non nuls', fontsize=12, fontweight='bold') axes[1].axhline(y=0, color='red', linestyle='--', alpha=0.5) axes[1].grid(alpha=0.3) # LASSO axes[2].stem(range(n_features), lasso_cv.coef_, basefmt=' ') axes[2].set_xlabel('Index de la feature', fontsize=11) axes[2].set_ylabel('Coefficient', fontsize=11) axes[2].set_title(f'LASSO (λ={lasso_cv.alpha_:.3f})\n{np.sum(lasso_cv.coef_ != 0)} coef non nuls', fontsize=12, fontweight='bold') axes[2].axhline(y=0, color='red', linestyle='--', alpha=0.5) axes[2].grid(alpha=0.3) plt.tight_layout() plt.show() # Analyser la sélection de variables par LASSO print("\n" + "="*70) print("SÉLECTION DE VARIABLES PAR LASSO") print("="*70) selected_features = np.where(lasso_cv.coef_ != 0)[0] print(f"Features sélectionnées: {list(selected_features)}") print(f"Nombre de features sélectionnées: {len(selected_features)}") # Vérifier le recall true_features = np.where(true_coef != 0)[0] correctly_selected = np.intersect1d(selected_features, true_features) recall = len(correctly_selected) / len(true_features) precision = len(correctly_selected) / len(selected_features) if len(selected_features) > 0 else 0 print(f"\nRecall (sensibilité): {recall:.2%}") print(f"Precision: {precision:.2%}")
[!example]- Chemin de régularisation (Regularization path) Visualiser l’évolution des coefficients en fonction de λ :
from sklearn.linear_model import lasso_path, ridge_path # Calculer les chemins de régularisation alphas_path = np.logspace(-3, 1, 100) # LASSO path alphas_lasso, coefs_lasso, _ = lasso_path(X_train, y_train, alphas=alphas_path) # Ridge path (via calcul manuel) coefs_ridge = [] for alpha in alphas_path: ridge_temp = Ridge(alpha=alpha) ridge_temp.fit(X_train, y_train) coefs_ridge.append(ridge_temp.coef_) coefs_ridge = np.array(coefs_ridge).T # Visualisation fig, axes = plt.subplots(1, 2, figsize=(16, 6)) # LASSO path for i in range(10): # Montrer seulement les 10 premières features axes[0].plot(np.log10(alphas_lasso), coefs_lasso[i, :], linewidth=2) axes[0].axvline(x=np.log10(lasso_cv.alpha_), color='red', linestyle='--', linewidth=2, label=f'λ optimal = {lasso_cv.alpha_:.3f}') axes[0].set_xlabel('log₁₀(λ)', fontsize=12) axes[0].set_ylabel('Coefficients', fontsize=12) axes[0].set_title('LASSO: Chemin de Régularisation', fontsize=13, fontweight='bold') axes[0].legend() axes[0].grid(alpha=0.3) # Ridge path for i in range(10): axes[1].plot(np.log10(alphas_path), coefs_ridge[i, :], linewidth=2) axes[1].axvline(x=np.log10(ridge_cv.alpha_), color='red', linestyle='--', linewidth=2, label=f'λ optimal = {ridge_cv.alpha_:.3f}') axes[1].set_xlabel('log₁₀(λ)', fontsize=12) axes[1].set_ylabel('Coefficients', fontsize=12) axes[1].set_title('Ridge: Chemin de Régularisation', fontsize=13, fontweight='bold') axes[1].legend() axes[1].grid(alpha=0.3) plt.tight_layout() plt.show() print("\nObservation :") print("- LASSO: Les coefficients tombent exactement à zéro → sélection de variables") print("- Ridge: Les coefficients se rapprochent de zéro mais ne l'atteignent jamais")
Applications avancées en Chimie : plans d’expériences et surfaces de réponse
[!NOTE] Définition Les plans d’expériences (Design of Experiments, DoE) optimisent le nombre et la disposition des expériences pour modéliser une surface de réponse.
Plans factoriels :
- Complet $2^k$ : Toutes les combinaisons de $k$ facteurs à 2 niveaux
- Fractionnaire : Sous-ensemble du plan complet (économie d’expériences)
- Central composite : Ajoute des points centraux et étoiles pour modéliser la courbure
Modèle de surface de réponse (ordre 2) :
\[Y = \beta_0 + \sum_{i=1}^{k}\beta_i X_i + \sum_{i=1}^{k}\beta_{ii}X_i^2 + \sum_{i<j}\beta_{ij}X_iX_j + \varepsilon\]Permet de capturer effets principaux, quadratiques, et interactions.
[!example]- Optimisation d’une synthèse photochimique On cherche à maximiser le rendement d’une cycloaddition [2+2] photocatalysée en fonction de trois facteurs :
- $X_1$ : [Photocatalyseur] (0.1-1.0 mM)
- $X_2$ : Puissance d’irradiation (10-100 mW/cm²)
- $X_3$ : Température (0-40°C)
Plan central composite : 20 expériences (8 factorielles, 6 axiales, 6 points centraux).
Modèle ajusté :
\[\text{Rendement} = 45.2 + 12.3X_1 + 8.7X_2 - 3.1X_3 - 2.8X_1^2 - 1.9X_2^2 + 0.4X_3^2 + 4.2X_1X_2 - 1.5X_1X_3\]
- $R^2 = 0.94$ : excellent ajustement
- Interaction $X_1 \times X_2$ positive : l’effet du catalyseur est amplifié par une forte irradiation
- Termes quadratiques négatifs pour $X_1$ et $X_2$ : maxima internes
Optimisation numérique :
Conditions optimales prédites : [Cat] = 0.67 mM, $I$ = 78 mW/cm², $T$ = 15°C
Rendement prédit : 68.4%
Rendement expérimental (validation) : 66.9% → excellente prédiction !
Résumé
| Concept | Formule |
|---|---|
| Régression linéaire simple | $Y = \beta_0 + \beta_1 X + \varepsilon$ |
| Pente OLS | $\beta_1 = \frac{\text{Cov}(X,Y)}{\text{Var}(X)}$ |
| Intercept OLS | $\beta_0 = \bar{y} - \beta_1 \bar{x}$ |
| SSR | $\sum_{i=1}^{n}(y_i - \hat{y}_i)^2$ |
| $R^2$ | $1 - \frac{\text{SSR}}{\text{SST}}$ |
| $R^2$ ajusté | $1 - \frac{(1-R^2)(n-1)}{n-p-1}$ |
| Test t pour $\beta_j$ | $t = \frac{\hat{\beta}_j}{SE(\hat{\beta}_j)}$ |
| Test F global | $F = \frac{R^2/p}{(1-R^2)/(n-p-1)}$ |
| VIF | $\text{VIF}(X_j) = \frac{1}{1-R_j^2}$ |
| Régression logistique | $\ln\left(\frac{p}{1-p}\right) = \boldsymbol{\beta}^T\mathbf{X}$ |
| Ridge | $\min {|\mathbf{Y}-\mathbf{X}\boldsymbol{\beta}|^2 + \lambda|\boldsymbol{\beta}|^2}$ |
| LASSO | $\min {|\mathbf{Y}-\mathbf{X}\boldsymbol{\beta}|^2 + \lambda|\boldsymbol{\beta}|_1}$ |
- Exploration des données
- Visualiser les relations (scatterplots, matrices de corrélation)
- Détecter les outliers
- Vérifier la distribution des variables
- Préparation des données
- Traiter les valeurs manquantes
- Encoder les variables catégorielles
- Centrer et réduire si nécessaire
- Transformer les variables si non-linéarité détectée
- Construction du modèle
- Commencer simple (régression simple ou multiple avec peu de variables)
- Ajouter de la complexité si justifié (interactions, polynômes)
- Vérifier la multicolinéarité (VIF)
- Évaluation
- Analyser les résidus (linéarité, homoscédasticité, normalité, indépendance)
- Tester la significativité (tests F, t, p-values)
- Calculer $R^2$, $R^2_{adj}$, AIC, BIC
- Validation
- Validation croisée (K-fold, LOO)
- Test sur ensemble indépendant si possible
- Vérifier la stabilité des coefficients
- Interprétation
- Interpréter les coefficients dans leur contexte physico-chimique
- Distinguer corrélation et causalité
- Identifier les limites du modèle (domaine de validité, hypothèses)
Références
-
Cours Le Wagon ![[Statistics & Probabilities.pdf]]
- Livres de référence :
- Kutner et al. (2004). Applied Linear Statistical Models. McGraw-Hill.
- James et al. (2021). An Introduction to Statistical Learning. Springer.
- Hastie et al. (2009). The Elements of Statistical Learning. Springer.
- Ressources en ligne :
- Applications en chimie :
- Consonni et al. (2019). “Comments on the definition of the Q² parameter for QSAR validation”. J. Chem. Inf. Model. 59(5), 1669-1673.
- Gramatica (2007). “Principles of QSAR models validation”. QSAR Comb. Sci. 26(5), 694-701.