Modèle linéaire gaussien

Didier Chauveau et Jean-François Delmas
(last modification date: October 10, 2017)
Version pdf de ce document
Version sans bandeaux

1 Peut-on augmenter la quantité de pluie?
2 (Facultatif) Analyse de la variance à un facteur
3 Régression linéaire

Contents

1 Peut-on augmenter la quantité de pluie?
 1.1 La problématique
 1.2 Les données
 1.3 Le modèle
 1.4 Comparaison des moyennes
 1.5 (Facultatif) Comparaison des variances des 2 échantillons gaussiens
2 (Facultatif) Analyse de la variance à un facteur
 2.1 Les données
 2.2 Le modèle
3 Régression linéaire
 3.1 Les données
 3.2 Régression linéaire simple
 3.3 Régression linéaire multiple

1 Peut-on augmenter la quantité de pluie?

1.1 La problématique

Il existe deux types de nuages qui donnent lieu à des précipitations: les nuages chauds et les nuages froids. Ces derniers possèdent une température maximale de l’ordre de -10°C à -25°C. Ils sont composés de cristaux de glace et de gouttelettes d’eau. Ces gouttelettes d’eau subsistent alors que la température ambiante est inférieure à la température de fusion. On parle d’eau surfondue. Leur état est instable. De fait, quand une particule de glace rencontre une gouttelette d’eau, elles s’aggrègent pour ne former qu’une seule particule de glace. Les particules de glace, plus lourdes que les gouttelettes, tombent sous l’action de la gravité. Enfin si les températures des couches d’air inférieures sont suffisamment élevées, les particules de glace fondent au cours de leur chute formant ainsi de la pluie.

En l’absence d’un nombre suffisant de cristaux de glace pour initier le phénomène décrit ci-dessus, on peut injecter dans le nuage froid des particules qui ont une structure cristalline proche de la glace, par exemple de l’iodure d’argent (environ 100 à 1000 grammes par nuage). Autour de ces particules, on observe la formation de cristaux de glace, ce qui permet, on l’espère, de déclencher ou d’augmenter les précipitations. Il s’agit de l’ensemencement des nuages. Signalons que cette méthode est également utilisée pour limiter le risque de grêle.

Il est évident que la possibilité de modifier ainsi les précipitations présente un grand intérêt pour l’agriculture. De nombreuses études ont été et sont encore consacrées à l’étude de l’efficacité de l’ensemencement des nuages dans divers pays. L’étude de cette efficacité est cruciale et délicate. Le débat est encore d’actualité.

1.2 Les données

On dispose des données concernant l’ensemencement par iodure d’argent des nuages en Floride (1975)1 .

Les volumes de pluie déversées en 107m3 (cf. les deux tableaux ci-dessous) concernent 23 jours similaires dont n1 = 11 jours avec ensemencement correspondant aux réalisations des variables aléatoires X1,,Xn1 et n2 = 12 jours sans ensemencement correspondant aux réalisations des variables aléatoires Y 1,,Y n2.














i 1 2 3 4 5 6 7 8 9 10 11












X i7.454.707.304.054.469.7015.108.518.132.202.16













Table 1: Volume de pluie en 107m3 déversée avec ensemencement















j 1 2 3 4 5 6 7 8 9 10 11 12













Y j15.5310.394.503.445.708.246.736.217.584.171.093.50














Table 2: Volume de pluie en 107m3 déversée sans ensemencement

X=[7.45, 4.70, 7.30, 4.05, 4.46, 9.70, 15.10, 8.51, 8.13, 2.20, 2.16];  
Y=[15.53, 10.39, 4.50, 3.44, 5.70, 8.24, 6.73, 6.21, 7.58, 4.17, 1.09, 3.50];

1.3 Le modèle

On suppose que les volumes de pluies suivent une loi gaussienne et que l’effet de l’ensemencement ne modifie pas la variance des lois gaussiennes.

Ceci revient à dire dans ce contexte paramétrique que les observations sont constituées de deux échantillons gaussiens indépendants

X ,...,X     i.i.d.  𝒩 (μ  ,σ2),
 1       n1             1  2
 Y1,...,Yn2  i.i.d.  𝒩 (μ2, σ ).
On souhaite tester l’hypothèse nulle H0: “le procédé n’a pas d’effet” contre “le procédé produit une augmentation significative de la quantité de pluie”, autrement dit
H   = {μ  =  μ }  contre   H  = {μ  >  μ }.
  0      1    2             1      1    2

Il est indispensable de s’assurer d’abord que les hypothèses de modèle que l’on a faites sont raisonnables:

  1. Les données peuvent-elles être considérées comme des réalisations de lois gaussienne pour les paramètres (moyenne, variance) appropriés?
  2. Peut-on considérer que les variances sont les mêmes?

