Re: [R] Variance Component/ICC Confidence Intervals via Bootstrap or Jackknife

From: Spencer Graves <spencer.graves_at_pdf.com>
Date: Sun 29 Oct 2006 - 19:06:23 GMT

      I can think of two ways to get confidence intervals on intraclass correlations (ICCs) and more accurate intervals for variance components: (1) modifying 'simulate.lme' to store the estimated variance components as well as "logLik" and (2) using 'lmer' and 'mcmcsamp' in library(lme4).

      The difficulty with (1) is that you have to make a local copy of 'simulate.lme', then figure out where and how to modify it. I've just looked at the code, and it looks like the required modifications should be fairly straightforward. The problem with (2) is that you would have to learn the syntax for a nested model for 'lmer'. It's different from that for 'lme' but not difficult. The 'lmer' function is newer and better in many ways but is not as well documented and does not have as many helper functions yet. The best documentation available for it may be the 'MlmSoftRev' vignette in the 'mlmRev' package plus the "R News 5/1" article from May 2005. If you are not familiar with vignettes, RSiteSearch("graves vignette") produced 74 hits for me just now. Find 'vignette' in the first hit led me to an earlier description of vignettes.

      If it were my problem, I'd probably try the second, though I might try both and compare. If you try them both, I'd be interested in the comparison. If the answers were substantially different, I'd worry.

      Hope this helps. 
      Spencer Graves

Rick Bilonick wrote:
> I'm using the lme function in nmle to estimate the variance components
> of a fully nested two-level model:
>
> Y_ijk = mu + a_i + b_j(i) + e_k(j(i))
>
> lme computes estimates of the variances for a, b, and e, call them v_a,
> v_b, and v_e, and I can use the intervals function to get confidence
> intervals. My understanding is that these intervals are probably not
> that robust plus I need intervals on the intraclass correlation
> coefficients:
>
> v_a/(v_a + v_b + v_e)
>
> and
>
> (v_a + v_b)/(v_a + v_b + v_e).
>
> I would like to use a bootstrap or a jackknife estimate of these
> confidence intervals. Is there any package available?
>
> I've tried to write the R code myself but I'm not sure if I understand
> exactly how to do it. The way I've tried to do a bootstrap estimate is
> to sample with replacement for i units, then sample with replacement the
> j(i) units, and then finally sample with replacement the k(j(i)) units.
>
> But the few references I've been able to track down (Arvesen, Biometrcs,
> 1970 is one), seem to say that I should just sample with replacement the
> i units. Plus they seem to indicate that a log transform is needed. The
> Arvesen reference used something like using log(v_a/v_e) as an estimator
> for sigma_a^2/sigma_e^2 and getting an interval and then transforming to
> get to an interval for the ICC (although it's not clear to me how to get
> the other ICC in a two-level nested design).
>
> Any insights would be appreciated.
>
> Rick B.
>
> ______________________________________________
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>



R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Mon Oct 30 06:11:56 2006

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.1.8, at Mon 30 Oct 2006 - 09:30:14 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.