1 Comparando múltiplas médias

É recomendado que as variáveis de grupo estejam no formato factor.

2 ANOVA - Análise de Variância

É usada para comparar a média de múltiplos grupos (três ou mais).

3 One-way ANOVA (um fator)

É uma extensão do teste t para amostras independentes em situações onde há mais de dois grupos.
A ideia por trás da ANOVA é: se a variação média entre grupos é grande o suficiente se comparada à variação média dentro dos grupos, então você pode concluir que pelo menos um grupo possui média diferente dos outros.

A matemática por trás da ANOVA é a seguinte:
- Computa-se a variância dentro do grupo (variância residual).
- Computa a variância entre as médias dos grupos.
- Computa-se a estatística F como a razão entre variancia.entre.grupos/variancia.dentro.dos.grupos
- F < 1 indica que não há diferença significativa entre as médias dos grupos a serem comparados.
- F ≥ 1, implica que a variação entre grupos é muito maior que a variação de observações individuais em cada grupo.

3.1 Suposições da ANOVA de um fator

Independência das observações
Sem ouliers extremos
Normalidade
Homogeneidade de variâncias

Caso as suposições não sejam satisfeitas, há uma alternativa não paramétrica, o teste de Kruskal-Wallis.

3.2 Preparação do ambiente

library(tidyverse)
library(ggpubr)
library(rstatix)

3.3 Base de dados

data("PlantGrowth")
PlantGrowth %>% 
  sample_n(5)

glimpse(PlantGrowth)

## Rows: 30
## Columns: 2
## $ weight <dbl> 4.17, 5.58, 5.18, 6.11, 4.50, 4.61, 5.17, 4.53, 5.33, 5.14, 4.8…
## $ group  <fct> ctrl, ctrl, ctrl, ctrl, ctrl, ctrl, ctrl, ctrl, ctrl, ctrl, trt…

table(PlantGrowth$group)

## 
## ctrl trt1 trt2 
##   10   10   10

class(PlantGrowth$group)

## [1] "factor"

Queremos saber se o crescimento médio das plantas, medido através de seu peso, apresenta diferença significativa para três condições diferentes.

3.4 Estatística Descritiva

PlantGrowth %>% 
  dplyr::group_by(group) %>% 
  get_summary_stats(weight, type = "mean_sd")

PlantGrowth %>% 
  dplyr::group_by(group) %>% 
  get_summary_stats(weight, type = "full")

3.5 Visualização inicial

boxplot_plantas <- PlantGrowth %>% 
  ggboxplot(x = "group", y = "weight")
boxplot_plantas

3.6 Checando suposições

3.6.1 Independência

Sabe-se que as amostras são independentes, uma vez que cada planta só recebe um tratamento.

3.6.2 Outliers

PlantGrowth %>% 
  dplyr::group_by(group) %>% 
  identify_outliers(weight)

Não há outiliers extremos.

3.6.3 Normalidade

Há duas maneiras de testar normalidade para ANOVA:

Analisando os resíduos da ANOVA: checa a normalidade para todos os grupos juntos. Esta abordagem é mais simples e muito útil quando se tem muitos grupos ou se há poucos dados por grupo.
Checando a normalidade para cada grupo de maneira separada: Esta abordagem deve ser utilizada quando há poucos grupos e muitos dados por grupo.

No nosso caso temos poucos grupos (3), e também poucos dados. Portanto, vamos ver as duas abordagens.

3.6.3.1 Checando normalidade pelos resíduos do modelo

Construindo um modelo linear:

modelo_linear <- lm(weight ~ group, data = PlantGrowth)

Criando um qqplot para os resíduos do modelo:

ggqqplot(residuals(modelo_linear))

Fazendo o teste de Shapiro-Wilk nos resíduos do modelo:

shapiro_test(residuals(modelo_linear))

No qq plot, apenas dois pontos saem do intervalo de confiança, mas a grande maioria cai na reta de referência, podemos assumir normalidade, o que é suportado pelo teste de Shapiro-Wilk (p = 0.44).

3.6.3.2 Checando normalidade por grupos

PlantGrowth %>% 
  ggqqplot("weight", facet.by = "group")

