R-formules
ECUE 1 – Modélisation linéaire

Formules R du cours

Toutes les commandes R utilisées dans le cours de régression linéaire simple, régression multiple et analyse de la variance (ANOVA), regroupées par thème.

Étape 1

Chargement et affichage des données

Importer les données, les afficher et les explorer avant toute modélisation.

Import Lire un fichier CSV (séparateur ;)

Importe un fichier CSV dont les colonnes sont séparées par des point-virgules (format français).

# Changer le répertoire si nécessaire : Fichier > Changer le répertoire courant
mesdonnees <- read.csv2("exemple_1.csv")
Import Saisir des données au clavier

Pour de petits jeux de données, on peut créer directement les vecteurs x et y.

surface_xi <- c(28, 50, 55, 60, 105, 106, 117, 190, 196)
prix_yi    <- c(132, 280, 268, 320, 375, 570, 670, 790, 800)
mesdonnees <- data.frame(surface_xi = surface_xi, prix_yi = prix_yi)
Exploration Afficher et résumer les données

Afficher le tableau complet, les premières/dernières lignes, et un résumé statistique.

mesdonnees              # affiche tout le tableau
head(mesdonnees)         # 6 premières lignes
tail(mesdonnees, 4)     # 4 dernières lignes
names(mesdonnees)        # noms des colonnes
summary(mesdonnees)      # résumé statistique
Visualisation Nuage de points

Représenter le nuage de points (xi, yi) pour visualiser la liaison entre les variables.

plot(mesdonnees)
Chapitre I

Régression linéaire simple

Estimer les coefficients, tester leur significativité, construire des intervalles de confiance et faire des prévisions.

Estimation Ajuster un modèle de régression simple

Calcule les estimateurs des moindres carrés â et b̂ du modèle Y = ax + b + ε.

reg <- lm(prix_yi ~ surface_xi, data = mesdonnees)
Formule : variable_Y ~ variable_X
Résultats Résumé complet du modèle (summary)

Affiche â, b̂, leurs écarts-types, les statistiques de Student, les p-valeurs et R².

summary(reg)
Coefficients Extraire les coefficients individuellement

Obtenir â et b̂ séparément depuis l'objet de régression.

coef(reg)     # tous les coefficients : b̂ et â
coef(reg)[1]  # b̂ (Intercept)
coef(reg)[2]  # â (pente)
ANOVA Tableau d'analyse de la variance

Affiche les sommes des carrés SCreg, SCres et les degrés de liberté.
Permet d'obtenir SCres = Σei².

aov(reg)
# Pour plus de chiffres significatifs :
options("digits" = 12)
aov(reg)
ANOVA Test de Fisher (anova)

Tableau ANOVA complet avec la statistique F pour tester H₀ : a = 0.
Équivalent au test de Student au carré.

anova(reg)
Tests Quantile de Student (valeur critique)

Obtenir la valeur critique t_{1-α/2}(n-2) pour un test de Student à n-2 degrés de liberté.

# Pour α = 0.05 et n = 9 (ddl = n-2 = 7) :
qt(0.975, 7)
# Résultat : 2.365 (valeur critique bilatérale)
Tests Quantile de Fisher (valeur critique)

Obtenir la valeur critique F(α, k, n-k-1) pour un test de Fisher.

# Pour α = 0.05, k = 1 variable, n-k-1 = 7 ddl :
qf(0.95, 1, 7)
# Résultat : 5.59 (valeur critique)
Intervalles Intervalles de confiance pour â et b̂

