1
$\begingroup$

I am in the last stage of my paper review and it is very close to getting submitted. However, there is a comment I am highly struggling with:

The reviewer wants to know what the ICC and the s.e. of the ICC is.

The model I have is a linear mixed effect model looking like:

mixedmodel <- lmer(y ~ x1 + x2 + x3 + (1 | countries), weights=weight, Alldata)

Before I had never heard about the ICC. When looking for documentation, I do find information on ICC, but not on how to calculate the s.e. of the ICC.

Can anyone please help me with this? I am using R.

$\endgroup$
3
$\begingroup$

You can easily get the ICC from the sjstats-package, which also offers bootstrapping. I would try bootstrapping to obtain se- and p-values (you can also take a look at this discussion in the r-sig-mailinglist for other approaches, but bootstrapping is one recommendation, two).

library(sjstats) # for ICC and bootstrapping
library(dplyr)   # for the pipe
library(lme4)

# compute ICC
fit <- lmer(Reaction ~ Days + (Days | Subject), sleepstudy)
icc(fit)
> Linear mixed model
>  Family: gaussian (identity)
> Formula: Reaction ~ Days + (Days | Subject)
>
>   ICC (Subject): 0.483090

# compute SE of ICC
dummy <- sleepstudy %>% 
  # generate 100 bootstrap replicates of dataset
  bootstrap(100) %>% 
  # run mixed effects regression on each bootstrap replicate
  mutate(models = lapply(.$strap, function(x) {
    lmer(Reaction ~ Days + (Days | Subject), data = x)
  })) %>% 
  # compute ICC for each "bootstrapped" regression
  mutate(icc = unlist(lapply(.$models, icc)))

# now compute SE and p-values for the bootstrapped ICC
boot_se(dummy, icc)
>   term   std.err
> 1  icc 0.1024704

boot_p(dummy, icc)
>   term      p.value
> 1  icc 1.936526e-08
  1. The above code takes the sleepstudy dataset and creates 100 bootstrap-samples, which are saved as "list-variable" named strap, in a new data frame.

  2. This list-variable (i.e. the data frame with this variable), containing the 100 bootstrap-samples, is then "piped" into the lapply-command to run 100 lmer-commands (which are saved in another list-variable named models, see mutate-command).

  3. This list-variable is then used to compute the ICC with the sjstats::icc function, so you have 100 ICC-values (for all models applied to the 100 bootstrap-samples) in the variable named icc.

Now boot_se and boot_p compute the standard error and p-value of any bootstrapped values, in this case the ICC-values.

$\endgroup$
  • $\begingroup$ Depending on your sample size, the computation of the mixed models for all bootstrap-samples may take very much time. I have extended the se()-function in my sjstats-package, which now a) also accepts objects from the icc()-function, and b) prints a progress bar, so you can guess when the calculation is finished. I will commit the changes to GitHub (github.com/sjPlot/sjstats) tonight (in about 3-4 hours), so you may want to install the GitHub-version of sjstats and then, after fitting your model, simply run se(icc(your_model)). $\endgroup$ – Daniel Aug 30 '16 at 13:10
  • $\begingroup$ I have over 60000 observations and the regression has 50 variables in Total. Nevertheless, it 'only' took about 30 minutes so I was very happy with it! $\endgroup$ – user33125 Aug 30 '16 at 13:13

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.