Modèle linéaire

class: center, middle, inverse, title-slide

# Modèle linéaire
## généralisé
### Marie-Pierre Etienne
### <a href="https://github.com/marieetienne" class="uri">https://github.com/marieetienne</a>
### 2020/09/11 (updated: 2023-10-02)

---

name: intro
# Introduction

---
template: intro
## Etude de l'aire de répartition de la grenouille à pattes rouges

237 points d'eau références

Pour chacun
- `Source` source de l'information (Museum, Literature, PersCom ou Field Note),
- `Source2` précision sur la source,
- `Status` présence ou l'absence de grenouilles
- `Latitude` latitude
- `Longitude` longitude
--

```r
data("grenouilles")
head(grenouilles)
```

```
  Source Source2 Status Latitude Longitude
1 Museum    AMNH      A 37.88222  -119.855
2 Museum     CAS      P 38.10111  -122.856
3 Museum     CAS      P 34.03269  -119.702
4 Museum     CAS      P 35.39003  -120.608
5 Museum     CAS      P 35.30360  -120.376
6 Museum     CAS      A 35.10289  -120.605
```

---
template: intro
## Etude de l'aire de répartition de la grenouille à pattes rouges

---
template: intro
## Etude de l'aire de répartition de la grenouille à pattes rouges

---
template: intro
## Etude de l'aire de répartition de la grenouille à pattes rouges

<a class=question> Peut on expliquer la présence de la grenouille à pattes rouges par la latitude ou la longitude ? </a>

---
template: intro
## Cadre général de la régression logistique

On cherche à expliquer la 
- variable Y, qualitative (2 modalités) par
- des variables quantitatives

--
## Cadre du modèle linéaire généralisé

On cherche à expliquer la 
- variable Y,qui est  mal remprésentée par une loi normale
- des variables quantitatives et ou qualitatives.

---
template: intro
## Sur l'exemple Grenouille

La variable à expliquer `$Y$` est la présence ou l'absence de grenouille

Les variables explicatives sont la latitude et la longitude

---
name: model
# Le modèle linéaire généralisé

---
template: model
## Rappel sur le modèle linéaire

`$$\begin{align}
Y_k &\overset{ind}{\sim}\mathcal{N}( \mathbb{E}(Y_k) , \sigma^2)\\
& \mathbb{E}(Y_k) = X_i \theta,
\end{align}$$`

et `$X_i$` désigne la ligne `$i$` de la matrice `$X$`.

--
## Le modèle linéaire généralisé

`$$\begin{align}
Y_k &\overset{ind}{\sim}\mathcal{L}( \mathbb{E}(Y_k) , \phi )\\
& g(\mathbb{E}(Y_k)) = X_i \theta,
\end{align}$$`

--
## Le modèlé linéaire est un modèle linéaire généralisé
`$$\mathcal{L} =\mathcal{N},\quad g=Id$$`

---
template: model
## Sur l'exemple Grenouille

- `$Y_k$` présence ou absence de la grenouille au site `$k$`

`$$Y_k \overset{ind}{\sim}\mathcal{B}(p_k)$$`

- `$p_k$` est la probabilité de présence au site `$k$`.

On aimerait écrire
`$$p_k = \beta_0 + \beta_1 x_{k}^{(1)} + \beta_2 x_{k}^{(2)}.$$`

Mais `$p_k$` est  un nombre entre 0 et 1 !!

---
count: false
template: model
## Sur l'exemple Grenouille

- `$Y_k$` présence ou absence de la grenouille au site `$k$`

`$$Y_k \overset{ind}{\sim}\mathcal{B}(p_k)$$`

- `$p_k$` est la probabilité de présence au site `$k$`.

`$$logit(p_k) = \beta_0 + \beta_1 x_{k}^{(1)} + \beta_2 x_{k}^{(2)}$$`

`$$logit(p) =\log{ \left (p /  (1-p) \right))}.$$`

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

---
template: estimation
## Rappel modèle linéaire

- Estimation des paramètres en minimisant RSS. 
- Ce qui revient à maximiser la vraisemblance i.e
`$$L(\theta) = \prod_{k} \mathbb{P}_{\theta}(y_k)$$`
--

- ou en encore la log vraisemblance
`$$\ell(\theta) = \sum_{k} \log\left \{  \mathbb{P}_{\theta}(y_k) \right \} = \frac{n}{2} \log{\sigma^2} - \frac{1}{2} \sum_{k=1}^n \frac{(y_k - \hat{y}_k)^2}{\sigma^2}$$`
--

- On obtient une formule explicite

---
template: estimation
## Modèle linéaire généralisé

- `$\hat{\theta}$` qui maximise  la log vraisemblance

`$$\ell(\theta) = \sum_{k} \log\left \{  \mathbb{P}_{\theta}(y_k) \right \}.$$`
--

Malheureusement dans le cas général pas de formule explicite.

On utilise une optimsation numérique.

