ANÁLISE DE DADOS CATEGORIZADOS STA13833

class: center, middle, inverse, title-slide

.title[
# ANÁLISE DE DADOS CATEGORIZADOS STA13833
]
.subtitle[
## Modelos Loglineares
]
.author[
### Nátaly A. Jiménez Monroy
]
.institute[
### LECON/DEST - UFES
]
.date[
### Vitória. ES
]

---

class: animated, fadeIn

#Características

* Modelam os padrões de associação e de interação entre variáveis categóricas.

* São adequados para esquemas amostrais de Poisson, Multinomiais e Produto de Multinomiais.

* Apropriados para situações em que não há distinção clara entre resposta e variáveis explicativas ou quando há mais de uma resposta.

>**Essa é a principal diferença entre os modelos logísticos e os loglineares**.

* Assumem que as variáveis discretas são nominais, mas é possível fazer ajustes para considerar ordinalidade e com dados pareados.

* São mais gerais que os modelos logit, mas alguns modelos loglineares têm correspondência direta com os modelos logit.

---
class: animated, fadeIn

#Modelos Loglineares para tabelas de duas vias - I

Sejam `$\mu_{ij}$` as contagens esperadas, `$E(n_{ij})$`, em uma tabela `$I\times J$` formada pelo cruzamento das variáveis `$A$` e `$B$`. O objetivo é modelar as contagens `$\mu_{ij}=n\pi_{ij}$`.

Sejam `$N=I\times J$`, as contagens das caselas, observações independentes de uma variável aleatória Poisson, `$n_{ij}=Poisson(\mu_{ij})$`. Um modelo loglinear análogo à ANOVA de duas vias com interação é dado por

`$$\log(\mu_{ij})=\lambda+\lambda_i^A+\lambda_j^B+\lambda_{ij}^{AB},$$`
onde `$i=1,\ldots,I; j=1,\ldots,J$` são os níveis das variáveis categóricas `$A$` e `$B$`. Com o objetivo de garantir a unicidade das estimativas, impõem-se as restrições `$\sum_i \lambda_i=\sum_j \lambda_j=\sum_i\sum_j \lambda_{ij}=0$`.

>Esse modelo é superparametrizado porque o termo `$\lambda_{ij}$` já tem `$I\times J$` parâmetros correspondentes às médias das caselas `$\mu_{ij}$`. A constante `$\lambda$` e os efeitos principais `$\lambda_i$` e `$\lambda_j$` correspondem a `$1+I+J$` parâmetros adicionais.

Nesse contexto, há dois tipos de modelos a serem considerados

* Modelo de independência (A, B).
 * Modelo saturado (AB).
 
 
 
 
---
class: animated, fadeIn

##Modelo de Independência - I

A independência pode ser escrita em termos de probabilidades das caselas como o produto das probabilidades marginais

`$$\pi_{ij}=\pi_{i+}\pi_{+j}\quad i=1,\ldots,I, j=1,\ldots,J$$`
e em termos das frequências das caselas

`$$\mu_{ij}=n\pi_{ij}=n\pi_{i+}\pi_{+j}\quad i=1,\ldots,I, j=1,\ldots,J.$$`

Tomando logaritmos obtemos o modelo loglinear de independência

`$$\log(\mu_{ij})=\lambda + \lambda_i^A + \lambda_j^B,$$`
 onde `$A$` e `$B$` denotam as duas variáveis categóricas.
 
 
---
class: animated, fadeIn

##Modelo de Independência - II

**Observações:**
 
  *  `$\lambda$` representa o efeito global, ou média geral dos logaritmos das contagens esperadas, garante que `$\sum_i\sum_j \mu_{ij}=n$`, isto é, as contagens esperadas sob o modelo ajustado totalizam o tamanho amostral `$n$`.
 
  * `$\lambda_i^A$` representa o efeito principal da variável A. Garante que `$\sum_j \mu_{ij}=n_{i+}$`, isto é, os totais marginais sob o modelo ajustado totalizam as contagens marginais observadas. Representa o efeito da classificação na linha `$i$`.
 
  * `$\lambda_j^B$` representa o efeito principal da variável B. Garante que `$\sum_i \mu_{ij}=n_{+j}$`. Representa o efeito da classificação na coluna `$j$`.
 
  * Restrições: `$\lambda_I^A=\lambda_J^B=0$`, ou alternativamente, `$\sum_i \lambda_i^A=\sum_j \lambda_j^B=0$`.
 
--

Os valores estimados por Máxima Verossimilhança para as contagens são os mesmos valores esperados sob o teste de independência em tabelas de duas vias, i.e., `$E(\mu_{ij})=n_{i+}n_{+j}/n$`. Portanto, as estatísticas `$\chi^2$` e `$G^2$` (deviance) para independência são testes de bondade do ajuste para o modelo loglinear de independência, onde as hipóteses são: o modelo de independência é o correto vs. o modelo saturado é o correto. Esse modelo também implica que todos os odds ratios são iguais a 1.

---
class: animated, fadeIn

##Modelo de Independência - III

