Analyse de la variance

# Analyse de la variance
## à 2 facteurs
### Marie-Pierre Etienne
### <a href="https://github.com/marieetienne" class="uri">https://github.com/marieetienne</a>
### 2020/09/11 (updated: 2023-10-02)

---

---
template: intro
## Notre exemple de fréquence cardiaque

On a mesuré la fréquence cardiaque de 20 femmes et 20 hommes.

```r
library(coursesdata)
data(freqdata)
```

--
<p class="question"> Les hommes et les femmes ont-ils la même fréquence cardiaque au repos ?</p>

Nous avions vu dans le cours sur l'analyse de la variance à 1 facteur qu'il n'y avait pas de différence entre homme et femmes.

---
template: intro
## Notre exemple de fréquence cardiaque

En plus du sexe, le jeu de données contient des informations sur le niveau de pratique sportive

<img src="anova2_files/figure-html/sport-1.png" width="60%" height="30%" style="display: block; margin: auto;" />
--
<p class=question> Y a t il un effet du sexe ou du sport sur la fréquence cardiaque ? </p>
<p class=question> Y a t il un effet du sexe  sur la fréquence cardiaque ? </p>
<p class=question> Y a t il un effet du sport sur la fréquence cardiaque ? </p>

---
template: intro
## Cadre général du modèle d'analyse de la variance à 2 facteurs

On étudie le lien entre  
- une variable quantitative notée `$Y$` (la fréquence cardiaque),
- et deux variables facteurs pouvant prendre respectivement `$I$` modalités (I=2 pour le Sexe)  et J modalités (J=5 pour le Sport)

Les données peuvent être visualisées à l'aide d'un boxplot.

--
<p class="question"> Y a t il un effet de l'un au moins des facteurs ?</p>
<p class="question"> Y a t il un effet du facteur 1 ?</p>
<p class="question"> Y a t il un effet du facteur 2 ?</p>

---
name: model
# Le modèle d'analyse de la variance à 2 facteurs

---
template: model

## Graphiquement

Une visualisation graphique du modèle d'analyse de la variance.

Comment imagine-t-on le processus aléatoire qui a conduit à nos données ?

---

---
count: false

---
count: false

---
count:false

---

---
count: false

---
count: false

---
template: model

Lequel de ces mécanismes est le plus crédible au vu des donées ?

---
template: model

## Version régulière du modèle Mcomp

`$$Y_{ijk} = \mu_{ij} +E_{ijk},\quad E_{ijk}\overset{ind}{\sim}\mathcal{N}(0, \sigma^2),$$`
avec 
- `$i=1,\ldots,I$` le niveau du facteur 1,
- `$j=1,\ldots,J$` le niveau du facteur 2,
- `$k= 1,\ldots, n_{ij}$` le numéro de l'individu dans le groupe `$(i,j)$`, 
- `$n_{ij}$` le nombre d'individus dans le groupe `$(i,j)$` et `$n=\sum_ j\sum_i n_{ij}$` le nombre total d'individus,
- `$\mu_{ij}$` le comportement moyen du groupe `$(i,j)$`,
- `$\sigma^2$` la variance commune à tous les groupes.

### Une écriture équivalente

`$$Y_{ik} \overset{ind}{\sim}\mathcal{N}(\mu_{ij}, \sigma^2).$$`

### Nombre de paramètres du modèle

- `$I*J$` paramètres de moyenne  `$(\mu_{11}, \mu_{12}, \ldots,  \mu_{1J}, \ldots \mu_{IJ})$`; 
- 1 paramètre de variance `$\sigma^2$`

---
template: model

## Version régulière du modèle Mcomp sur l'exemple 1

`$$Y_{ijk} = \mu_{ij} +E_{ijk},\quad E_{ijk}\overset{ind}{\sim}\mathcal{N}(0, \sigma^2),$$`
avec 
- `$i=1,\ldots,2$` le niveau du facteur Sexe,
- `$j=1,\ldots,5$` le niveau du facteur Sport,
- `$k= 1,\ldots, 4$` le numéro de l'individu dans le groupe `$(i,j)$`, 
- `$4$` le nombre d'individus dans le groupe `$(i,j)$` et `$n= 40=\sum_ j\sum_i n_{ij}$` le nombre total d'individus,
- `$\mu_{ij}$` le comportement moyen du groupe `$(i,j)$`,
- `$\sigma^2$` la variance commune à tous les groupes.