---
template: estimation
## Loi des estimateurs

Grâce à des résultats mathématiques, on a quand même

`$${T} \underset{n\to\infty}{\sim} \mathcal{N}(\theta, \Sigma),$$`
`$\Sigma$` connu.

On peut faire des tests asymptotiques (si `$n$` grand).

---
template: estimation
## Sur l'exemple des grenouilles

```

Call:
glm(formula = Status ~ Latitude + Longitude, family = binomial(link = "logit"), 
    data = grenouilles)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-2.0575  -0.4404  -0.2393   0.3011   2.6435

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept) 265.0912    32.1226   8.252  < 2e-16 ***
Latitude      2.1157     0.3171   6.671 2.53e-11 ***
Longitude     2.8382     0.3523   8.056 7.87e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

(Dispersion parameter for binomial family taken to be 1)

Null deviance: 326.74  on 235  degrees of freedom
Residual deviance: 146.81  on 233  degrees of freedom
AIC: 152.81

Number of Fisher Scoring iterations: 6
```

---
name: test
# Test de comparaison de modèles
---
template: test
## Rappel modèle linéaire

Les modèles sont comparés à l'aide de différences de RSS. Par exemple comparaison entre `$M_0$` et `$M_1$`:

`$$RSS_{M_0} - RSS_{M_1}$$`

## Modèle linéaire généralisé

On définit la déviance d'un modèle `$D_{M1}$` et l'idée sera de mesurer l'écart entre deux modèles par l'écart de déviance.

La déviance représente, comme RSS, la variabilité non capturée par le modèle

---
template: test
## Modèle linéaire généralisé

### Modèle saturé

On appelle modèle saturé `$M_s$`, le modèle ayant autant de paramètres que d'observations. C'est le modèle le plus ajusté possible. On note `$\ell_{Sat}$` la vraisemblance de ce modèle.

### Déviance

La déviance pour un modèle `$M_1$` est donné par

`$$D_{M_1}= -2(\ell(\hat{\theta}) - \ell_{Sat} )$$`

### Comparaison de modèles

Sous `$H_0 : \left\lbrace M_1 \mbox{ et } M_0 \mbox{ sont équivalents}\right\rbrace,$`

`$$D_{M0} - D_{M1} \underset{n\to\infty}{\sim} \chi^2(p),$$`
où le nombre de degrés de liberté `$p$` est le nombre de paramètres d'écart entre les deux modèles.

---
template: test
## Graphiquement

---
template: test
## Test du modèle complet sur l'exemple des grenouilles

```r
glm0 <- glm(Status ~ 1, family = binomial(link = logit), data = grenouilles)
anova(glm0, glm_comp, test = 'Chisq')
```

```
Analysis of Deviance Table

Model 1: Status ~ 1
Model 2: Status ~ Latitude + Longitude
  Resid. Df Resid. Dev Df Deviance  Pr(>Chi)    
1       235     326.74                          
2       233     146.81  2   179.93 < 2.2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```

---
template: test
## Décomposition de type 1

```r
anova(glm_comp, test = 'Chisq')
```

```
Analysis of Deviance Table

Model: binomial, link: logit

Response: Status

Terms added sequentially (first to last)

Df Deviance Resid. Df Resid. Dev  Pr(>Chi)    
NULL                        235     326.74              
Latitude   1    8.929       234     317.81  0.002806 ** 
Longitude  1  171.002       233     146.81 < 2.2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```

```r
glm1 <- glm(Status ~ Latitude , family = binomial(link = "logit"), data = grenouilles)
glm2 <- glm(Status ~ Longitude, family = binomial(link = "logit"), data = grenouilles)
anova(glm0, glm1, test = 'Chisq')
```

```
Analysis of Deviance Table

Model 1: Status ~ 1
Model 2: Status ~ Latitude
  Resid. Df Resid. Dev Df Deviance Pr(>Chi)   
1       235     326.74                        
2       234     317.81  1   8.9292 0.002806 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```

---
template: test
## Décomposition de type 1

---
template: test
## Décomposition de type 2

```r
library(car)
Anova(glm_comp, test = 'LR')
```

```
Analysis of Deviance Table (Type II tests)

Response: Status
          LR Chisq Df Pr(>Chisq)    
Latitude    100.04  1  < 2.2e-16 ***
Longitude   171.00  1  < 2.2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```

```r
anova(glm2, glm_comp, test = 'Chisq')
```

```
Analysis of Deviance Table

Model 1: Status ~ Longitude
Model 2: Status ~ Latitude + Longitude
  Resid. Df Resid. Dev Df Deviance  Pr(>Chi)    
1       234     246.85                          
2       233     146.81  1   100.04 < 2.2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```

---
template: test
## Exercice

Reprendre le dessin de décomposition du type 1 et l'adapter au type 2.

---
# Fin du cours

## Des questions ?

Ecrivez moi !!!

## Bon travail