Dependendo do tipo de restrição imposta, podemos obter diferentes estimativas dos parâmetros (i.e, diferentes valores dos `$\lambda s$`). Entretanto, as diferenças (log-odds) são únicos: `$\lambda_i^A - \lambda_{i'}^A$` e `$\lambda_j^B - \lambda_{j'}^B$`, onde o subíndice `$i$` denota um nível da variável categórica `$A$` e `$i'$` denota outro nível da mesma varíável; similarmente para `$B$`.

Dessa forma, os odds ratios também são únicos.

`$$\begin{align*}
   \log(odds)&=\log\left(\frac{\mu_{i1}}{\mu_{i2}}\right)=\log(\mu_{i1})-\log(\mu_{i2})\\
   &=(\lambda+\lambda_i^A+\lambda_1^B)-(\lambda+\lambda_i^A+\lambda_2^B)\\
   &=(\lambda_1^B-\lambda_2^B)\\
   \end{align*}$$`
   
`$$\Longrightarrow odds=\exp(\lambda_1^B-\lambda_2^B)$$`

---
class: animated, fadeIn

##Modelo de Independência - IV

Sob o modelo de independência, espera-se que o log-odds seja 0, isto é, que o OR=1. Assim,

`$$\begin{align*}
\log(odds)&=\log\left(\frac{\mu_{11}\mu_{22}}{\mu_{12}\mu_{21}}\right)\\
&=\log(\mu_{11})+\log(\mu_{22})-\log(\mu_{12})-\log(\mu_{21})\\
&=\lambda+\lambda_1^A+\lambda_1^B+\lambda+\lambda_2^A+\lambda_2^B-\lambda-\lambda_1^A-\lambda_2^B-\lambda-\lambda_2^A-\lambda_1^B\\
&=0
\end{align*}$$`

>**Observação:** 
O odds ratio mede a força da associação e depende apenas do termo de interação `${\lambda_{ij}^{AB}}$`, que claramente não aparece neste modelo.

---
class: animated, fadeIn

##Modelo Saturado

No modelo saturado as `$N=I\times J$` contagens nas caselas ainda são assumidas como observações independentes de uma variável aleatória Poisson, mas nãos e assume independência entre `$A$` e `$B$`. O modelo é representado como

`$$\log(\mu_{ij})=\lambda+\lambda_i^A+\lambda_j^B+\lambda_{ij}^{AB}.$$`
Nesse modelo também se aplicam as restrições semelhantes às da ANOVA para contornar a superparametrização.

**Estimativas dos parâmetros:**

O odds ratio está diretamente relacionado com o termo de interação. Por exemplo, para uma tabela `$2\times 2$`:

`$$\begin{align*}
\log(OR)&=\log\left(\frac{\mu_{11}\mu_{22}}{\mu_{12}\mu_{21}}\right)\\
&=\lambda+\lambda_1^A+\lambda_1^B+\lambda_{11}^{AB}+\lambda+\lambda_2^A+\lambda_2^B+\lambda_{22}^{AB}-\lambda-\lambda_1^A-\lambda_2^B-\lambda_{12}^{AB}-\lambda-\lambda_2^A-\lambda_1^B-\lambda_{21}^{AB}\\
&=\lambda_{11}^{AB}+\lambda_{22}^{AB}-\lambda_{12}^{AB}-\lambda_{21}^{AB}.
\end{align*}$$`

---
class: animated, fadeIn

##Exemplo - I

Um estudo duplo-cego investigou o efeito terapeutico da vitamina C para tratar resfriados comuns. O estudo foi conduzido durante um período de duas semanas em uma amostra de 280 esquiadores franceses. Uma observação foi perdida por inconsistências.

| Tratamento || Resfriado | Sem resfriado || Total |
|:----------:||----------:|:-------------:||:-----:|
|            ||           |               ||       |
| Placebo    ||    31     |      109      ||  140  |
| Vit. C     ||    17     |      122      ||  139  |
|  Total     ||    48     |      231      ||  279  |

```r
dados<-matrix(c(31, 17, 109, 122), ncol=2, dimnames=list(Tratamento=c("Placebo", "VitC"), Resfriado=c("Com", "Sem")))
dados
```

```
##           Resfriado
## Tratamento Com Sem
##    Placebo  31 109
##    VitC     17 122
```

---
class: animated, fadeIn

##Exemplo - II

```r
library(vcd)
assocstats(dados)
```

```
##                     X^2 df P(> X^2)
## Likelihood Ratio 4.8717  1 0.027301
## Pearson          4.8114  1 0.028272
## 
## Phi-Coefficient   : 0.131 
## Contingency Coeff.: 0.13 
## Cramer's V        : 0.131
```

---
class: animated, fadeIn

##Exemplo - III

```r
OR<-oddsratio(dados, log=FALSE)
OR
```

```
##  odds ratios for Tratamento and Resfriado 
## 
## [1] 2.041015
```

```r
confint(OR)
```

```
##                         2.5 %  97.5 %
## Placebo:VitC/Com:Sem 1.070353 3.89193
```

>A chance de pegar resfriado é aproximadamente o dobro para aqueles que tomaram o placebo, quando comparados com os que tomaram vitamina C.

---
class: animated, fadeIn