`$I=2, J= 5$` et  la convention `$i=1$` pour les femmes et `$i=2$` pour les hommes.

### Nombre de paramètres
- 10 paramètres de moyenne
- 1 paramètre de variance

---
template: model
## Version singulière du modèle du modèle Mcomp

`$$Y_{ijk} = \mu + \alpha_i +\beta_j +\gamma_{ij} +E_{ijk},\quad E_{ijk}\overset{ind}{\sim}\mathcal{N}(0, \sigma^2),$$`
avec 
- `$i=1,\ldots,I$` le niveau du facteur 1,
- `$j=1,\ldots,J$` le niveau du facteur 2,
- `$k= 1,\ldots, n_{ij}$` le numéro de l'individu dans le groupe `$(i,j)$`, 
- `$\mu$`  le comportement moyen de référence
- `$\alpha_i$`  l'effet différentiel du niveau `$i$` 
- `$\beta_j$`  l'effet différentiel du niveau `$j$` 
- `$\gamma_{ij}$`  l'effet différentiel de la rencontre des niveaux `$(i,j),$`
- `$\sigma^2$` la variance commune à tous les groupes.

### Nombre de paramètres du modèle

- `$1 + I + J + I*J$` paramètres de moyenne  
- 1 paramètre de variance `$\sigma^2$`

#### La version dans les logiciels et qui permet de séparer les effets des différens facteurs.

---
template: model
## Lien entre les deux versions du même modèle

<table style="width:100%">
  <tr>
    <th>Groupe</th>
    <th>V. régulière</th>
    <th>V. singulière</th>
  </tr>
  <tr>
    <td>1</td>
    <td> `$\mu_{11}$` </td>
    <td> `$\mu +\alpha_1 + \beta_1 + \gamma_{11}$` </td>
  </tr>
  <tr>
    <td>2</td>
    <td> `$\mu_{12}$` </td>
    <td> `$\mu +\alpha_1 + \beta_2 + \gamma_{12}$` </td>
  </tr>
    <tr>
    <td> </td>
    <td>  </td>
    <td>   </td>
  </tr>
  </tr>
    <tr>
    <td> </td>
    <td>  </td>
    <td>   </td>
  </tr>
  </tr>
    <tr>
    <td> </td>
    <td>  </td>
    <td>   </td>
  </tr>
    <tr>
    <td>I</td>
    <td> `$\mu_{IJ}$` </td>
    <td> `$\mu +\alpha_I +\beta_J +\gamma_{IJ}$` </td>
  </tr>
</table>

<a class=care> Problème </a>  Le problème est identique au modèle d'analyse de la variance, le modèle sous cette forme est  <a style="font-weight:400;"> indéterminé</a>.

#### Solution : ajouter des contraintes
Par défaut dans R :
`$\alpha_1=0,$`  `$\beta_1=0,$` et `$\gamma_{1j}=0$` pour tous les `$j$` et `$\gamma_{i1}=0$` pour tous les `$i$`.

Ce qui correspond à `$1 + 1 + J + (I-1)=I  + J +1$` contraintes.

#### Nombre de paramètres 
- On a donc `$1 + I + J + IJ$` paramètres dont `$I  + J +1$` sont contraints, il reste donc
`$IJ$` à estimer pour le comportement moyen
- 1 paramètre de variance `$\sigma^2$`.

---
template: model
## Que signifient ces contraintes (à partir de l'exemple)

Le groupe `$(1,1)$`, c'est à dire les femmes très sportives ont un comportement moyen dans le modèle égal à `$\mu +\alpha_1 +\beta_1 + \gamma_{11},$` ce qui vaut `$\mu$` étant données les contraintes.

Le groupe `$(1,j)$` désigne les femmes à la pratique sportive identifiée `$j$`. En moyenne leur fréquence cardiaque dans le modèle est `$\mu + \alpha_1 +\beta_j +\gamma_{1j}$`, soit étant données les contraintes `$\mu +\beta_j.$` `$\beta_j$` est donc l'effet différentiel du sport `$j$` par rapport au sport `$1$`  <a class=care> pour les femmes</a>.

