Como conduzir uma ancova em r

Este tutorial fornece um exemplo de como realizar uma ANCOVA em R.

Exemplo: ANCOVA em R

Realizaremos uma ANCOVA para testar se a técnica de estudo tem ou não impacto nos resultados dos exames usando as seguintes variáveis:

• Estudo técnico : A variável independente que queremos analisar
• Nota atual do aluno : A covariável que desejamos levar em consideração
• Pontuação de revisão : as variáveis de resposta que queremos analisar

O conjunto de dados a seguir contém informações sobre 90 alunos divididos aleatoriamente em três grupos de 30.

O conjunto de dados mostra a técnica de estudo que cada aluno usou (A, B ou C) , sua nota atual na turma quando começou a usar a técnica de estudo e a nota que recebeu no exame depois de usar a técnica de estudo por um mês para se preparar para o exame. exame:

 #make this example reproducible
set.seed(10)

#create dataset
data <- data.frame(technique = rep(c("A", "B", "C"), each = 30),
                   current_grade = runif(90, 65, 95),
                   exam = c(runif(30, 80, 95), runif(30, 70, 95), runif(30, 70, 90)))

#view first six lines of dataset
head(data)

# technical current_grade exam
#1 A 80.22435 87.32759
#2 A 74.20306 90.67114
#3 A 77.80723 88.87902
#4 A 85.79306 87.75735
#5 A 67.55408 85.72442
#6 A 71.76310 92.52167

Etapa 1: explore os dados

Antes de ajustar o modelo ANCOVA, primeiro precisamos explorar os dados para melhor compreendê-los e verificar se não existem valores extremos que possam distorcer os resultados.

Primeiro, podemos exibir um resumo de cada variável no conjunto de dados:

 summary(data)

# technical current_grade exam      
#A:30 Min. :65.43 Min. :71.17  
# B:30 1st Qu.:71.79 1st Qu.:77.27  
# C:30 Median:77.84 Median:84.69  
# Mean:78.15 Mean:83.38  
# 3rd Qu.:83.65 3rd Qu.:89.22  
# Max. :93.84 Max. :94.76  

Podemos observar que cada valor da técnica de estudo ( A, B e C) aparece 30 vezes nos dados.

Também podemos ver como as pontuações atuais dos alunos estavam distribuídas no início do estudo. A nota mínima da turma foi 65,43, a máxima foi 93,84 e a média foi 78,15.

Da mesma forma, podemos verificar que a nota mínima obtida no exame foi 71,17, a nota máxima foi 94,76 e a média foi 83,38.

Em seguida, podemos usar o pacote dplyr para encontrar facilmente a média e o desvio padrão das notas atuais e dos resultados dos exames para cada técnica de estudo:

 #load dplyr
library(dplyr)

data %>%
  group_by (technical) %>%
  summarize (mean_grade = mean(current_grade),
            sd_grade = sd(current_grade),
            mean_exam = mean(exam),
            sd_exam = sd(exam))

# A tibble: 3 x 5
# technique mean_grade sd_grade mean_exam sd_exam                      
#1 A 79.0 7.00 88.5 3.88
#2 B 78.5 8.33 81.8 7.62
#3 C 76.9 8.24 79.9 5.71

Podemos ver que a média e os desvios padrão da nota atual dos alunos que utilizam cada técnica de estudo são aproximadamente semelhantes.

Podemos observar também que a nota média do exame é significativamente maior para os alunos que utilizaram a técnica de estudo A em comparação às técnicas B e C.

Também podemos visualizar a distribuição dos resultados dos exames com base na técnica de estudo usando boxplots :

 boxplot(exam ~ technique,
data = data,
main = "Exam Score by Studying Technique",
xlab = "Studying Technique",
ylab = "Exam Score",
col = "steelblue",
border = "black"
)

Da mesma forma, também podemos usar boxplots para visualizar a distribuição das notas atuais com base na técnica de estudo:

 boxplot(current_grade ~ technical,
data = data,
main = "Current Grade by Studying Technique",
xlab = "Studying Technique",
ylab = "Current Grade",
col = "steelblue",
border = "black"
)

Passo 2: Verifique as suposições do modelo

Depois de realizarmos alguma exploração básica de dados e estarmos familiarizados com os dados, precisamos verificar se as seguintes premissas para ANCOVA são atendidas:

• A covariável e o tratamento são independentes – é necessário verificar se a covariável ( nota atual) e o tratamento (técnica de estudo) são independentes entre si, pois adicionar um termo de covariável no modelo só tem sentido se a covariável e o tratamento atua independentemente da variável resposta ( exame ).
• Homogeneidade de variância – precisamos verificar se as variâncias entre os grupos são iguais