##Exemplo - IV

###Modelo de Independência

```r
dados<-as.table(dados)
mod.indep<-loglin(dados, list(1,2), fit=TRUE, param=TRUE)
```

```
## 2 iterations: deviation 0
```

```r
mod.indep[-3]#eliminando a impressão do terceiro item da lista (graus de liberdade)
```

```
## $lrt
## [1] 4.871697
## 
## $pearson
## [1] 4.811413
## 
## $margin
## $margin[[1]]
## [1] "Tratamento"
## 
## $margin[[2]]
## [1] "Resfriado"
## 
## 
## $fit
##           Resfriado
## Tratamento       Com       Sem
##    Placebo  24.08602 115.91398
##    VitC     23.91398 115.08602
## 
## $param
## $param$`(Intercept)`
## [1] 3.963656
## 
## $param$Tratamento
##      Placebo         VitC 
##  0.003584245 -0.003584245 
## 
## $param$Resfriado
##        Com        Sem 
## -0.7856083  0.7856083
```

---
class: animated, fadeIn

##Exemplo - V

**log-odds estimados:**

`$$\begin{align*}
\log(\mu_{11})&=3.963656+0.003584245-0.7856083\\
\cdots\\
\log(\mu_{22})&=3.963656-0.003584245+0.7856083\\
\end{align*}$$`

**odds de resfriado:**

`$$\begin{align*}
\exp(\lambda_1^R-\lambda_2^R)&=\exp(-0.7856083-0.7856083)\\
&=\exp(-1.5712)\\
&=0.208
\end{align*}$$`

---
class: animated, fadeIn

##Exemplo - VI

###Modelo Saturado

```r
mod.sat<-loglin(dados, list(c(1,2)), fit=TRUE, param=TRUE)
```

```
## 2 iterations: deviation 0
```

```r
mod.sat
```

```
## $lrt
## [1] 0
## 
## $pearson
## [1] 0
## 
## $df
## [1] 0
## 
## $margin
## $margin[[1]]
## [1] "Tratamento" "Resfriado" 
## 
## 
## $fit
##           Resfriado
## Tratamento Com Sem
##    Placebo  31 109
##    VitC     17 122
## 
## $param
## $param$`(Intercept)`
## [1] 3.940642
## 
## $param$Tratamento
##    Placebo       VitC 
##  0.1220252 -0.1220252 
## 
## $param$Resfriado
##        Com        Sem 
## -0.8070421  0.8070421 
## 
## $param$Tratamento.Resfriado
##           Resfriado
## Tratamento        Com        Sem
##    Placebo  0.1783618 -0.1783618
##    VitC    -0.1783618  0.1783618
```

---
class: animated, fadeIn

##Exemplo - VII

**Contagens esperadas sob o modelo ajustado**

```
##           Resfriado
## Tratamento Com Sem
##    Placebo  31 109
##    VitC     17 122
```

**Parâmetros estimados (efeitos principais)**

```
## $`(Intercept)`
## [1] 3.940642
## 
## $Tratamento
##    Placebo       VitC 
##  0.1220252 -0.1220252 
## 
## $Resfriado
##        Com        Sem 
## -0.8070421  0.8070421
```

---
class: animated, fadeIn

##Exemplo - VIII

**Parâmetros estimados (Interação)**

```
## $Tratamento.Resfriado
##           Resfriado
## Tratamento        Com        Sem
##    Placebo  0.1783618 -0.1783618
##    VitC    -0.1783618  0.1783618
```

>Temos então:
`$$\begin{align*}
\log(\mu_{11})&=3.940642+0.1220252-0.8070421+0.1783618\\
\cdots\\
\log(\mu_{22})&=3.940642-0.1220252+0.8070421+0.1783618\\
\end{align*}$$`

---
class: animated, fadeIn

##Exemplo - IX

###Alternativa:

```r
dados1<-as.data.frame(dados)
dados1
```

```
##   Tratamento Resfriado Freq
## 1    Placebo       Com   31
## 2       VitC       Com   17
## 3    Placebo       Sem  109
## 4       VitC       Sem  122
```

```r
modelo.ind<-glm(dados1$Freq~dados1$Tratamento+dados1$Resfriado, family=poisson())
modelo.ind[1]
```

```
## $coefficients
##           (Intercept) dados1$TratamentoVitC   dados1$ResfriadoSem 
##           3.181631652          -0.007168489           1.571216700
```

---
class: animated, fadeIn

##Exemplo - X

```r
summary(modelo.ind)
```

```
## 
## Call:
## glm(formula = dados1$Freq ~ dados1$Tratamento + dados1$Resfriado, 
## family = poisson())
## 
## Coefficients:
## Estimate Std. Error z value Pr(>|z|) 
## (Intercept) 3.181632 0.156179 20.372 <2e-16 ***
## dados1$TratamentoVitC -0.007168 0.119738 -0.060 0.952 
## dados1$ResfriadoSem 1.571217 0.158626 9.905 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
## Null deviance: 135.4675 on 3 degrees of freedom
## Residual deviance: 4.8717 on 1 degrees of freedom
## AIC: 34.004
## 
## Number of Fisher Scoring iterations: 4
```

