@Wolfgang ya dio una gran respuesta. Quiero ampliarlo un poco para mostrar que también puede llegar al ICC estimado de 0,75 en su conjunto de datos de ejemplo al implementar literalmente el algoritmo intuitivo de seleccionar aleatoriamente muchos pares de valores , de donde provienen los miembros de cada par mismo grupo, y luego simplemente calculando su correlación. Y luego, este mismo procedimiento puede aplicarse fácilmente a conjuntos de datos con grupos de cualquier tamaño, como también mostraré.y
Primero cargamos el conjunto de datos de @ Wolfgang (no se muestra aquí). Ahora definamos una función R simple que tome un data.frame y devuelva un solo par de observaciones seleccionadas al azar del mismo grupo:
get_random_pair <- function(df){
# select a random row
i <- sample(nrow(df), 1)
# select a random other row from the same group
# (the call to rep() here is admittedly odd, but it's to avoid unwanted
# behavior when the first argument to sample() has length 1)
j <- sample(rep(setdiff(which(dat$group==dat[i,"group"]), i), 2), 1)
# return the pair of y-values
c(df[i,"y"], df[j,"y"])
}
Aquí hay un ejemplo de lo que obtenemos si llamamos a esta función 10 veces en el conjunto de datos de @ Wolfgang:
test <- replicate(10, get_random_pair(dat))
t(test)
# [,1] [,2]
# [1,] 9 6
# [2,] 2 2
# [3,] 2 4
# [4,] 3 5
# [5,] 3 2
# [6,] 2 4
# [7,] 7 9
# [8,] 5 3
# [9,] 5 3
# [10,] 3 2
Ahora, para estimar el ICC, simplemente llamamos a esta función un gran número de veces y luego calculamos la correlación entre las dos columnas.
random_pairs <- replicate(100000, get_random_pair(dat))
cor(t(random_pairs))
# [,1] [,2]
# [1,] 1.0000000 0.7493072
# [2,] 0.7493072 1.0000000
Este mismo procedimiento se puede aplicar, sin modificaciones en absoluto, a conjuntos de datos con grupos de cualquier tamaño. Por ejemplo, creemos un conjunto de datos que consta de 100 grupos de 100 observaciones cada uno, con el ICC verdadero establecido en 0,75 como en el ejemplo de @ Wolfgang.
set.seed(12345)
group_effects <- scale(rnorm(100))*sqrt(4.5)
errors <- scale(rnorm(100*100))*sqrt(1.5)
dat <- data.frame(group = rep(1:100, each=100),
person = rep(1:100, times=100),
y = rep(group_effects, each=100) + errors)
stripchart(y ~ group, data=dat, pch=20, col=rgb(0,0,0,.1), ylab="group")
Estimando el ICC basado en los componentes de varianza de un modelo mixto, obtenemos:
library("lme4")
mod <- lmer(y ~ 1 + (1|group), data=dat, REML=FALSE)
summary(mod)
# Random effects:
# Groups Name Variance Std.Dev.
# group (Intercept) 4.502 2.122
# Residual 1.497 1.223
# Number of obs: 10000, groups: group, 100
4.502/(4.502 + 1.497)
# 0.7504584
Y si aplicamos el procedimiento de emparejamiento aleatorio, obtenemos
random_pairs <- replicate(100000, get_random_pair(dat))
cor(t(random_pairs))
# [,1] [,2]
# [1,] 1.0000000 0.7503004
# [2,] 0.7503004 1.0000000
que concuerda estrechamente con la estimación del componente de varianza.
Tenga en cuenta que si bien el procedimiento de emparejamiento aleatorio es intuitivo y didácticamente útil, el método ilustrado por @Wolfgang es en realidad mucho más inteligente. Para un conjunto de datos como este de tamaño 100 * 100, el número de emparejamientos únicos dentro del grupo (sin incluir los emparejamientos propios) es de 505,000, un número grande pero no astronómico, por lo que es totalmente posible para nosotros calcular la correlación del conjunto completamente agotado de todos los emparejamientos posibles, en lugar de necesitar muestrear aleatoriamente desde el conjunto de datos. Aquí hay una función para recuperar todos los emparejamientos posibles para el caso general con grupos de cualquier tamaño:
get_all_pairs <- function(df){
# do this for every group and combine the results into a matrix
do.call(rbind, by(df, df$group, function(group_df){
# get all possible pairs of indices
i <- expand.grid(seq(nrow(group_df)), seq(nrow(group_df)))
# remove self-pairings
i <- i[i[,1] != i[,2],]
# return a 2-column matrix of the corresponding y-values
cbind(group_df[i[,1], "y"], group_df[i[,2], "y"])
}))
}
Ahora, si aplicamos esta función al conjunto de datos 100 * 100 y calculamos la correlación, obtenemos:
cor(get_all_pairs(dat))
# [,1] [,2]
# [1,] 1.0000000 0.7504817
# [2,] 0.7504817 1.0000000
Lo que concuerda bien con las otras dos estimaciones, y en comparación con el procedimiento de emparejamiento aleatorio, es mucho más rápido de calcular, y también debería ser una estimación más eficiente en el sentido de tener menos varianza.
