如何计算 R 中的类内相关性 (ICC)？

问题描述

我有一个包含 200 个变量、94 个主题的长格式数据集，每个主题对每个变量都有 1 到 3 个测量值。

例如：

ID measurement var1 var2 . . .
1  1           2    6
1  2           3    8 
1  3           6    12
2  1           3    9
2  2           4    4
2  3           5    3 
3  1           1    11
3  2           1    4
.  .           .    .
.  .           .    .
.  .           .    . 

但是，某些变量的三个测量值之一具有缺失值。有人建议我，在用受试者的均值插补缺失值之前，我应该使用重复测量方差分析或混合模型，以确认测量的可重复性。

我发现计算 ICC 的第一件事是 psych 包中的 ICC() 函数。但是，据我所知，这要求数据每个主题有一行，每个测量有一行，由于我需要单独计算 ICC 的 200 个变量，这会变得更加复杂。我确实继续计算了单个变量的 ICC，并获得了以下输出：

Intraclass correlation coefficients 
                         type  ICC   F df1 df2             p lower bound upper bound
Single_raters_absolute   ICC1 0.38 2.8  93 188 0.00000000067        0.27        0.49
Single_random_raters     ICC2 0.38 2.8  93 186 0.00000000068        0.27        0.49
Single_fixed_raters      ICC3 0.38 2.8  93 186 0.00000000068        0.27        0.49
Average_raters_absolute ICC1k 0.65 2.8  93 188 0.00000000067        0.53        0.74
Average_random_raters   ICC2k 0.65 2.8  93 186 0.00000000068        0.53        0.74
Average_fixed_raters    ICC3k 0.65 2.8  93 186 0.00000000068        0.53        0.74

 Number of subjects = 94     Number of Judges =  3

接下来，我尝试使用混合模型计算 ICC。使用此代码：

m1 <- lme(var1 ~ measurement,random=~1|ID,data=mydata,na.action=na.omit)
summary(m1)

输出如下：

Linear mixed-effects model fit by REML
 Data: mydata
        AIC       BIC   logLik
  -1917.113 -1902.948 962.5564

Random effects:
 Formula: ~1 | ORIGINAL_ID
        (Intercept)    Residual
StdDev: 0.003568426 0.004550419

Fixed effects: var1 ~ measurement
                          Value    Std.Error  DF  t-value p-value
(Intercept)         0.003998953 0.0008388997 162 4.766902  0.0000
measurement         0.000473053 0.0003593452 162 1.316429  0.1899
 Correlation: 
                    (Intr)
measurement         -0.83 

Standardized Within-Group Residuals:
        Min          Q1         Med          Q3         Max 
-3.35050264 -0.30417725 -0.03383329  0.25106803 12.15267443 

Number of Observations: 257
Number of Groups: 94 

这是用于评估 ICC 的正确模型吗？我不清楚相关性 (Intr) 测量的是什么，它与使用 ICC() 获得的 ICC 不同。

这是我第一次计算和使用类内相关性，因此感谢您的帮助！

解决方法

使用模拟数据集...

set.seed(42)  
n <- 6
dat <- data.frame(id=rep(1:n,2),group= as.factor(rep(LETTERS[1:2],n/2)),V1 = rnorm(n),V2 = runif(n*2,min=0,max=100),V3 = runif(n*2,V4 = runif(n*2,V5 = runif(n*2,max=100))

正在加载一些库...

library(lme4)
library(purrr)
library(tidyr)
# Add list of variable names to the vector below...
var_list <- c("V1","V2","V3","V4","V5")

map_dfr() 来自 purrr 库。我使用 lme4::VarCorr() 来获取每个级别的差异。

map_dfr(var_list,function(x){
                    formula_mlm = as.formula(paste0(x,"~ group + (1|id)"));
                    model_fit = lmer(formula_mlm,data=dat);
                    re_variances = VarCorr(model_fit,comp="Variance") %>% 
                      data.frame() %>% 
                      dplyr::mutate(variable = x);
                    return(re_variances)
                  }) %>% 
  dplyr::select(variable,grp,vcov) %>% 
  pivot_wider(names_from="grp",values_from="vcov") %>% 
  dplyr::mutate(icc = id/(id+Residual))