Para verificar se a covariável e o tratamento são independentes, podemos realizar uma ANOVA usando a nota atual como variável resposta e a técnica de estudo como variável preditora:

 #fit anova model
anova_model <- aov(current_grade ~ technique, data = data)
#view summary of anova model
summary(anova_model)

# Df Sum Sq Mean Sq F value Pr(>F)
#technical 2 74 37.21 0.599 0.552
#Residuals 87 5406 62.14    

O valor p é maior que 0,05, portanto a covariável ( nota atual) e o tratamento ( técnica de estudo ) parecem independentes.

Então, para verificar se há homogeneidade de variância entre os grupos, podemos realizar o teste de Levene:

 #load car library to conduct Levene's Test
libary(car)

#conduct Levene's Test
leveneTest(technical exam, data = data)

#Levene's Test for Homogeneity of Variance (center = median)
# Df F value Pr(>F)    
#group 2 9.4324 0.0001961 ***
#87   

O valor p do teste é igual a 0,0001961, o que indica que as variâncias entre os grupos não são iguais. Embora pudéssemos tentar uma transformação nos dados para corrigir esse problema, não nos preocuparemos muito com diferenças de variância neste momento.

Passo 3: Ajustar o modelo ANCOVA

A seguir, ajustaremos o modelo ANCOVA usando a pontuação do exame como variável de resposta, a técnica de estudo como variável preditora (ou “tratamento”) e a nota atual como covariável.

Usaremos a função Anova() no pacote porque para fazer isso, apenas para poder especificar que queremos usar a soma dos quadrados do tipo III para o modelo, pois a soma dos quadrados do tipo I depende da ordem em que os preditores são inseridos no modelo:

 #load car library
library(car)

#fit ANCOVA model
ancova_model <- aov(exam ~ technique + current_grade, data = data)

#view summary of model
Anova(ancova_model, type="III") 

#Answer: exam
# Sum Sq Df F value Pr(>F)    
#(Intercept) 7161.2 1 201.4621 < 2.2e-16 ***
#technical 1242.9 2 17.4830 4.255e-07 ***
#current_grade 12.3 1 0.3467 0.5576    
#Residuals 3057.0 86         

Podemos observar que o valor p da técnica é extremamente baixo, indicando que a técnica de estudo tem um efeito estatisticamente significativo nas notas dos exames, mesmo após o controle da nota atual.

Etapa 4: teste post-hoc

Embora os resultados da ANCOVA nos tenham dito que a técnica de estudo teve um efeito estatisticamente significativo nas pontuações dos exames, precisamos realizar testes post hoc para determinar quais técnicas de estudo diferem umas das outras.

Para fazer isso, podemos usar a função glht() no pacote multcomp em R para realizar o teste de Tukey para comparações múltiplas:

 #load the multcomp library
library(multicomp)

#fit the ANCOVA model
ancova_model <- aov(exam ~ technique + current_grade, data = data)

#define the post hoc comparisons to make
postHocs <- glht(ancova_model, linfct = mcp(technique = "Tukey"))

#view a summary of the post hoc comparisons
summary(postHocs)

#Multiple Comparisons of Means: Tukey Contrasts
#
#Fit: aov(formula = exam ~ technique + current_grade, data = data)
#
#Linear Assumptions:
#Estimate Std. Error t value Pr(>|t|)    
#B - A == 0 -6.711 1.540 -4.358 0.000109 ***
#C - A == 0 -8.736 1.549 -5.640 < 1e-04 ***
#C - B == 0 -2.025 1.545 -1.311 0.393089    

#view the confidence intervals associated with the multiple comparisons
confint(postHocs)

# Simultaneous Confidence Intervals
#
#Multiple Comparisons of Means: Tukey Contrasts
#
#Fit: aov(formula = exam ~ technique + current_grade, data = data)
#
#Quantile = 2.3845
#95% family-wise confidence level
#
#Linear Assumptions:
# Estimate lwr upr     
#B - A == 0 -6.7112 -10.3832 -3.0392
#C - A == 0 -8.7364 -12.4302 -5.0426
#C - B == 0 -2.0252 -5.7091 1.6588

Pelo resultado, podemos perceber que há uma diferença estatisticamente significativa (em α = 0,05) nos resultados dos exames entre estudar a técnica A e estudar a técnica B (valor p: 0,000109), bem como entre a técnica A e a técnica C (valor de p: <1e-04).

Podemos também ver que não há diferença estatisticamente significativa (em α = 0,05) entre as técnicas B e C. Os intervalos de confiança entre as técnicas também apoiam estes resultados.

Assim, podemos concluir que o uso da técnica de estudo A leva a uma nota de exame estatisticamente significativamente maior para os alunos em comparação às técnicas B e C , mesmo após controlar a nota atual do aluno na turma.