# BALANCED DATA (how = "aov" and "lme"/"lmer" do YES provide the same value)
str(InsectSprays)
icc_11(x = InsectSprays$"count", grp = InsectSprays$"spray", how = "aov")
icc_11(x = InsectSprays$"count", grp = InsectSprays$"spray", how = "lme")
icc_11(x = InsectSprays$"count", grp = InsectSprays$"spray", how = "lmer")
icc_11(x = InsectSprays$"count", grp = InsectSprays$"spray",
how = "raw") # biased estimator and not recommended. Only available for teaching purposes.
# UN-BALANCED DATA (how = "aov" and "lme"/"lmer" do NOT provide the same value)
dat <- as.data.frame(lmeInfo::Bryant2016)
icc_11(x = dat$"outcome", grp = dat$"case", how = "aov")
icc_11(x = dat$"outcome", grp = dat$"case", how = "lme")
icc_11(x = dat$"outcome", grp = dat$"case", how = "lmer")
icc_11(x = dat$"outcome", grp = dat$"case", how = "lme", REML = FALSE)
icc_11(x = dat$"outcome", grp = dat$"case", how = "lmer", REML = FALSE)
# how = "lme" does not account for any correlation structure
lme_obj <- nlme::lme(outcome ~ 1, random = ~ 1 | case,
data = dat, na.action = na.exclude,
correlation = nlme::corAR1(form = ~ 1 | case), method = "REML")
var_corr <- nlme::VarCorr(lme_obj) # VarCorr.lme
vars <- as.double(var_corr[, "Variance"])
btw <- vars[1]
wth <- vars[2]
btw / (btw + wth)
Run the code above in your browser using DataLab