Des éléments de réponses au point (1) seront abordés ultérieurement dans le cours. Le point (2) est abordé dans le paragraphe facultatif 1.5.

Nous admettons donc dans un premier temps ces deux hypothèses.

1.4 Comparaison des moyennes

Afin de se faire une idée des données et de l’écart entre les deux populations, on peut commencer par calculer les moyennes empiriques

         n                n
¯    -1-∑ 1       ¯   -1-∑ 2
X =  n      Xi,  Y  = n      Yi,
      1 i=1             2 i=1
ainsi que les variances empiriques
                                       [                  ]
            1   ∑n1               n1     1  n∑1       (  )2
S2X  =   -------   (Xi − X¯)2 = -------  ---   X2i −  X¯    ,
         n1 − 1 i=1             n1 − 1   n1 i=1
                ∑n2                   [    ∑n2      (  ) ]
 S2  =   ---1---   (Yi − ¯Y )2 =  --n2--- -1-   Y 2 −  ¯Y  2 .
  Y      n2 − 1 i=1             n2 − 1  n2 i=1  i
// taille des echantillons  
n1 = length(X);  
n2 = length(Y);  
// calcul des moyennes  
Xbar=sum(X)/n1  
Ybar=sum(Y)/n2  
// calcul de la somme des carres  
SSX=sum(X.^2);  
SSY=sum(Y.^2);  
// calcul des variances empiriques sans biais  
SX2 = SSX/(n1-1)- n1/(n1-1)*Xbar^2  
SY2 = SSY/(n2-1)- n2/(n2-1)*Ybar^2

On rappelle que les variables aléatoires X (resp. Y ) et (n1 1)SX2∕σ2 (resp. (n 2 1)SY 2∕σ2) sont indépendantes et de loi 𝒩(μ12∕n 1) (resp. 𝒩(μ22∕n 2)) et χ2(n 1 1) (resp. χ2(n 2 1)). Comme SX2 et S Y 2 sont indépendants, on en déduit que la loi de (n  − 1)S2
---1--2---X
     σ + (n  − 1)S2
--2---2---Y
    σ est la loi du χ2 à n 1 + n2 2 degrés de liberté.

Sous H0, la loi de ∘ --------
  -n1n2---
  n1 + n2( ¯X-−-¯Y)-
   σ est la loi gaussienne centrée réduite. On en déduit que sous H0,

    ∘  --------
       -n1n2--------( ¯X-−-Y¯)-----
T =    n1 + n2 ∘ (n1−-1)S2+-(n2−1)S2-
                 -----nX1+n2−-2---Y

suit une loi de Student de paramètre n1 + n2 2.

Comme sous H1, X Y converge p.s. vers μ1 μ2 > 0 quand min(n1,n2) →∞, et que lim min(n1,n2)→∞∘ --------
   -n1n2---
   n1 + n2 = +, on en déduit que sous H1, la statistique T diverge vers + quand min(n1,n2) →∞.

Question 1 Construire, à partir de la statistique T, un test d’égalité des moyennes en déterminant sa région critique, la valeur critique au niveau α et la p-valeur.

Question 2 En utilisant la fonction ttest2 définie ci-dessous, calculer numériquement la valeur de la statistique de test, la valeur critique au niveau 5%, et la p-valeur pour les données de pluie. Conclusion.

La fonction ttest2 (ttest2(n1,Xbar,SX2,n2,Ybar,SY2,alpha)) retourne la valeur de la statistique T observée, la valeur critique associée au niveau α et la p-valeur du test suggéré dans la question précédente. Cliquer sur le lien ci-dessus pour obtenir le code de la fonction, le sauvegarder, sous le nom ttest2.sce, dans le répertoire où vous utilisez scilab. Pour charger la fonction, utiliser la commande: getf  'ttest2.sce'. Pour appliquer la fonction, utiliser la commande:

ttest2(n1,Xbar,SX2,n2,Ybar,SY2,alpha)

On donne quelques indications pour la compréhension du code Scilab de ttest2 qui utilise la fonction de répartion (“cumulative distribution function” en anglais) de la loi de Student:

Pour avoir plus d’information sur la fonction cdft, on peut consulter le manuel en utilisant la commande help cdft.

Le but de ce qui suit est d’observer les résultats du tests quand les valeurs observées ou la taille de l’échantillon varient.

