[R-sig-ME] Partial R2 in mixed models
Joan Molibo
joanmolibo at gmail.com
Thu Nov 20 17:00:54 CET 2014
Good afternoon;
First, I am not a statistician although I am in the way (I am a medical
doctor studying the grade in statistic, still in the first course). I would
like to compute de partial R-squared of the fixed effects of a model. I
have found a function from the LMERConvenienceFunctions package, but it
computes these from the lme4 anova extraction function, which gives a
sequential anova.
I have created an ad hoc function to compute the R-squared for each term
conditionally to the other terms in the model (based on the pamer.fnc). For
other hand, I have done something similar based with the recommedations
given by Snijders in his book (2nd edition, pages 111-113) to compute de R2
in two levels models.
I am not very sure of what I have done, but I think that the function works
so I would appreciate some light. For other hand, could I call the
calculated value as partial R-squared value of the fixed effects of a mixed
model?
Thank you very much.
As example:
partialR2 <- function(model){
# Based on SS
term <- attr(terms(model), "term.labels")
dv <- gsub(" ", "", gsub("(.*)~.*", "\\1", as.character(model at call)[2]))
ss.tot <- sum((model at frame[, dv] - mean(model at frame[, dv]))^2)
n <- length(term)
ss.var <- numeric(n)
form <- unlist(lapply(term, function(x) paste(paste(".~. -", x, sep =
""), x, sep = "+")))
for (i in 1:(n - 1)){
ss.var[i] <- as.data.frame(anova(update(model,
as.formula(form[i]))))[n, 2]
}
ss.var[n] <- as.data.frame(anova(model))[n, 2]
names(ss.var[n]) <- term[n]
out <- cbind(round(100 * ss.var/ss.tot, 5))
rownames(out) <- term
colnames(out) <- "Partial R2"
#Snijder
form <- paste(paste(".~ 1 + (1 |", names(ngrps(model)), sep = ""), ")",
sep = "")
m.null <- update(model, as.formula(form))
var.g.null <- VarCorr(m.null)[[1]][1]
var.r.null <- sigma(m.null)^2
var.null <- var.g.null + var.r.null
var.g.full <- VarCorr(model)[[1]][1]
var.r.full <- sigma(model)^2
var.full <- var.g.full + var.r.full
form <- unlist(lapply(term, function(x) paste(".~. -", x, sep = " ")))
var.red <- numeric(n)
for (i in n:1){
var.g.red <- VarCorr(update(model, as.formula(form[i])))[[1]][1]
var.r.red <- sigma(update(model, as.formula(form[i])))^2
var.red[i] <- var.g.red + var.r.red
}
out2 <- round(100 * (var.red - var.full)/var.null, 5)
return(cbind(out, Partial_R2_Snijders = out2))
}
##################################################
##################################################
library(LMERConvenienceFunctions)
library(foreign)
library(lme4)
dd <- read.dta("
http://www.ats.ucla.edu/stat/stata/examples/mlm_ma_snijders/mlbook1.dta")
str(dd)
m1 <- lmer(langpost ~ sex + ses + iq_perf + langpret + (1|schoolnr), data =
dd)
summary(m1)
anova(m1)
pamer.fnc(m1)
partialR2(m1)
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models
mailing list