Le groupe `$(i,1)$` désigne les individus de sexe `$i$`  à la pratique sportive identifiée `$1$`. En moyenne leur fréquence cardiaque dans le modèle est `$\mu + \alpha_i +\beta_1 +\gamma_{i1}$`, soit étant données les contraintes `$\mu +\alpha_i.$` `$\alpha_i$` est donc l'effet différentiel du sexe `$i$`  par rapport au sexe de référence `$1$`   <a class=care> pour les sportifs de niveau 1</a>.

Les termes `$\gamma_{ij}$` sont difficiles à interpréter et seront discuter en TP.

---
template: model
## Sous forme matricielle
 `$$\bf{Y = X\theta + E}$$`
### Forme régulière

`$$Y=\begin{pmatrix}
Y_{111}\\
Y_{112}\\
\vdots\\
Y_{11n_{11}}\\
Y_{121}\\
\vdots\\
Y_{12n_{12}}\\
\vdots\\
Y_{ijk}\\
\vdots\\
Y_{IJn_{IJ}}\end{pmatrix},
 \quad
{\bf{X}} =\overset{\color{gray}{\begin{matrix}\mu_{11} & \mu_{12} & \ldots & \ldots &\mu_{IJ}\end{matrix}}}{\begin{pmatrix}
1 & 0 & 0 & \ldots 0\\
1 & 0 & 0 & \ldots 0\\
\vdots\\
1 & 0 & 0 & \ldots 0\\
0 & 1 & 0 & \ldots 0 \\
\vdots\\
0 & 1 & 0 & \ldots 0 \\
\vdots\\
0 & \ldots  & 1 & \ldots 0 \\\vdots\\
0 & 0 & 0 & \ldots 1 \end{pmatrix}},\quad
{\bf{\theta}} =\begin{pmatrix}
\mu_{11}\\
\mu_{12}\\
\vdots\\
\mu_{1J}\\
\mu_{21}\\
\vdots\\
\mu_{IJ}\end{pmatrix}, \quad{\bf{E}} = \overset{}{\begin{pmatrix}
E_{111}\\
E_{112}\\
\vdots\\
E_{11n_{11}}\\
E_{121}\\
\vdots\\
E_{12n_{12}}\\
\vdots\\
E_{IJn_{IJ}}\end{pmatrix}}$$`
---
template: model

## Sous forme matricielle
 `$$\bf{Y = X\theta + E}$$`
### Forme singulière

<p class=question> Ecrire le modèle en version singulière sous forme matricielle en intégrant les contraintes.</p>

---
template: model

## sur l'exemple de la fréquence cardiaque

```r
Mcomp <- lm(freqC ~ Sexe + Sport + Sexe:Sport, data = freqdata)
model.matrix(Mcomp)
```

```
   (Intercept) SexeM Sport SexeM:Sport
1            1     1     1           1
2            1     0     1           0
3            1     1     1           1
4            1     0     1           0
5            1     1     1           1
6            1     0     1           0
7            1     1     1           1
8            1     0     1           0
9            1     1     2           2
10           1     0     2           0
11           1     1     2           2
12           1     0     2           0
13           1     1     2           2
14           1     0     2           0
15           1     1     2           2
16           1     0     2           0
17           1     1     3           3
18           1     0     3           0
19           1     1     3           3
20           1     0     3           0
21           1     1     3           3
22           1     0     3           0
23           1     1     3           3
24           1     0     3           0
25           1     1     4           4
26           1     0     4           0
27           1     1     4           4
28           1     0     4           0
29           1     1     4           4
30           1     0     4           0
31           1     1     4           4
32           1     0     4           0
33           1     1     5           5
34           1     0     5           0
35           1     1     5           5
36           1     0     5           0
37           1     1     5           5
38           1     0     5           0
39           1     1     5           5
40           1     0     5           0
attr(,"assign")
[1] 0 1 2 3
attr(,"contrasts")
attr(,"contrasts")$Sexe
[1] "contr.treatment"
```

---
template: model

## sur l'exemple de la fréquence cardiaque
*Attention, il manque une colonne*

```r
Mcomp <- lm(freqC ~ Sexe + Sport_Fact + Sexe:Sport_Fact, data = freqdata)
model.matrix(Mcomp)
```