---
class: animated, fadeIn

##Exemplo - XI

###Modelo Saturado

```r
mod.sat<-loglin(dados, list(c(1,2)), fit=TRUE, param=TRUE)
```

```
## 2 iterations: deviation 0
```

```r
mod.sat
```

---
class: animated, fadeIn

##Exemplo - XII

**Contagens esperadas sob o modelo ajustado**

```
##           Resfriado
## Tratamento Com Sem
##    Placebo  31 109
##    VitC     17 122
```

**Parâmetros estimados (efeitos principais)**

```
## $`(Intercept)`
## [1] 3.940642
## 
## $Tratamento
##    Placebo       VitC 
##  0.1220252 -0.1220252 
## 
## $Resfriado
##        Com        Sem 
## -0.8070421  0.8070421
```

---
class: animated, fadeIn

##Exemplo - XIII

```r
anova(modelo.ind)
```

```
## Analysis of Deviance Table
## 
## Model: poisson, link: log
## 
## Response: dados1$Freq
## 
## Terms added sequentially (first to last)
## 
## 
##                   Df Deviance Resid. Df Resid. Dev
## NULL                                  3    135.468
## dados1$Tratamento  1    0.004         2    135.464
## dados1$Resfriado   1  130.592         1      4.872
```

---
class: animated, fadeIn

##Exemplo - XIV

```r
modelo.sat<-glm(dados1$Freq~dados1$Tratamento*dados1$Resfriado, family=poisson())
modelo.sat
```

```
## 
## Call:  glm(formula = dados1$Freq ~ dados1$Tratamento * dados1$Resfriado, 
##     family = poisson())
## 
## Coefficients:
##                               (Intercept)  
##                                    3.4340  
##                     dados1$TratamentoVitC  
##                                   -0.6008  
##                       dados1$ResfriadoSem  
##                                    1.2574  
## dados1$TratamentoVitC:dados1$ResfriadoSem  
##                                    0.7134  
## 
## Degrees of Freedom: 3 Total (i.e. Null);  0 Residual
## Null Deviance:	    135.5 
## Residual Deviance: -1.021e-14 	AIC: 31.13
```

---
class: animated, fadeIn

#Modelos Loglineares para tabelas de três vias

Seja `$\mu_{ijk}$` que representa a média dos níveis `$i$`, `$j$` e `$k$` das variáveis `$A$`, `$B$` e `$C$` respectivamente. O modelo é dado por

`$$\log(\mu_{ijk})=\lambda+\lambda_i^A+\lambda_j^B+\lambda_k^C+\lambda_{ij}^{AB}+\lambda_{ik}^{AC}+\lambda_{jk}^{BC}+\lambda_{ijk}^{ABC}.$$`

Há diversos modelos que podem ser testados e ajustados, são eles:

* Saturado

* Independência completa

* Independência conjunta

* Independência condicional

* Associação Homogênea

---
class: animated, fadeIn

##Modelo Saturado - I

É o modelo default e serve para testar a bondade do ajuste dos outros modelos.

Sejam `$N=I\times J\times K$` as contagens nas caselas, assumidas como observações independentes de uma variável aleatória Poisson. O modelo é dado por

`$$\log(\mu_{ijk})=\lambda+\lambda_i^A+\lambda_j^B+\lambda_k^C+\lambda_{ij}^{AB}+\lambda_{ik}^{AC}+\lambda_{jk}^{BC}+\lambda_{ijk}^{ABC}.$$`

Sujeito  a `$\lambda_I^A=\lambda_J^B=\lambda_K^C=\lambda_{Ij}^{AB}=\cdots=\lambda_{ijK}^{ABC}=0.$`

* `$\lambda$` representa o efeito global o grande média (na escala logarítmica) das contagens esperadas, garante que `$\sum_i\sum_j\sum_k \mu_{ijk}=n.$`

* `$\lambda_i^{A}$`, `$\lambda_j^{B}$` e `$\lambda_k^{C}$` representam os efeitos principais das variáveis `$A$`, `$B$` e `$C$`, ou desvios da média global.

* `$\lambda_{ij}^{AB}$`, `$\lambda_{ik}^{AC}$` e `$\lambda_{jk}^{BC}$` representam a interação/associação entre duas variáveis quando controlada a terceira (i.e., odds ratios condicionais, testes de associação parcial) e reflete o afastamento da independência.

* `$\lambda_{ijk}^{ABC}$` representa a associação/interação entre três variáveis e reflete o afastamento da independência.

---
class: animated, fadeIn

##Modelo Saturado - II

**Observações:**

* Se há termo de interação significante, normalmente não olhamos para os termos de ordem inferior e apenas interpretamos os termos de ordem superior.
 
 * O modelo saturado tem ajuste perfeito, com `$G^2=0$` e 0 graus de liberdade, pois o número de caselas é igual ao número de parâmetros únicos no modelo.

---
class: animated, fadeIn

##Exemplo 2 - I

Lembremos das informações sobre admissões nas pós-graduações nos seis maiores departamentos da U.C. Berkeley em 1973.

