[R] Help with ICC

Patrick E. McKnight
2001-11-29T06:36:47Z

Tom,

Perhaps the following code will do the trick :

	ICC <- lm(x~j + i,data=iccdata)
	edf <- df.residual(ICC)
	ems <- deviance(ICC)/edf
	bdf <- anova(ICC)[1,1]
	bms <- anova(ICC)[1,2]/bdf
	msb <- bms
	jdf <- anova(ICC)[2,1]
	jms <- anova(ICC)[2,2]/jdf
	k <- jdf + 1
	msw <- ((ems*edf)+(jms*jdf))/(edf+jdf)
	wms <- msw
	n <- bdf + 1
	theta <- (msb - msw)/(k*msw)
	winer11 <- theta / (1+theta)
	winer1k <- (k*theta)/(1+k*theta)
	sf11 <- (bms - wms)/(bms+(k-1)*wms)
	sf21 <- (bms - ems)/((bms)+((k-1)*ems)+((k*(jms-ems))/n))
	sf31 <- (bms - ems)/(bms+((k-1)*ems))
	sf1k <- (bms - wms)/bms
	sf1k <- (bms - ems)/(bms+((jms-ems)/n))
	sf3k <- (bms - ems)/bms

Invoke a few print commands for the output desired and voila, Shrout and
Fleiss ICC's right from R.  I might have a few typos in the above code
but the general picture works.  I tested the critical parts.  In
addition, the parceling of the anova output is rather clunky.  I hope
others will show us a more efficient method of extracting the sums of
squares. Let me know if the results are consistent with your hand
calculations.

Cheers,

Patrick

NB:  The general approach above was abstracted from Robert M. Hamer's
SAS macro.
__________________________________________________________________________________________
Patrick E. McKnight
Evaluation Group for Analysis of Data
University of Arizona
Department of Psychology
Tucson, AZ  85721
520-621-5463
pem at theriver.com

Original message....

Hello, R-folks:

Here is a statement I use to make a data frame:

iccdata <- data.frame(i=rep(1:10,rep(2,10)),j=rep(1:2,10),
x=c(0.35011,0.11989,0.13081,0.09919,0.16000,0.12000,0.00000,0.00000,
0.44023,0.32977,2.67081,2.63919,0.09050,0.03950,0.44019,0.30981,0.59000,
0.57000,4.03000,3.77000))

<snip>

Now, how do you use (lme in) R to compute the Intraclass correlation for
these data?  I have written brute-force code to get the "right"
answer, but
I cannot get the proper mean squares from any regression model I have tried.
It seems like it should be trivial, but I cannot do it!  I need this to work
in order to apply the technique to a much larger data set.  I could use my
silly code, in SAS, but I want to see and understand the variance
components.  After studying Pinheiro and Bates for a while, I naively tried
the following, which is not what I want.


>> ICC <- groupedData(x~j|i,data=iccdata)
>> ICC.lme <- lme(x~j,data=ICC)
>

Can you help me?  Thanks in advance!

Tom



-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at
stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._