```
   (Intercept) SexeM Sport_Fact2 Sport_Fact3 Sport_Fact4 Sport_Fact5
1            1     1           0           0           0           0
2            1     0           0           0           0           0
3            1     1           0           0           0           0
4            1     0           0           0           0           0
5            1     1           0           0           0           0
6            1     0           0           0           0           0
7            1     1           0           0           0           0
8            1     0           0           0           0           0
9            1     1           1           0           0           0
10           1     0           1           0           0           0
11           1     1           1           0           0           0
12           1     0           1           0           0           0
13           1     1           1           0           0           0
14           1     0           1           0           0           0
15           1     1           1           0           0           0
16           1     0           1           0           0           0
17           1     1           0           1           0           0
18           1     0           0           1           0           0
19           1     1           0           1           0           0
20           1     0           0           1           0           0
21           1     1           0           1           0           0
22           1     0           0           1           0           0
23           1     1           0           1           0           0
24           1     0           0           1           0           0
25           1     1           0           0           1           0
26           1     0           0           0           1           0
27           1     1           0           0           1           0
28           1     0           0           0           1           0
29           1     1           0           0           1           0
30           1     0           0           0           1           0
31           1     1           0           0           1           0
32           1     0           0           0           1           0
33           1     1           0           0           0           1
34           1     0           0           0           0           1
35           1     1           0           0           0           1
36           1     0           0           0           0           1
37           1     1           0           0           0           1
38           1     0           0           0           0           1
39           1     1           0           0           0           1
40           1     0           0           0           0           1
   SexeM:Sport_Fact2 SexeM:Sport_Fact3 SexeM:Sport_Fact4 SexeM:Sport_Fact5
1                  0                 0                 0                 0
2                  0                 0                 0                 0
3                  0                 0                 0                 0
4                  0                 0                 0                 0
5                  0                 0                 0                 0
6                  0                 0                 0                 0
7                  0                 0                 0                 0
8                  0                 0                 0                 0
9                  1                 0                 0                 0
10                 0                 0                 0                 0
11                 1                 0                 0                 0
12                 0                 0                 0                 0
13                 1                 0                 0                 0
14                 0                 0                 0                 0
15                 1                 0                 0                 0
16                 0                 0                 0                 0
17                 0                 1                 0                 0
18                 0                 0                 0                 0
19                 0                 1                 0                 0
20                 0                 0                 0                 0
21                 0                 1                 0                 0
22                 0                 0                 0                 0
23                 0                 1                 0                 0
24                 0                 0                 0                 0
25                 0                 0                 1                 0
26                 0                 0                 0                 0
27                 0                 0                 1                 0
28                 0                 0                 0                 0
29                 0                 0                 1                 0
30                 0                 0                 0                 0
31                 0                 0                 1                 0
32                 0                 0                 0                 0
33                 0                 0                 0                 1
34                 0                 0                 0                 0
35                 0                 0                 0                 1
36                 0                 0                 0                 0
37                 0                 0                 0                 1
38                 0                 0                 0                 0
39                 0                 0                 0                 1
40                 0                 0                 0                 0
attr(,"assign")
 [1] 0 1 2 2 2 2 3 3 3 3
attr(,"contrasts")
attr(,"contrasts")$Sexe
[1] "contr.treatment"

attr(,"contrasts")$Sport_Fact
[1] "contr.treatment"
```

---
class: inverse
name: pause
# Pause

<br><br><br><br>
<p style="color:#B40F20;font-size:35px;text-align:center;">Prenons une petite pause !!!</p>

---
name: parametre
# Estimation des paramètres

---
template: parametre
## Estimation des paramètres du modèle version matricielle

Le modèle sous forme matricielle s'écrit

`$$\bf{Y = X\theta + E}.$$`
--

### Estimation de `$\theta$`

`$$\hat{\theta} = (X^\intercal X )^{-1} X^\intercal Y_{obs}.$$`

### Estimateur de `$\theta$`

`$$T = (X^\intercal X )^{-1} X^\intercal Y.$$`
--

### Loi de l'estimateur de `$\theta$`

`$$T  \sim \mathcal{N}_{I}\left(\theta, \sigma^2 (X^\intercal X )^{-1}\right).$$`

---
template: parametre
## Le paramètre de variance

La somme des carrés résiduelles s'écrit sous la forme

`$$RSS = || Y- X \hat{\theta} ||^2$$`

### Estimateur de la variance

$$S^2 =\frac{1}{DF_{res}} RSS, $$
est un <a class=care> estimateur sans bias de  `$\sigma^2$` </a> .