Question 3

  1. Faire varier Y , les autres valeurs étant fixées. Pour quelles valeurs de Y rejetez vous H0?
  2. Avec la valeur de Y ainsi déterminée, qu’observez vous si SY 2 augmente, si S Y 2 diminue?
  3. En reprenant les valeurs numériques du problème, qu’observez vous si n1 et n2 augmentent, si n1 et n2 diminuent jusqu’aux valeurs limites?

Question 4 Reprendre les questions 1 et 2, si l’hypothèse nulle est {μ1 μ2}.

1.5 (Facultatif) Comparaison des variances des 2 échantillons gaussiens

Pour s’assurer du point (2), on peut tout d’abord supposer que ces variances sont différentes, i.e.

                2
Xi   ∼   𝒩 (μ1,σ1),  i = 1,...,n1,
 Yj  ∼   𝒩 (μ2,σ2),  j = 1,...,n2,
                2
puis construire un test d’égalité des variances, afin de voir si la différence sur les variances (ou les écart-types) observée est significative.

Question 5 (Facultatif) Quelle est la loi de la statistique F = SX2∕S Y 2 sous l’hypothèse nulle d’égalité des variances? Quel est son comportement sous H1, quand n1 → ∞ et n2 →∞?

Question 6 (Facultatif) Construire, à partir de la statistique F, un test d’égalité des variances. Déterminer la région critique du test, évaluer la valeur critique au niveau α = 5% et la p-valeur pour les données de pluie? Conclusion.

On donne les commandes Scilab suivantes concernant la fonction de répartion de la loi de Fisher:

Pour avoir plus d’information sur la fonction cdff, on peut consulter le manuel en utilisant la commande help cdff.

2 (Facultatif) Analyse de la variance à un facteur

2.1 Les données

On dispose de données (réelles) relatives à 6 jeux de notes pour un même contrôle, corrigé par un unique correcteur, mais concernant six petites classes (et donc six enseignants) différentes. On désire savoir si les résultats des petites classes sont significativement différents.

|------|-------|------|-------|------|-------|
|Note--|Classe-|Note--|Classe-|-Note-|-Classe-|
|15.50  |  1    |15.50 |   1   | 12.00 |   5   |
|16.00  |  3    |11.75 |   2   | 18.50 |   2   |
|14.00  |  5    |18.50 |   4   | 8.00 |   6   |
|13.75  |  2    | 8.00  |   5   | 10.50 |   3   |
|      |       |      |       |      |       |
|17.00  |  2    |15.00 |   4   | 7.00 |   4   |
|11.00  |  6    |12.50 |   2   | 11.00 |   6   |
|16.50  |  1    |15.50 |   4   | 12.50 |   5   |
|12.75  |  5    |19.00 |   1   | 15.00 |   5   |
|15.00  |  6    |19.00 |   1   | 12.50 |   4   |
|13.50  |  5    | 9.50  |   2   | 16.50 |   2   |
|14.00  |  2    | 8.00  |   5   | 16.50 |   1   |
|      |       |      |       |      |       |
|19.50  |  1    |16.50 |   1   | 18.00 |   4   |
|15.50  |  5    | 9.50  |   6   | 9.50 |   5   |
|14.50  |  6    |19.00 |   1   | 13.50 |   6   |
|13.00  |  2    |17.00 |   6   | 13.00 |   3   |
|18.00  |  4    |17.00 |   5   | 16.00 |   2   |
|15.00  |  5    | 9.50  |   5   | 16.00 |   4   |
|14.50  |  6    |15.50 |   2   | 12.00 |   3   |
|14.50  |  5    |15.00 |   4   | 18.50 |   4   |
|      |       |      |       |      |       |
|18.50  |  3    |17.00 |   2   | 12.00 |   3   |
|17.50  |  1    |19.50 |   3   | 13.00 |   6   |
|18.00  |  3    |16.50 |   6   | 14.50 |   3   |
|11.00  |  3    | 7.00  |   4   | 12.00 |   3   |
|10.50  |  5    | 8.50  |   2   | 19.50 |   1   |
|10.00  |  1    |12.50 |   4   | 12.00 |   4   |
|16.00  |  3    |18.75 |   6   | 10.50 |   1   |
|      |       |      |       |      |       |
|12.25  |  4    |16.00 |   4   | 15.00 |   6   |
|9.50  |  1    |19.50 |   3   | 15.50 |   2   |
|12.00  |  5    |17.50 |   1   | 14.25 |   3   |
|13.50  |  3    |14.50 |   4   | 14.00 |   2   |
|19.50  |  1    | 9.00  |   6   | 15.50 |   3   |
|18.50  |  6    |      |       |      |       |
|------|--------------------------------------
|      |

2.2 Le modèle