Calcule les IC à 95 % pour les paramètres a (pente) et b (ordonnée à l'origine).

confint(reg, level = 0.95)
#                 2.5 %     97.5 %
# (Intercept) -52.39     200.98
# surface_xi    2.80       4.99
Prévision Intervalles de confiance pour E(Y|X=x)

Pour chaque xi, calcule l'IC pour la valeur moyenne prédite E(Yx) = ax + b. L'intervalle s'élargit quand x s'éloigne de x̄.

# IC pour tous les points du jeu de données :
conf <- predict.lm(reg, interval = "confidence")

# IC pour une valeur particulière x = 120 :
predict(reg, data.frame(surface_xi = c(120)),
        interval = "conf")
#      fit      lwr      upr
# 1 542.17  476.55   607.79
Astuce : conf[,2] donne toutes les bornes inférieures, conf[,3] les bornes supérieures.
Prévision Intervalle de prédiction pour Yx

Plus large que l'IC de confiance : contient une observation future Y avec probabilité 1-α.

# Intervalle de prédiction pour tous les points :
pred <- predict.lm(reg, interval = "prediction")

# Intervalle de prédiction pour x = 120 :
predict.lm(reg, data.frame(surface_xi = c(120)),
           interval = "pred")
#      fit      lwr      upr
# 1 542.17  344.52   739.82
Visualisation Tracer les intervalles de confiance/prédiction

Superposer les courbes des intervalles de confiance (bleu) et de prédiction (vert) au nuage de points.

plot(mesdonnees)
# Intervalles de confiance (bleu, pointillés) :
lines(mesdonnees$surface_xi, conf[,2], lwd=2, col="blue", lty=2)
lines(mesdonnees$surface_xi, conf[,3], lwd=2, col="blue", lty=2)
# Intervalles de prédiction (vert, tirets) :
lines(mesdonnees$surface_xi, pred[,2], lwd=2, col="green", lty=6)
lines(mesdonnees$surface_xi, pred[,3], lwd=2, col="green", lty=6)
Chapitre II

Régression linéaire multiple

Extension au cas de k variables explicatives : y = Xβ + ε. Estimation, tests et intervalles de confiance.

Estimation Ajuster un modèle avec plusieurs variables

Estime β̂ = (X'X)⁻¹X'y pour le modèle y = β₀ + β₁x₁ + β₂x₂ + ε.

mesdonnees <- read.csv2("exemple_2.csv")
reg <- lm(Production ~ Travail + Capital, data = mesdonnees)
summary(reg)
Syntaxe générale : lm(Y ~ X1 + X2 + … + Xk, data = …)
Intervalles Intervalles de confiance (régression multiple)

IC à 95 % pour chaque βⱼ. Utilise la loi de Student avec n-k-1 degrés de liberté.

confint(reg, level = 0.95)
Variance Matrice de variance-covariance des estimateurs

Retourne σ̂² (X'X)⁻¹. La diagonale donne Var(β̂ⱼ) ; les autres termes donnent Cov(β̂ⱼ, β̂ₗ).

vcov(reg)
Rappel : ŝⱼ = √(terme diagonal jj) est l'écart-type estimé de β̂ⱼ.
ANOVA Tableau ANOVA (régression multiple)

Décompose SCtot = SCreg + SCres avec les degrés de liberté k et n-k-1. Teste H₀ : β₁ = … = βk = 0.

anova(reg)
# Puis le résumé global :
summary(reg)
Tests Quantile de Fisher (régression multiple)

Valeur critique F(α, k, n-k-1) pour le test de Fisher global en régression multiple.

# Pour α = 0.05, k = 2 var., n-k-1 = 6 ddl :
qf(0.95, 2, 6)
# Valeur critique pour le test de Student :
# Pour n-k-1 = 6 ddl :
qt(0.975, 6)   # donne t_{0.025}(6) = 2.447
Chapitre II – C

Sélection des variables

Choisir les variables explicatives les plus pertinentes par méthodes ascendante ou pas-à-pas.

Exploration Diagramme de dispersion de toutes les paires

Visualiser simultanément toutes les relations entre variables avant de construire le modèle.

mesdonnees <- read.csv2("exemple_2.csv")
pairs(mesdonnees)
Méthode ascendante Ajout pas à pas (add1)

Teste quelles variables à ajouter au modèle courant en minimisant l'AIC et en regardant Pr(>F). Ajouter celle avec le Pr(>F) le plus faible (si < α).

mesdonnees <- read.csv2("ozone.txt")
# Étape 1 : modèle vide, chercher 1ère variable
add1(lm(maxO3 ~ 1, data=mesdonnees),
     maxO3 ~ T9+T12+T15+Ne9+Ne12+Ne15+maxO3v,
     test="F")
# Étape 2 : après ajout de T12
add1(lm(maxO3 ~ T12, data=mesdonnees),
     maxO3 ~ T9+T12+T15+Ne9+Ne12+Ne15+maxO3v,
     test="F")
# Continuer jusqu'à ce qu'aucun Pr(>F) < seuil
AIC : AIC = n·[log(2πσ̂) + log(SCres/n)] + 2k – 2. Minimiser l'AIC.
Stepwise Méthode pas-à-pas bidirectionnelle (step)

Algorithme automatique : ajoute et retire des variables à chaque étape jusqu'à convergence. Minimise l'AIC.

step(lm(maxO3 ~ 1, data=mesdonnees),
     maxO3 ~ T9+T12+T15+Ne9+Ne12+Ne15+maxO3v,
     direction="both")
direction : "both" (bidirectionnel), "forward" (ascendant), "backward" (descendant)
Chapitre II – C

Diagnostics du modèle

Vérifier les hypothèses du modèle : homoscédasticité, normalité des résidus, absence de valeurs aberrantes.

Résidus Graphique des résidus vs valeurs estimées

Si le modèle est adéquat (homoscédasticité), les résidus eᵢ = yᵢ - ŷᵢ doivent être répartis uniformément dans une bande horizontale.

plot(reg$residuals)
# Ou avec les résidus studentisés (bornes ±2) :
plot(rstudent(reg))
Normalité QQ-Plot des résidus

Vérifie la normalité des erreurs : les points doivent être alignés sur la droite gaussienne. Ajouter qqline pour la droite de référence.

qqnorm(reg$residuals)
qqline(reg$residuals)
# Pour afficher les labels d'identification :
coords <- qqnorm(reg$residuals)
qqline(reg$residuals)
text(coords$x, coords$y,
     labels=mesdonnees$obs, pos=4, cex=0.7)
Visualisation Droites de régression par groupe (ggplot2)

Représenter les droites de régression par modalité d'une variable qualitative (ex. type de cidre).

library(ggplot2)
# install.packages("ggplot2") si nécessaire

ggplot(mesdonnees,
       aes(y=S.Sucree, x=S.Amere, col=Type)) +
  geom_point() +
  geom_smooth(method="lm", se=FALSE)
# se=TRUE pour afficher les IC autour des droites
Chapitre III

Analyse de la variance (ANOVA)

Tester l'égalité des moyennes entre groupes définis par une ou plusieurs variables qualitatives.

ANOVA 1 facteur Ajuster un modèle ANOVA à un facteur

Teste H₀ : μ₁ = μ₂ = … = μI (égalité des moyennes par groupe). Variable réponse quantitative, facteur qualitatif.

mesdonnees <- read.csv2("notes.csv")
reg <- lm(Note ~ Examinateur, data=mesdonnees)
anova(reg)
# Valeur critique F(α, I-1, n-I) :
qf(0.95, 2, 18)   # ici I=3 examinateurs, n=21
ANOVA 2 facteurs ANOVA à deux facteurs avec interaction

Teste l'effet de chaque facteur et leur interaction. L'opérateur * inclut les deux effets principaux et l'interaction.

mesdonnees <- read.csv2("ble.txt")
# Modèle avec interaction (A*B = A + B + A:B) :
reg <- lm(rendement ~ phyto * variete, data=mesdonnees)
anova(reg)
# Modèle sans interaction :
reg2 <- lm(rendement ~ phyto + variete, data=mesdonnees)
anova(reg2)
Opérateurs : * pour facteur A × facteur B (effets + interaction) ; + pour effets additifs sans interaction.
ANCOVA Analyse de covariance (ANCOVA)

Modèle mixte : variable quantitative Y expliquée par une variable quantitative X et un facteur qualitatif. Permet de tester l'égalité des pentes et des ordonnées à l'origine.

mesdonnees <- read.csv2("cidre.csv")
# Modèle avec pentes différentes (interaction) :
reg <- lm(S.Sucree ~ S.Amere * Type, data=mesdonnees)
anova(reg)
# Modèle avec pentes égales (sans interaction) :
reg2 <- lm(S.Sucree ~ S.Amere + Type, data=mesdonnees)
anova(reg2)
# Modèle avec seul effet de S.Amere :
reg3 <- lm(S.Sucree ~ S.Amere, data=mesdonnees)
anova(reg3)
Résumé Résumé global du modèle ANOVA

Affiche les coefficients estimés (effets de chaque modalité), les statistiques de Student associées, et R².

summary(reg)
Lecture Comment lire summary(reg) — régression simple
Exemple de sortie R
Call:
lm(formula = prix_yi ~ surface_xi, data = mesdonnees)

Residuals:
    Min      1Q  Median      3Q     Max
-119.07  -46.68   15.90   37.16  109.34

Coefficients:
             Estimate Std. Error t value Pr(>|t|)
(Intercept)  74.2961   53.5743    1.387   0.2074
surface_xi    3.8989    0.4632    8.417   5.36e-05 ***

Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1

Residual standard error: 78.85 on 7 degrees of freedom
Multiple R-squared:  0.9102,	Adjusted R-squared:  0.8974
F-statistic: 70.84 on 1 and 7 DF,  p-value: 5.358e-05
Estimate : valeurs de b̂ (Intercept) et â (pente)
t value : statistique de Student = Estimate / Std. Error
Pr(>|t|) : p-valeur du test H₀ : coefficient = 0
Codes de significativité : indiquent à quel point on peut rejeter H₀ : coefficient = 0
***  p < 0.001 — résultat très hautement significatif : risque de se tromper < 0.1 %. On rejette H₀ avec une très grande certitude.
**   p < 0.01  — résultat hautement significatif : risque < 1 %. On rejette H₀ avec une grande certitude.
*    p < 0.05  — résultat significatif : risque < 5 %. Seuil classique en statistique — on rejette H₀.
.    p < 0.10  — résultat marginalement significatif : tendance, mais insuffisant au seuil 5 %.
     p ≥ 0.10  — résultat non significatif : on ne peut pas rejeter H₀, le coefficient pourrait être nul.
R-squared : R² = part de variance expliquée (ici 91 %)
F-statistic : test global H₀ : a = 0 (ici F = 70.84 > F critique)
Lecture Comment lire anova(reg)
Exemple de sortie R
Analysis of Variance Table

Response: prix_yi
           Df  Sum Sq Mean Sq F value    Pr(>F)
surface_xi  1 440381  440381  70.839 5.358e-05 ***
Residuals   7  43517    6217
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05
Sum Sq : SCreg = 440 381 (variance expliquée) et SCres = 43 517 (résiduelle)
Mean Sq : MCreg = SCreg/1 et MCres = SCres/7 = σ̂² = 6 217
F value : F = MCreg / MCres — à comparer à qf(0.95, 1, 7) = 5.59
Pr(>F) : p-valeur du test de Fisher — si < α on rejette H₀ : a = 0
*** : p < 0.001 — la variable surface_xi est très hautement significative dans le modèle. On rejette H₀ : a = 0 avec un risque < 0.1 %.
Lecture Comment lire confint(reg, level=0.95)
Exemple de sortie R
                  2.5 %      97.5 %
(Intercept) -52.386951  200.979174
surface_xi    2.803539    4.994333
2.5 % / 97.5 % : bornes inférieure et supérieure de l'IC à 95 %
(Intercept) : IC pour b — contient 0, donc b n'est pas significativement différent de 0
surface_xi : IC pour a = [2.80 ; 4.99] — ne contient pas 0, donc a est significatif
Lecture Comment lire summary(reg) — régression multiple
Exemple de sortie R (Production ~ Travail + Capital)
Call:
lm(formula = Production ~ Travail + Capital, data = mesdonnees)

Coefficients:
             Estimate Std. Error t value Pr(>|t|)
(Intercept) -437.710   57.933   -7.556  0.000283 ***
Travail        0.336    0.090    3.748  0.009627  **
Capital        0.410    0.196    2.092  0.081584    .

Residual standard error: 23.07 on 6 degrees of freedom
Multiple R-squared:  0.9784,	Adjusted R-squared:  0.9712
F-statistic: 135.9 on 2 and 6 DF,  p-value: 4.384e-06
Estimate : β̂₀ = −437.71, β̂₁ = 0.336, β̂₂ = 0.41
Pr(>|t|) Capital = 0.082 : > 0.05 donc Capital non significatif au seuil 5 % → envisager de le retirer
Rappel des codes de significativité :
*** p < 0.001 — très hautement significatif (β₀ ici)
**  p < 0.01  — hautement significatif (Travail ici)
.   p < 0.10  — marginalement significatif (Capital ici) : tendance à la limite du seuil 5 %
    p ≥ 0.10  — non significatif : envisager de retirer la variable du modèle
R² ajusté : pénalise les variables superflues — à préférer au R² brut pour comparer des modèles
F global = 135.9 : au moins une variable est significative dans le modèle
Lecture Comment lire anova(reg) — régression multiple
Exemple de sortie R
Analysis of Variance Table

Response: Production
          Df  Sum Sq Mean Sq  F value    Pr(>F)
Travail    1 142369  142369  267.369   1.29e-06 ***
Capital    1   2326    2326    4.376    0.08158    .
Residuals  6   3194     532
---
SCtot = 142369 + 2326 + 3194 = 147889
Sum Sq : contribution de chaque variable à la variance expliquée (séquentielle)
Residuals Sum Sq = 3 194 : SCres — variance non expliquée par le modèle complet
F value : test d'apport de chaque variable ajoutée séquentiellement
Pr(>F) Travail *** : p < 0.001 — Travail est très hautement significatif, son apport au modèle est indiscutable.
Pr(>F) Capital = 0.082 (.) : p < 0.10 mais > 0.05 — Capital est marginalement significatif ; son apport est discutable une fois Travail inclus. On peut envisager un modèle sans Capital.