# NOT RUN {
## Predict effects (direct, total)
m <- Shipley.SEM
e <- Shipley.SEM.Eff
dir <- dirEff(e); tot <- totEff(e)
f.dir <- predEff(m, effects = dir, type = "response")
f.tot <- predEff(m, effects = tot, type = "response")
## Using new data for predictors
d <- na.omit(Shipley)
xn <- c("lat", "DD", "Date", "Growth")
seq100 <- function(x) seq(min(x), max(x), length = 100)
nd <- data.frame(sapply(d[xn], seq100))
f.dir <- predEff(m, nd, dir, type = "response")
f.tot <- predEff(m, nd, tot, type = "response")
## Add CIs
# dir.b <- dirEff(e, "boot"); tot.b <- totEff(e, "boot")
# f.dir <- predEff(m, nd, dir, dir.b, type = "response")
# f.tot <- predEff(m, nd, tot, tot.b, type = "response")
## Predict an interactive effect (e.g. Live ~ Growth * DD)
xn <- c("Growth", "DD")
d[xn] <- scale(d[xn]) # scale predictors (improves fit)
m <- lme4::glmer(Live ~ Growth * DD + (1 | site) + (1 | tree),
family = binomial, data = d)
nd <- with(d, expand.grid(
Growth = seq100(Growth),
DD = mean(DD) + c(-sd(DD), sd(DD)) # two levels for DD
))
f <- predEff(m, nd, type = "response", interaction = "Growth:DD")
## Add CIs (need to bootstrap model...)
# system.time(B <- bootEff(m, ran.eff = "site", R = 1000))
# f <- predEff(m, nd, B, type = "response", interaction = "Growth:DD")
## Model-averaged predictions (several approaches)
m <- Shipley.Growth # candidate models (list)
w <- runif(length(m), 0, 1) # weights
e <- stdEff(m, w) # averaged effects
f1 <- predEff(m[[1]], effects = e) # pass avg. effects
f2 <- predEff(m, weights = w) # pass weights argument
f3 <- avgEst(predEff(m), w) # use avgEst function
stopifnot(all.equal(f1, f2))
stopifnot(all.equal(f2, f3))
## Compare model fitted values: predEff() vs. fitted()
m <- Shipley.SEM$Live
f1 <- predEff(m, unique.x = FALSE, re.form = NULL, type = "response")
f2 <- fitted(m)
stopifnot(all.equal(f1, f2))
## Compare predictions for standardised vs. raw effects
f1 <- predEff(m)
f2 <- predEff(m, use.raw = TRUE)
stopifnot(all.equal(f1, f2))
# }
Run the code above in your browser using DataLab