|  Departamento || Homens admitidos ||  Homens rejeitados  ||  Mulheres admitidas  ||  Mulheres rejeitadas  |
 |:-------------:||:----------------:||:-------------------:||:--------------------:||:---------------------:|
 |    A          ||     512          ||         313         ||      89              ||          19           |
 |    B          ||     353          ||         207         ||      17              ||          8            |
 |    C          ||     120          ||         205         ||      202             ||        391            |
 |    D          ||     139          ||         279         ||      131             ||        244            |
 |    E          ||      53          ||         138         ||       94             ||        299            |
 |    F          ||      22          ||         351         ||       24             ||        317            |

Denotemos por `$D=$` Departamento, `$S=$` Sexo e `$A=$` Admissão (admitido ou rejeitado). Seja `$Y$` a frequência observada em uma casela particular da tabela de três vias.

---
class: animated, fadeIn

##Exemplo 2 - II

```r
UCBAdmissions
```

```
## , , Dept = A
## 
##           Gender
## Admit      Male Female
##   Admitted  512     89
##   Rejected  313     19
## 
## , , Dept = B
## 
##           Gender
## Admit      Male Female
##   Admitted  353     17
##   Rejected  207      8
## 
## , , Dept = C
## 
##           Gender
## Admit      Male Female
##   Admitted  120    202
##   Rejected  205    391
## 
## , , Dept = D
## 
##           Gender
## Admit      Male Female
##   Admitted  138    131
##   Rejected  279    244
## 
## , , Dept = E
## 
##           Gender
## Admit      Male Female
##   Admitted   53     94
##   Rejected  138    299
## 
## , , Dept = F
## 
##           Gender
## Admit      Male Female
##   Admitted   22     24
##   Rejected  351    317
```

---
class: animated, fadeIn

##Exemplo 2 - III

```r
berkeley<-as.data.frame(UCBAdmissions)
berkeley$Gender = relevel(berkeley$Gender, ref='Female')
berkeley$Dept = relevel(berkeley$Dept, ref='F')
berkeley.sat = glm(Freq~Admit*Gender*Dept, family=poisson(), data=berkeley)
summary(berkeley.sat)
```

```
## 
## Call:
## glm(formula = Freq ~ Admit * Gender * Dept, family = poisson(), 
## data = berkeley)
## 
## Coefficients:
## Estimate Std. Error z value Pr(>|z|) 
## (Intercept) 3.17805 0.20412 15.569 < 2e-16 ***
## AdmitRejected 2.58085 0.21171 12.190 < 2e-16 ***
## GenderMale -0.08701 0.29516 -0.295 0.7682 
## DeptA 1.31058 0.23001 5.698 1.21e-08 ***
## DeptB -0.34484 0.31700 -1.088 0.2767 
## DeptC 2.13021 0.21591 9.866 < 2e-16 ***
## DeptD 1.69714 0.22204 7.644 2.11e-14 ***
## DeptE 1.36524 0.22870 5.969 2.38e-09 ***
## AdmitRejected:GenderMale 0.18890 0.30516 0.619 0.5359 
## AdmitRejected:DeptA -4.12505 0.32968 -12.512 < 2e-16 ***
## AdmitRejected:DeptB -3.33462 0.47817 -6.974 3.09e-12 ***
## AdmitRejected:DeptC -1.92041 0.22876 -8.395 < 2e-16 ***
## AdmitRejected:DeptD -1.95888 0.23781 -8.237 < 2e-16 ***
## AdmitRejected:DeptE -1.42370 0.24250 -5.871 4.33e-09 ***
## GenderMale:DeptA 1.83670 0.31672 5.799 6.66e-09 ***
## GenderMale:DeptB 3.12027 0.38572 8.090 5.99e-16 ***
## GenderMale:DeptC -0.43376 0.31687 -1.369 0.1710 
## GenderMale:DeptD 0.13907 0.31938 0.435 0.6632 
## GenderMale:DeptE -0.48599 0.34151 -1.423 0.1547 
## AdmitRejected:GenderMale:DeptA 0.86318 0.40267 2.144 0.0321 * 
## AdmitRejected:GenderMale:DeptB 0.03113 0.53349 0.058 0.9535 
## AdmitRejected:GenderMale:DeptC -0.31382 0.33741 -0.930 0.3523 
## AdmitRejected:GenderMale:DeptD -0.10691 0.34013 -0.314 0.7533 
## AdmitRejected:GenderMale:DeptE -0.38908 0.36500 -1.066 0.2864 
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
## Null deviance: 2.6501e+03 on 23 degrees of freedom
## Residual deviance: 1.2612e-13 on 0 degrees of freedom
## AIC: 207.06
## 
## Number of Fisher Scoring iterations: 3
```

---
class: animated, fadeIn

##Exemplo 2 - IV

**Observações:**

* Os efeitos principais de `$D$`, `$S$` e `$A$` são difíceis de interpretar e não muito úteis, devido a que há associações de duas e três vias significantes

* Por exemplo, o coeficiente estimado da associação SA é `$0.1889$`. Dessa forma, o odds ratio estimado para o Departamento `$F$` (departamento de referência) é `$\exp(0.1889)=1.208$`. A categoria de referência para S é "Mulher" e para A é "admitido". A tabela S `$\times$` A no Departamento F, i.e., a tabela parcial é dada por