On fait l’hypothèse que les notes sont distribuées suivant une loi gaussienne de même variance, et de moyenne dépendant de l’enseignant. Le modèle s’écrit donc

Xi,j = mi + 𝜀i,j,  i = 1,...,6,  j =  1,...,ni,
mi est la note moyenne (inconnue) liée à l’enseignant i, et les 𝜀i,j sont des v.a. indépendantes de loi 𝒩(02), où σ2 > 0 est inconnu également.
// notes  
X=[  
15.50; 15.50; 12.00; 16.00; 11.75; 18.50; 14.00; 18.50;  8.00; 13.75;  
 8.00; 10.50; 17.00; 15.00;  7.00; 11.00; 12.50; 11.00; 16.50; 15.50;  
12.50; 12.75; 19.00; 15.00; 15.00; 19.00; 12.50; 13.50;  9.50; 16.50;  
14.00;  8.00; 16.50; 19.50; 16.50; 18.00; 15.50;  9.50;  9.50; 14.50;  
19.00; 13.50; 13.00; 17.00; 13.00; 18.00; 17.00; 16.00; 15.00;  9.50;  
16.00; 14.50; 15.50; 12.00; 14.50; 15.00; 18.50; 18.50; 17.00; 12.00;  
17.50; 19.50; 13.00; 18.00; 16.50; 14.50; 11.00;  7.00; 12.00; 10.50;  
 8.50; 19.50; 10.00; 12.50; 12.00; 16.00; 18.75; 10.50; 12.25; 16.00;  
15.00;  9.50; 19.50; 15.50; 12.00; 17.50; 14.25; 13.50; 14.50; 14.00;  
19.50;  9.00; 15.50; 18.50];  
// groupe  
Grp=[  
1; 1; 5; 3; 2; 2; 5; 4; 6; 2; 5; 3; 2; 4; 4; 6; 2; 6; 1; 4; 5; 5; 1; 5;  
6; 1; 4; 5; 2; 2; 2; 5; 1; 1; 1; 4; 5; 6; 5; 6; 1; 6; 2; 6; 3; 4; 5; 2;  
5; 5; 4; 6; 2; 3; 5; 4; 4; 3; 2; 3; 1; 3; 6; 3; 6; 3; 3; 4; 3; 5; 2; 1;  
1; 4; 4; 3; 6; 1; 4; 4; 6; 1; 3; 2; 5; 1; 3; 3; 4; 2; 1; 6; 3; 6];  
 
 
 
 

On peut résumer les données.

// nombre de groupes  
k = max(Grp);  
// calcul de l'effectif n(i), de la moyenne moy(i)  
// et de la variance s2(i) du groupe i  
for i=1:k,  
 n(i) = length(X(Grp==i));  
 moy(i) = sum(X(Grp==i))/n(i);  
 ss(i)=sum(X(Grp==i).^2);  
 s2(i)=ss(i)/(n(i)-1)- n(i)/(n(i)-1)*moy(i)^2;  
end;  
[n,moy,s2] // effectif, moyenne et variance par groupe

Question 7 Quelles sont les estimations des mi, i = 1,, 6? Quelle est l’estimation de σ2?

Question 8 Peut-on considérer qu’il existe une différence significative entre les notes des classes, dues aux qualités pédagogiques des enseignants? Effectuer le test de l’hypothèse nulle H0 : “les enseignants ne présentent pas de différence” contre “c’est faux”. Écrire la table d’analyse de la variance et conclure.

On pourra utiliser les commandes suivantes:

// calcul de la moyenne generale  
mg=sum(X)/length(X);  
// calcul des sommes de carres  
SSM = sum(n.*(moy-mg).^2);  
SSE = sum(s2.*(n-1));

La fonction anova (anova(n,moy,s2,alpha)) retourne la table d’analyse de la variance associée. Les paramètres sont n, vecteur des tailles de chacun des groupes; moy, vecteur des moyennes de chacun des groupes; s2, vecteur des variances de chacun des groupes; alpha, niveau du test d’analyse de la variance. Cliquer sur le lien ci-dessus pour obtenir le code de la fonction, le sauvegarder sous le nom anova.sce dans le répertoire où vous utilisez scilab. Pour charger la fonction, utiliser la commande: getf  'anova.sce'. Pour appliquer la fonction, utiliser la commande:

anova(n,moy,s2, alpha)

Le but de ce qui suit est d’observer les résultats du tests quand les valeurs observées ou la taille de l’échantillon varient.

Question 9

  1. Faire varier le vecteur moy, les autres valeurs étant fixées. Pour quelles valeurs de moy rejetez vous H0?
  2. Qu’observez vous si le vecteur, s2, des variances augmente, s’il diminue?
  3. En reprenant les valeurs numériques du problème, qu’observez vous si le vecteur, n, des tailles des groupes augmente ou diminue?