PlantGrowth %>%
  dplyr::group_by(group) %>% 
  shapiro_test(weight)

O peso segue uma distribuição normal em todos os grupos (p > 0.05, teste de Shapiro-Wilk). O QQ plot mostra que os pontos caem aproximadamente na linha de referência e intervalo de confiança.

Se sua amostra é maior que 50, o QQ plot é preferível, dado que o teste de Shapiro-Wilk é muito sensível, mesmo à um mínimo desvio na normalidade.

3.6.4 Homogeneidade de variância

O gráfico de resíduos vs. ajuste pode ser checado para a homogeneidade de variâncias

plot(modelo_linear, 1)

Não há nenhma relação evidente entre os resíduos e os valores ajustados (as médias de cada grupo), o que é bom. Então, podemos assumir homogeneidade de variâncias.

Também é possível utilizar o teste de Levene:

PlantGrowth %>% 
  levene_test(weight ~ group)

p = 0.341, não há diferença significativa entre ad variâncias para os grupos.

No caso em que a suposição de homogeneidade de variâncias não for satisfeita, usas-se o teste Welch para ANOVA de uma via (welch_anova_test() do pacote rstatix). Este teste não requer que as variâncias sejam homogêneas.

3.7 Computando o teste

resultado_anova <- PlantGrowth %>% 
  anova_test(weight ~ group)

## Coefficient covariances computed by hccm()

resultado_anova

ges: eta quadrado generalizado (tamanho do efeito). O valor significa que 26.4% da variabilidade no peso pode ser explicado pelo tratamento adotado.
F: indica a estatística de teste com 2 e 27 graus de liberdade.
p: p-valor

Pelo resultado da ANOVA, podemos ver que há diferença estatisticamente significativa entre a média de peso das plantas para diferentes tratamentos (p = 0.016, F(2, 27) = 4.85, eta2[g] = 0.26).

3.8 Teste post-hoc

Sabemos que há diferença entre os grupos, agora vamos realizar comparações dois a dois para ter mais detalhes sobre as diferenças entre grupos.
Vamos usar o teste post-hoc de Tukey.

comp_dois_dois <- PlantGrowth %>% 
  tukey_hsd(weight ~ group)
comp_dois_dois

O resultado nos mostra:

estimate: estimativa da diferença entre as médias dos dois grupos.
conf.low, conf.high: limites inferior e superior do intervalo de confiança de 95% (padrão).
p.adj: p-valor depois do ajuste pelas comparações dois a dois.

Podemos ver que apenas tr1 e tr2 apresentam diferença significativa entre si (p = 0.012).

3.9 Reportando o resultado

O teste ANOVA de um fator foi executado para avaliar se o crescimentos de plantas foi diferente para 3 grupos de tratamento: ctr (n = 10), trt1 (n = 10) e trt2 (n = 10).

Os dados são apresentados como média ± desvio-padrão. O cresciemnto das plantas foi estatisticamente diferente entre os grupos de tratamento considerados (F(2, 27) = 4.85, p = 0.016, eta^2 generalizado = 0.26).

O crescimento da planta foi menor com o tratamento 1 trt1 (4.66 ± 0.79) quando comparado com o controle ctr (5.03 ± 0.58). O crescimento foi maior com o tratamento 2 trt2 (5.53 ± 0.44) quando se compara tanto com o trt1 quanto com o ctr.

O teste post-hoc de Tukey revelou que o aumento no crescimento das plantas comparando os tratamentos trt1 e trt2 (0.87, IC 95%: (0.17-1.56)) foi estatisticamente significativo (p = 0.012), mas nenhum outro grupo demonstrou diferenças estatisticamente significativas.

Gráfico final

comp_dois_dois <- comp_dois_dois %>% 
  add_xy_position(x = "group") 

boxplot_plantas +
  stat_pvalue_manual(comp_dois_dois, hide.ns = TRUE) +
  labs(subtitle = get_test_label(resultado_anova, detailed = TRUE),
       caption = get_pwc_label(comp_dois_dois),
       x = "Tratamento", y = "Peso (kg)")

ggsave("box_plot_anova.png",
       width = 10, height = 10,
       units = "cm")

Aula 6 - ANOVA

by Edneide Ramalho

31 Julho, 2021