Chuan Tang - 3 months ago 10

R Question

The working data looks like:

`set.seed(1234)`

df <- data.frame(y = rnorm(1:30),

fac1 = as.factor(sample(c("A","B","C","D","E"),30, replace = T)),

fac2 = as.factor(sample(c("NY","NC","CA"),30,replace = T)),

x = rnorm(1:30))

The

`lme`

`library(lme4)`

mixed <- lmer(y ~ x + (1|fac1) + (1|fac2), data = df)

I used

`bootMer`

`mixed_boot_sum <- function(data){s <- sigma(data)`

c(beta = getME(data, "fixef"), theta = getME(data, "theta"), sigma = s)}

mixed_boot <- bootMer(mixed, FUN = mixed_boot_sum, nsim = 100, type = "parametric", use.u = FALSE)

My first question is how to obtain the coefficients(slope) of each individual levels of the two random effects from the bootstrapping results

`mixed_boot`

I have no problem extracting the coefficients(slope) from

`mixed`

`augment`

`broom`

`library(broom)`

mixed.coef <- augment(mixed, df)

However, it seems like

`broom`

`boot`

`mixed_boot`

I also tried to modify the

`mixed_boot_sum`

`mmList`

`Error in bootMer(mixed, FUN = mixed_boot_sum, nsim = 100, type = "parametric", :`

bootMer currently only handles functions that return numeric vectors

Furthermore, is it possible to obtain CI of both fixed&random effects by specifying

`FUN`

Now, I am very confused about the correct specifications for the

`FUN`

Answer

My first question is how to obtain the coefficients(slope) of each individual levels of the two random effects from the bootstrapping results mixed_boot ?

I'm not sure what you mean by "coefficients(slope) of each individual level". `broom::augment(mixed, df)`

gives the predictions (residuals, etc.) for *every observation*. If you want the predicted coefficients at each level I would try

```
mixed_boot_coefs <- function(fit){
unlist(coef(fit))
}
```

which for the original model gives

```
mixed_boot_coefs(mixed)
## fac1.(Intercept)1 fac1.(Intercept)2 fac1.(Intercept)3 fac1.(Intercept)4
## -0.4973925 -0.1210432 -0.3260958 0.2645979
## fac1.(Intercept)5 fac1.x1 fac1.x2 fac1.x3
## -0.6288728 0.2187408 0.2187408 0.2187408
## fac1.x4 fac1.x5 fac2.(Intercept)1 fac2.(Intercept)2
## 0.2187408 0.2187408 -0.2617613 -0.2617613
## ...
```

If you want the resulting object to be more clearly named you can use:

```
flatten <- function(cc) setNames(unlist(cc),
outer(rownames(cc),colnames(cc),
function(x,y) paste0(y,x)))
mixed_boot_coefs <- function(fit){
unlist(lapply(coef(fit),flatten))
}
```

When run through `bootMer`

/`confint`

/`boot::boot.ci`

these functions will give confidence intervals for each of these values (note that all of the slopes `fac`

are identical across groups because the model assumes random variation in the intercept only). In other words, whatever information you know how to extract from a fitted model (conditional modes/BLUPs [**W**.x**Z**`ranef`

], predicted intercepts and slopes for each level of the grouping variable [`coef`

], parameter estimates [`fixef`

, `getME`

], random-effects variances [`VarCorr`

], predictions under specific conditions [`predict`

] ...) can be used in `bootMer`

's `FUN`

argument, *as long as you can flatten its structure into a simple numeric vector*.