Question 10 Si l’on rejette H0, il est naturel de déterminer les classes qui sont significativement différentes.

  1. Proposer des intervalles de confiance pour les notes moyennes pour chaque classe.
  2. Proposer un test de l’hypothèse nulle
    H0  = {mi =  mj }  contre  H1  = {mi ⁄=  mj },
    fondé sur la loi de Student, et l’utiliser afin de comparer quelques couples de moyennes.

3 Régression linéaire

3.1 Les données

L’exemple utilisé a déjà été proposé dans les exercices sur le chapitre 2. On dispose de données concernant l’âge (X1), le kilométrage en milliers de kms (X2), et le prix en milliers d’euros (Y ) pour un échantillon de voitures d’occasion d’un même type:

|----|-----|-----|
|X1--|-X2--|-Y---|
| 5  | 92  | 7.8  |
| 4  | 64  | 9.5  |
|    |     |     |
| 6  |124  | 6.4  |
| 5  | 97  | 7.5  |
| 5  | 79  | 8.1  |
| 5  | 76  | 9.0  |
| 6  | 93  | 6.1  |
|    |     |     |
| 6  | 63  | 8.7  |
| 2  | 13  |15.4 |
| 7  |111  | 6.4  |
| 7  |143  | 4.4  |
|----|-----------
|    |

Si l’on veut étudier la liaison entre les variables potentiellement explicatives et la variable à expliquer (le prix Y ), on peut commencer par visualiser les nuages de points (X1,Y ) et (X2,Y ):

donnees = [  
5  92  7.8; 4  64  9.5; 6 124  6.4; 5  97  7.5; 5  79  8.1; 5  76  9.0;  
6  93  6.1; 6  63  8.7; 2  13 15.4; 7 111  6.4; 7 143  4.4];  
// definition des variables  
X1=donnees(:,1); X2=donnees(:,2); Y=donnees(:,3);  
// allures des nuages  
xbasc();  
subplot(121), plot2d(X1,Y,-2,"111","Age",[1,0,8,16]);  
subplot(122), plot2d(X2,Y,-3,"111","Kms",[0,0,150,16]);

Question 11 Que suggèrent les deux nuages?

3.2 Régression linéaire simple

On propose tout d’abord le modèle M1:

Y  = β + α1X1  + 𝜀,
𝜀 est un échantillon de 𝒩(02).

Question 12 L’âge a-t-il une influence significative sur le prix?.

La fonction regression (regression(M,Y,alpha)) calcule, à partir de la donnée du modèle sous la forme Y = M𝜃 + 𝜀, où Y est la variable à expliquer et M la matrice définissant le modèle, les éléments usuels:

Cliquer sur le lien ci-dessus pour obtenir le code de la fonction, le sauvegarder sous le nom regression.sce dans le répertoire où vous utilisez scilab. Pour charger la fonction, utiliser la commande: getf  'regression.sce'. Pour appliquer la fonction, utiliser la commande:

regression(M,Y,alpha)

Pour le modèle, où on explique le prix en fonction de l’âge, on utilise les commaneds suivantes:

// MODELE (M1)  Y = theta1(0) + theta1(1) X1 + eps  
M1 = [ones(X1) X1]; // matrice des regresseurs  
// Pour obtenir les estimations de theta, sigma2,  
// et la p-valeur  
[theta1,sigma21,pval1]=regression(M1,Y,alpha);

On représente enfin la droite de régression sur le nuage:

xbasc();  
// trace de la droite  
z=1:0.1:8; t=theta1(1)+theta1(2)*z; xbasc();  
plot2d(X1,Y,-2,"111","Age",[1,0,8,16]);  
plot2d(z,t,[1],"000");

3.3 Régression linéaire multiple

On propose d’essayer à présent le modèle “complet”

Y = γ + α1X1  + α2X2  + 𝜀.

Question 13 Donner, à l’aide de la fonction regression des estimations sans biais des paramètres (γ,α12) et σ2.

Question 14 Test de l’utilité d’un régresseur: Tester, au niveau 5%, l’hypothèse nulle “le kilométrage n’a pas d’effet sur le prix” contre “c’est faux”.

La fonction regression n’effectue pas le test demandé (quel test effectue-t-elle?). Ecrire les commandes Scilab permettant de calculer la statistique et la p-valeur du test demandé, et conclure.

Question 15 A l’aide de la question précédente, déterminer quel sous-modèle, parmi les deux sous-modèles à un seul régresseur, explique le mieux le prix?