# Course: EDF 6403 - Quantitative Foundations of Educational Research
# Stats Analysis 9: Repeated measures ANVOA
# Author: Nguyet Hoang | hoangt@ufl.edu
# Step 0: Load library
library(psych)
library(haven)
require(reshape2)
## Loading required package: reshape2
require(afex)
## Loading required package: afex
## Loading required package: lme4
## Loading required package: Matrix
## ************
## Welcome to afex. For support visit: http://afex.singmann.science/
## - Functions for ANOVAs: aov_car(), aov_ez(), and aov_4()
## - Methods for calculating p-values with mixed(): 'S', 'KR', 'LRT', and 'PB'
## - 'afex_aov' and 'mixed' objects can be passed to emmeans() for follow-up tests
## - Get and set global package options with: afex_options()
## - Set sum-to-zero contrasts globally: set_sum_contrasts()
## - For example analyses see: browseVignettes("afex")
## ************
##
## Attaching package: 'afex'
## The following object is masked from 'package:lme4':
##
## lmer
library(emmeans)
## Welcome to emmeans.
## Caution: You lose important information if you filter this package's results.
## See '? untidy'
library(ggplot2)
##
## Attaching package: 'ggplot2'
## The following objects are masked from 'package:psych':
##
## %+%, alpha
library(effectsize)
##
## Attaching package: 'effectsize'
## The following object is masked from 'package:psych':
##
## phi
library(moments)
library(scales)
##
## Attaching package: 'scales'
## The following objects are masked from 'package:psych':
##
## alpha, rescale
library(nortest)
# Step 1: Read the data
data <- read_sav("src/SA9_Baked_beans.sav")
head(data)
## # A tibble: 6 × 5
## ID week1 week2 week8 week9
## <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 26.6 12.4 91.9 74.3
## 2 2 16.6 11.6 93.2 82.8
## 3 3 17.2 23.3 77.6 88.3
## 4 4 7.00 17.6 86.5 99.0
## 5 5 26.7 22.1 86.7 77.2
## 6 6 19.6 28.7 93.0 93.0
# Step 2: Convert to the long format
data_long = melt(data=data, id.vars="ID", measure.vars=c("week1", "week2", "week8", "week9"),
variable.name="Week", value.name="Score")
head(data_long)
## ID Week Score
## 1 1 week1 26.614524
## 2 2 week1 16.616344
## 3 3 week1 17.156468
## 4 4 week1 6.995224
## 5 5 week1 26.731173
## 6 6 week1 19.566529
# Set things to factors that should be factors but aren't
data_long$ID = as.factor(data_long$ID)
# Step 3: Assess assumptions (which is skipped for this assignment)
# Step 4: Run the ANOVA
(model <- aov_car(Score ~ Week + Error(ID/Week), data=data_long))
## Anova Table (Type 3 tests)
##
## Response: Score
## Effect df MSE F ges p.value
## 1 Week 2.14, 40.59 89.02 504.58 *** .958 <.001
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
##
## Sphericity correction method: GG
# Obtain partial effect sizes
eta_squared(model, partial = TRUE)
## # Effect Size for ANOVA (Type III)
##
## Parameter | Eta2 (partial) | 95% CI
## -----------------------------------------
## Week | 0.96 | [0.95, 1.00]
##
## - One-sided CIs: upper bound fixed at [1.00].
omega_squared(model, partial = TRUE)
## # Effect Size for ANOVA (Type III)
##
## Parameter | Omega2 (partial) | 95% CI
## -------------------------------------------
## Week | 0.96 | [0.94, 1.00]
##
## - One-sided CIs: upper bound fixed at [1.00].
# Step 5: Run the follow-up tests
# Obtain emmeans across times
em <- emmeans(model, ~ Week)
# Assess polynomial trends
contrast(em, method = "poly")
## contrast estimate SE df t.ratio p.value
## linear 282.724 8.94 19 31.613 <0.0001
## quadratic -0.767 2.45 19 -0.313 0.7574
## cubic -126.559 8.96 19 -14.126 <0.0001
# Post hoc pairwise comparisons (which is also skipped being reported for this assignment)
contrast(em, method = "pairwise", adjust = "bonferroni")
## contrast estimate SE df t.ratio p.value
## week1 - week2 -3.34 1.97 19 -1.697 0.6361
## week1 - week8 -69.58 1.91 19 -36.403 <0.0001
## week1 - week9 -72.16 3.11 19 -23.185 <0.0001
## week2 - week8 -66.24 2.52 19 -26.311 <0.0001
## week2 - week9 -68.82 2.26 19 -30.493 <0.0001
## week8 - week9 -2.58 3.06 19 -0.842 1.0000
##
## P value adjustment: bonferroni method for 6 tests
# Step 6: Plot the means
# Add the means for each time
week_averages <- colMeans(data[, c("week1", "week2", "week8", "week9")], na.rm = TRUE)
week_df <- data.frame(
Times = names(week_averages),
Score = as.numeric(week_averages)
)
# A plot
ggplot(week_df, aes(x = Times, y = Score)) +
geom_line(data = week_df, aes(group = 1), color = "salmon", size = 1.2) +
geom_point(data = week_df, color = "red", size = 3) +
theme_minimal() + ylim(0, 100) +
labs(title = "Average Scores Across Weeks", x = "Week", y = "Average Score")
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