|  Dep. F  || Rejeitado  |  Admitido  |
|:--------:||:----------:|:----------:|
|  Homem   ||    351     |     22     |
|  Mulher  ||    317     |     24     |
 
De onde o odds ratio estimado é
 
 `$$OR=\frac{351\times 24}{317\times 22}=1.208.$$`

---
class: animated, fadeIn

##Exemplo 2 - V

* A estatística de Wald para o coeficiente SA é dada por
`$$z=0.1889/0.3052,$$`
que corresponde à estatística `$\chi^2$` `$(0.62^2=0.38)$` com valor `$p=0.5359$` e indica que  o odds ratio SA para o departamento F não é significativamente diferente de 1 (log-odds ratio não é significativamente diferente de 0).

* Para obter o odds ratio SA para qualquer outro departamento, devemos combinar os coeficientes SA com um dos coeficientes DSA. Por exemplo, os odds ratio SA para o departamento A é
`$$\exp(0.1889+0.8632)=2.864.$$`

* A estatística `$z$` indica que o odds ratio SA para o departamento A é significativamente diferente do odds ratio SA no departamento F.

* O modelo saturado se torna mais complicado na medida em que o número de variáves aumenta.

---
class: animated, fadeIn

##Modelo de Independência Completa

* É o modelo mais restritivo. Todas as variáveis são assumidas conjuntamente independentes, independentemente de qualquer condicionamento.

* Assume-se que as `$N=J\times J\times K$` contagens nas caselas são observações independentes de uma variável aleatória Poisson.

* Assume-se também que não há interações parciais `$\lambda_{ij}^{AB}=\lambda_{ik}^{AC}=\lambda_{jk}^{BC}$`, para todo `$i,j,k$`; nem interações triplas `$\lambda_{ijk}^{ABC}=0$`, para todo `$i,j,k.$`

* Restrições devem ser impostas para evitar a superparametrização.

* A estrutura do modelo é dada por

`$$\log(\mu_{ijk})=\lambda + \lambda_i^A + \lambda_j^B + \lambda_k^C.$$`

---
class: animated, fadeIn

##Exemplo 2 - VI

```r
berkeley.ind = glm(Freq~Admit + Gender + Dept, family=poisson(), data=berkeley)
summary(berkeley.ind)
```

```
## 
## Call:
## glm(formula = Freq ~ Admit + Gender + Dept, family = poisson(), 
## data = berkeley)
## 
## Coefficients:
## Estimate Std. Error z value Pr(>|z|) 
## (Intercept) 4.72072 0.04553 103.673 < 2e-16 ***
## AdmitRejected 0.45674 0.03051 14.972 < 2e-16 ***
## GenderMale 0.38287 0.03027 12.647 < 2e-16 ***
## DeptA 0.26752 0.04972 5.380 7.44e-08 ***
## DeptB -0.19927 0.05577 -3.573 0.000352 ***
## DeptC 0.25131 0.04990 5.036 4.74e-07 ***
## DeptD 0.10368 0.05161 2.009 0.044533 * 
## DeptE -0.20098 0.05579 -3.602 0.000315 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
## Null deviance: 2650.1 on 23 degrees of freedom
## Residual deviance: 2097.7 on 16 degrees of freedom
## AIC: 2272.7
## 
## Number of Fisher Scoring iterations: 5
```

---
class: animated, fadeIn

##Exemplo 2 - VII

```r
berkeley.ind$dev
```

```
## [1] 2097.671
```

O Deviance do modelo, com 16 graud de liberdade indica que o modelo não ajusta bem. Se o modelo ajustasse adequadamente, a razão "Deviance/gl" seria próxima de 1.

---
class: animated, fadeIn

##Independência Conjunta

* Com três variáveis, há três formas de obter independência conjunta. O pressuposto é que uma variável é independente das outras duas, mas essas duas podem ter qualquer associação arbitrária.

* Se A é cojuntamente independente de B e C, denotado por (A, BC), então A é independente de B e de C, e podemos fatorar a distribuição conjunta como o produto da distribuiçao marginal de A e da conjunta de (B, C).

* Sem perda de generalidade, para o caso (A, BC), assume-se que as `$N=I\times J\times K$` contagens são observações de uma variável aleatória Poisson.

* Assume-se também que não há interações envolvendo A: `$\lambda_{ij}^{AB}=\lambda_{ik}^{AC}=0$`, para todo `$i, j, k$`; nem ha interação tripla `$\lambda_{ijk}^{ABC}$`, para todo `$i,j,k.$`

* A estrutura do modelo é dada por

`$$\log(\mu_{ijk})=\lambda + \lambda_i^A + \lambda_j^B + \lambda_k^C + \lambda_{jk}^{BC}.$$`

---
class: animated, fadeIn

##Exemplo 2 - VIII

Consideremos o modelo onde a admissão é conjuntamente independente de departamento e sexo, denotado por (A, DS):

