0

Is there a way to modify (overwrite) random-effects within a lmer-model? For fixed effects there is a slot called my_lmer@beta and I could alter the fixed effects using:

my_lmer@beta[1] <- 0.5

Is there a similar way to do this for random effects? Does a lmer-object already contain the random effects, or are the calculated later by e.g. ranef().

winwin
  • 384
  • 6
  • 20

1 Answers1

0

It would indeed be good to find out what you're trying to do. Modifying the internals of reference class objects is particularly dangerous -- you could accidentally modify copies of objects or cause segmentation faults ... from here,

library(lme4)
fm1 <- lmer(Reaction~Days+(1|Subject),sleepstudy) ## just for example
fm1@pp$getRefClass()$methods()

will show you the methods ... however, you have to go a bit deeper than that ... it turns out (l. 85 of src/predModule.cpp) that what b actually does it to take the internal u

 VectorXd merPredD::b(const double& f) const {return d_Lambdat.adjoint() * u(f);}

which in turn calls

VectorXd merPredD::u(const double& f) const {return d_u0 + f * d_delu;}

which means that in order to change b you would need to change the corresponding values of u0; at present I don't think that's possible.

For reference, this is some code (from here) that evaluates the deviance when the random effects are displaced by a (vector) z from their estimated values ...

rr <- m@resp                        ## extract response module
u0 <- getME(m,"u")                  ## conditional modes
L <- getME(m,"L")
## sd <- 1/getME(pp,"L")@x
## filled elements of L matrix==diag for simple case
## for more general case need the following -- still efficient
sd <- sqrt(diag(chol2inv(L)))
## fixed-effects contribution to linear predictor
fc <- getME(m,"X") %*% getME(m,"beta")
ZL <- t(getME(m,"Lambdat") %*% getME(m,"Zt"))
## evaluate the unscaled conditional density on the deviance scale
dc <- function(z) {
    uu <- u0 + z * sd    ## displace conditional modes
    ##  should still work if z is a vector (by recycling, because u values
    ##  applying to each group are stored adjacent to each other)
    rr$updateMu(fc + ZL %*% uu)     ## update linear predictor
    drc <- unname(as.vector(tapply(rr$devResid(), ff, sum)))
    uuc <- colSums(matrix(uu * uu,nrow=nvar))
    (drc + uuc)[grps]
}
Ben Bolker
  • 211,554
  • 25
  • 370
  • 453
  • Thanks for your responses! The intention was to compare predicted values from to model to a hypothetical second model with given estimates. Basically I wanted to change to estimates to a given prior and see how the predicted values chance. I thought it is easier to edit the lmer-object and use the `predict` method, instead of multiplying it out by myself. – winwin Aug 31 '16 at 07:48
  • actually, I think it would be easier to multiply it out yourself; you can extract the `X` and (the transposed) `Z` matrices via `getME(fitted_model,c("X","Zt"))`; then `X %*% beta + t(Zt) %*% b` should give you the predictions ... – Ben Bolker Sep 01 '16 at 02:27
  • Thanks for the tipp! – winwin Sep 02 '16 at 04:31