Hoe voer je een ancova uit in r

Deze tutorial geeft een voorbeeld van hoe u een ANCOVA uitvoert in R.

Voorbeeld: ANCOVA in R

We voeren een ANCOVA uit om te testen of de studietechniek al dan niet invloed heeft op de examenresultaten, met behulp van de volgende variabelen:

• Technisch onderzoek : de onafhankelijke variabele die we willen analyseren
• Huidig leerlingcijfer : de covariabele waarmee we rekening willen houden
• Beoordelingsscore : de responsvariabelen die we willen analyseren

De volgende dataset bevat informatie over 90 studenten, willekeurig verdeeld in drie groepen van 30.

De dataset toont de studietechniek die elke leerling gebruikte (A, B of C) , hun huidige cijfer in de klas toen ze begonnen met het gebruik van de studietechniek, en het cijfer dat ze kregen op het examen nadat ze de studietechniek een maand lang ter voorbereiding hadden gebruikt. voor het examen. examen:

 #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

Stap 1: Verken de gegevens

Voordat we het ANCOVA-model passen, moeten we eerst de gegevens verkennen om deze beter te begrijpen en te verifiëren dat er geen extreme uitschieters zijn die de resultaten zouden kunnen vertekenen.

Ten eerste kunnen we een samenvatting van elke variabele in de dataset weergeven:

 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  

We kunnen zien dat elke studietechniekwaarde ( A, B en C) 30 keer in de gegevens voorkomt.

Ook kunnen we zien hoe de huidige scores van studenten bij aanvang van de studie verdeeld waren. De minimale score in de klas was 65,43, de maximale score 93,84 en het gemiddelde was 78,15.

Op dezelfde manier kunnen we zien dat de minimale score behaald op het examen 71,17 was, de maximale score 94,76 en het gemiddelde 83,38.

Dan kunnen we het dplyr- pakket gebruiken om eenvoudig de gemiddelde en standaardafwijking van huidige cijfers en examenresultaten voor elke studietechniek te vinden:

 #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

We kunnen zien dat de gemiddelde en standaardafwijkingen van het huidige leerjaar van studenten die elke studietechniek gebruiken ongeveer vergelijkbaar zijn.

Ook zien we dat de gemiddelde examenscore significant hoger is voor studenten die gebruik maken van studietechniek A vergeleken met technieken B en C.

We kunnen de verdeling van examenresultaten ook visualiseren op basis van studietechniek met behulp van boxplots :

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

Op dezelfde manier kunnen we ook boxplots gebruiken om de verdeling van de huidige cijfers te visualiseren op basis van de studietechniek:

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

Stap 2: Controleer de modelaannames

Zodra we wat basisgegevensverkenning hebben uitgevoerd en bekend zijn met de gegevens, moeten we verifiëren dat aan de volgende aannames voor ANCOVA wordt voldaan:

• De covariabele en de behandeling zijn onafhankelijk – het is noodzakelijk om te verifiëren dat de covariabele ( huidige graad) en de behandeling (onderzoekstechniek) onafhankelijk van elkaar zijn, omdat het toevoegen van een covariabele term in het model alleen zinvol is als de covariabele en de behandeling de behandeling werkt onafhankelijk op de responsvariabele ( onderzoek ).
• Homogeniteit van variantie – we moeten controleren of de varianties tussen groepen gelijk zijn

Om te verifiëren dat de covariabele en de behandeling onafhankelijk zijn, kunnen we een ANOVA uitvoeren met de huidige graad als responsvariabele en de onderzoekstechniek als voorspellende variabele:

 #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    

De p-waarde is groter dan 0,05, dus de covariabele ( huidige graad) en behandeling ( onderzoekstechniek ) lijken onafhankelijk.

Om vervolgens te controleren of er homogeniteit van variantie tussen groepen bestaat, kunnen we de Levene-test uitvoeren:

 #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   

De p-waarde van de test is gelijk aan 0,0001961, wat aangeeft dat de varianties tussen de groepen niet gelijk zijn. Hoewel we een transformatie van de gegevens zouden kunnen proberen om dit probleem te corrigeren, zullen we ons op dit moment niet al te veel zorgen maken over variantieverschillen.

Stap 3: Pas het ANCOVA-model aan

Vervolgens passen we het ANCOVA-model aan met behulp van de examenscore als de responsvariabele, de onderzoekstechniek als de voorspellende (of ‘behandelings’) variabele en het huidige cijfer als de covariabele.

We zullen de functie Anova() in het pakket gebruiken omdat we dit doen, gewoon om te kunnen specificeren dat we de kwadratensom van type III voor het model willen gebruiken, omdat de kwadratensom van type I afhangt van de volgorde waarin de Predictors worden in het model ingevoerd:

 #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         

We kunnen zien dat de p-waarde voor techniek extreem laag is, wat aangeeft dat studietechniek een statistisch significant effect heeft op examenscores, zelfs na controle voor het huidige cijfer.

Stap 4: Post-hoc testen

Hoewel de ANCOVA-resultaten ons vertelden dat de studietechniek een statistisch significant effect had op de examenscores, moeten we post-hoctesten uitvoeren om te bepalen welke studietechnieken van elkaar verschillen.

Om dit te doen, kunnen we de functie glht() in het multcomp- pakket in R gebruiken om de Tukey-test uit te voeren voor meerdere vergelijkingen:

 #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

Uit het resultaat kunnen we zien dat er een statistisch significant verschil is (bij α = 0,05) in de examenresultaten tussen het bestuderen van techniek A en het bestuderen van techniek B (p-waarde: 0,000109) en tussen techniek A en techniek C (p-waarde: <1e-04).

We kunnen ook zien dat er geen statistisch significant verschil is (bij α = 0,05) tussen technieken B en C. Betrouwbaarheidsintervallen tussen technieken ondersteunen deze bevindingen ook.

We kunnen dus concluderen dat het gebruik van studietechniek A leidt tot een statistisch significant hoger examencijfer voor studenten vergeleken met technieken B en C , zelfs na controle voor het huidige cijfer van de student in de klas.