Dans le cas du modèle d'analyse de la variance à 2 facteurs `$DF_res=n-IJ$` (n observations et IJ paramètres de moyennes à estimer, le nombre de composantes dans le vecteur `$\theta$`)

## Estimation de `$\sigma^2$`

`$$\hat{\sigma}^2 =\frac{1}{n-IJ} RSS_{obs}.$$`
---
template: parametre
## sur l'exemple de la fréquence cardiaque

```r
summary(Mcomp)
```

```

Call:
lm(formula = freqC ~ Sexe + Sport_Fact + Sexe:Sport_Fact, data = freqdata)

Residuals:
   Min     1Q Median     3Q    Max 
 -3.00  -1.00   0.00   1.00   3.25

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
(Intercept)        81.0000     0.8229  98.438  < 2e-16 ***
SexeM              -2.0000     1.1637  -1.719 0.095979 .  
Sport_Fact2         1.5000     1.1637   1.289 0.207247    
Sport_Fact3         1.5000     1.1637   1.289 0.207247    
Sport_Fact4         4.7500     1.1637   4.082 0.000305 ***
Sport_Fact5         8.7500     1.1637   7.519 2.21e-08 ***
SexeM:Sport_Fact2   0.5000     1.6457   0.304 0.763361    
SexeM:Sport_Fact3   6.5000     1.6457   3.950 0.000438 ***
SexeM:Sport_Fact4   1.5000     1.6457   0.911 0.369316    
SexeM:Sport_Fact5   2.7500     1.6457   1.671 0.105119    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.646 on 30 degrees of freedom
Multiple R-squared:  0.8689,	Adjusted R-squared:  0.8296 
F-statistic: 22.09 on 9 and 30 DF,  p-value: 6.689e-11
```

---
class: inverse
name: pause
# Pause

<br><br><br><br>
<p style="color:#B40F20;font-size:35px;text-align:center;">Prenons une petite pause !!!</p>

---
name: modcomp
# Test du modèle complet

---
template: modcomp
## Rappel exemple fréquence cardiaque (exemple 1)

On a mesuré la fréquence cardiaque de 20 femmes et 20 hommes et leur niveau de pratique sportive

--
<p class="question"> Le sexe ou le niveau de pratique sportive sont ils liés à la fréquence carsdiaque au repos ?</p>

---
template: modcomp
## Sous forme de comparaison de modèle

--
<p class="question"> Le modèle Mcomp est il plus pertinent que le modèle M0 ?</p>

---
template: modcomp
## Hypothèses du test

On va donc opposer une hypothèse de travail `$H_0$` contre une hypothèse alternative `$H_1$`. `$H_0$` peut donc prendre différentes formes:

`$$\begin{align} 
H_0 & =\left \lbrace \mbox{Pas de différence entre les différents groupes }\right\rbrace\\
    & =\left \lbrace  \mbox{pour tout }i, \alpha_i =0, \mbox{pour tout }j,\beta_j=0  \mbox{ et pour tout }(i,j),\gamma_{ij}=0   \right\rbrace\\
    & =\left \lbrace  M_{comp} \mbox{ est équivalent à } M0 \right\rbrace.
\end{align}$$`

`$H_1$` prend les formes équivalentes suivantes

`$$\begin{align} 
H_1 & =\left \lbrace \mbox{Au moins 1 groupe est différent des autres}\right\rbrace\\
    & =\left \lbrace  \mbox{Il existe un }i, \alpha_i  \ne 0  \mbox{ ou un }j, \beta_j \ne 0 \mbox{ ou un } (i,j), \gamma_{ij} \ne 0  \right\rbrace\\
    & =\left \lbrace  M_{comp} \mbox{ est préférable à } M0 \right\rbrace.
\end{align}$$`

Sous `$H_0$`, 
`$$F= \frac{\frac{SS_{M_{comp}}}{IJ-1}}{\frac{RSS}{n-IJ}} \underset{H_0}{\sim}\mathcal{F}(IJ-1, n-IJ)$$`

---
template: modcomp
## Loi de la statistique de test sous `$H_0$` - graphiquement

Sous `$H_0$` la loi de distribution de `$F$` est

---

---
count: false

---
count: false

---
name: variance_decomposition
# Decomposition de SSM

---
template: variance_decomposition
## Visualisation graphique de la décomposition de la variance

<br>