```r
berkeley.conj = glm(Freq~Admit + Gender + Dept + Gender*Dept, family=poisson(), data=berkeley)
berkeley.conj$coef
```

```
##      (Intercept)    AdmitRejected       GenderMale            DeptA 
##       4.88451279       0.45673941       0.08969594      -1.14975125 
##            DeptB            DeptC            DeptD            DeptE 
##      -2.61300665       0.55331192       0.09504355       0.14192713 
## GenderMale:DeptA GenderMale:DeptB GenderMale:DeptC GenderMale:DeptD 
##       1.94355622       3.01936502      -0.69106516       0.01646425 
## GenderMale:DeptE 
##      -0.81123213
```

```r
berkeley.conj$dev
```

```
## [1] 877.0564
```

---
class: animated, fadeIn

##Exemplo 2 - IX

* O modelo implica que a associação entre D e S não depende do nível da variável A. Dessa forma, a associação entre o departamento e sexo é independente da decisão de aceitar ou rejeitar.

* O primeiro coeficiente estimado para a associação DS indica que o odds ratio estimado entre sexo e departamento (especificamente A vs F) é  `$\exp(1.9436)=6,98$` com IC de 95% dado por 
`$$\left[\exp(1.695);\exp(2.192)\right]=(5,45; 8,95).$$`

* A estatística de bondade do ajuste (Deviance=877.06) indica que o modelo não ajusta bem, pois o valor "Deviance/gl" é muito maior que 1.

---
class: animated, fadeIn

##Modelo de independência Condicional

* Para três variáveis, há três formas de obter independência condicional. Assume-se que duas variáveis são independentes, dada a terceira variável. Por exemplo, se A e B são condicionalmente independentes, dado C, denotado por (AC, BC), então a distribuição  de (AB), dado C, pode ser fatorada como o produto das duas marginais, dado C.

* Sem perda de generalidade, para o caso (AC, BC), assume-se que as `$N=I\times J\times K$` contagens são observações de uma variável aleatória Poisson.

* Assume-se também que não há interação parcial entre A e B: `$\lambda_{ij}^{AB}=0$`, para todo `$i, j$`; nem ha interação de três vias `$\lambda_{ijk}^{ABC}$`, para todo `$i,j,k.$`

* A estrutura do modelo é dada por

`$$\log(\mu_{ijk})=\lambda + \lambda_i^A + \lambda_j^B + \lambda_k^C + \lambda_{ik}^{AC} + \lambda_{jk}^{BC}.$$`

---
class: animated, fadeIn

##Exemplo 2 - X

Consideremos o modelo onde a admissão e departamento são conjuntamente independentes, dado o sexo, denotado por (AS, DS):

```r
berkeley.cond = glm(Freq~Admit + Gender + Dept + Admit*Gender + Dept*Gender, family=poisson(), data=berkeley)
berkeley.cond$coef
```

```
##              (Intercept)            AdmitRejected               GenderMale 
##               4.63964796               0.83048640               0.47267109 
##                    DeptA                    DeptB                    DeptC 
##              -1.14975125              -2.61300665               0.55331192 
##                    DeptD                    DeptE AdmitRejected:GenderMale 
##               0.09504355               0.14192713              -0.61035238 
##         GenderMale:DeptA         GenderMale:DeptB         GenderMale:DeptC 
##               1.94355622               3.01936502              -0.69106516 
##         GenderMale:DeptD         GenderMale:DeptE 
##               0.01646425              -0.81123213
```

```r
berkeley.cond$dev
```

```
## [1] 783.607
```

---
class: animated, fadeIn

##Associação Homogênea

* O modelo de associações homogêneas também é conhecido como **Modelo sem interações triplas**. Denotado por (AB, AC, BC), a única restrição que esse modelo impõe é que nenhuma associação condicional por pares depende  do valor da terceira variável. Por exemplo, o odds ratio condicional entre A e B, dado que C está no primeiro nível, deve ser o mesmo odds ratio condicional entre A e B, dado C no segundo nível, e assim por diante.

* Assume-se que as `$N=I\times J\times K$` contagens são observações de uma variável aleatória Poisson.

* Não há interação tripla `$\lambda_{ijk}^{ABC}$`, para todo `$i,j,k.$`

* O modelo é dado por

`$$\log(\mu_{ijk})=\lambda + \lambda_i^A + \lambda_j^B + \lambda_k^C + \lambda_{ij}^{AB} + \lambda_{jk}^{BC} + \lambda_{ik}^{AC}.$$`

---
##Exemplo 2 - XI

```r
berkeley.homog = glm(Freq~(Admit + Gender + Dept)^2, family=poisson(), data=berkeley)
berkeley.homog$coef
```

```
##              (Intercept)            AdmitRejected               GenderMale 
##              3.137357882              2.624558572             -0.003730812 
##                    DeptA                    DeptB                    DeptC 
##              1.135552319             -0.342488552              2.222782088 
##                    DeptD                    DeptE AdmitRejected:GenderMale 
##              1.743871669              1.480917951              0.099870088 
##      AdmitRejected:DeptA      AdmitRejected:DeptB      AdmitRejected:DeptC 
##             -3.306480056             -3.263082125             -2.043882034 
##      AdmitRejected:DeptD      AdmitRejected:DeptE         GenderMale:DeptA 
##             -2.011873587             -1.567174318              2.002319157 
##         GenderMale:DeptB         GenderMale:DeptC         GenderMale:DeptD 
##              3.077139541             -0.662813561              0.043994838 
##         GenderMale:DeptE 
##             -0.792866728
```

```r
berkeley.homog$dev
```

```
## [1] 20.20428
```

---
##Exemplo 2 - XII

* O modelo implica que a associação condicional entre sexo  e departamento não depende o valor fixado de admissão. A associação condicional entre sexo e admissão não dependem do departamento fixado, e a associação entre departamento e admissão não dependem do sexo.

* O interesse está nas interações de ordem superior. Por exemplo, o primeiro coeficiente do cunjunto "AdmitRejected:GenderMale", é o log odds ratio condicional entre sexo e admissão. Nesse caso, para um departamento fixado, o odds de um homem ser rejeitado é `$\exp(0.0999)=1.1506$` vezes o odds de uma mulher ser rejeitada.

* Embora nessa interpretação o departamento esteja fixado (dessa forma a comparação é feita entre indivíduos no mesmo departamento), não interessa em qual departamento estamos focando, pois todos conduzem ao mesmo resultado sob esse modelo.

* O modelo ainda não parece ajustar adequadamente (Deviance = 20.204, com 5 gl). Mas comparado com os anteriores, parece ser melhor.

---
class: animated, fadeIn

##Seleção do Modelo

* Considerando que os modelos aqui estudados são hierárquicos, onde cada um é caso especial de outro, usamos a estatística de razão de verossimilhanças para medir a redução no ajuste do modelo menor (hipótese nula), com respeito ao modelo maior (hipótese alternativa). Os graus de liberdade desses testes são a diferença entre os números de parâmetros envolvidos nos dois modelos.

---
class: animated, fadeIn

##Exemplo 2 - XIII

Resumo de todos os possíveis modelos

| Modelo | gl | `$G^2$` | valor p |
|:-----------:|:----------:|:----------:|:----------:|
| (D, S, A) | 16 | 2097.671 | <0.01 |
| (DS, A) | 11 | 877.056 | <0.01 |
| (D, SA) | 15 | 2004.222 | <0.01 |
| (DA, S) | 11 | 1242.350 | <0.01 |
| (DS, SA) | 10 | 783.607 | <0.01 |
| (DS, DA) | 6 | 21.736 | <0.01 |
| (DA, SA) | 10 | 1148.901 | <0.01 |
| (DS, DA, SA) | 5 | 20.204 | <0.01 |
| (DSA) | 0 | 0.00 | |

---
class: animated, fadeIn

##Exemplo 2 - XIV

**Observações:**
 * O modelo saturado seria o melhor, considerando que o qualquer modelo reduzido tem ajuste significativamente pior.

* Por exemplo, para testar (DS, DA), que assume que sexo e admissão são condicionalmente independentes, dado o departamento, as hipóteses seriam
`$$\begin{align*}
H_0&: (DS, DA)\\
H_1&: (DS, DA, SA)
\end{align*}$$`

A estatística de teste pode ser calculada diretamente das estatísticas deviance dos modelos ajustados. Neste caso,
`$$G^2=783.607-20.204=763.4$$`
Com distribuição `$\chi^2$` com `$10-5=5$` graus de liberdade, que é altamente significante.

Caso o modelo de independência condicional não tivesse sido rejeitado, esse poderia ser colocado na hipótese alternativa para testar o modelo reducido de independência conjunta (DS, A), e assim por diante.

---
class: animated, hide-logo, bounceInDown
## Política de proteção aos direitos autorais

> O conteúdo disponível consiste em material protegido pela legislação brasileira, sendo certo que, por ser
o detentor dos direitos sobre o conteúdo disponível na plataforma, o **LECON** e o **NEAEST** detém direito
exclusivo de usar, fruir e dispor de sua obra, conforme Artigo 5o, inciso XXVII, da Constituição Federal
e os Artigos 7o e 28o, da Lei 9.610/98.
A divulgação e/ou veiculação do conteúdo em sites diferentes à plataforma e sem a devida autorização do
**LECON** e o **NEAEST**, pode configurar violação de direito autoral, nos termos da Lei 9.610/98, inclusive podendo
caracterizar conduta criminosa, conforme Artigo 184o, §1o a 3o, do Código Penal.
É considerada como contrafação a reprodução não autorizada, integral ou parcial, de todo e qualquer
conteúdo disponível na plataforma.

.pull-left[
<img src="images/logo_lecon.png" width="50%" style="display: block; margin: auto;" />
]
.pull-right[
<img src="images/logo_neaest.png" width="50%" style="display: block; margin: auto;" />
]
 
.center[
[https://lecon.ufes.br](https://lecon.ufes.br/)]

Material elaborado pela equipe LECON/NEAEST: 
Alessandro J. Q. Sarnaglia, Bartolomeu Zamprogno, Fabio A. Fajardo, Luciana G. de Godoi 
e Nátaly A. Jiménez.