級内相関係数, Case2, ICC(2, 1), ICC(2, k)

Rの関数でICC(2,1)をやってみましょう

Rの関数を使用してICC(2,1)を算出

irr::icc(
    dat[,2:4],
    model="twoway",
    type="agreement",
    unit="single"
)

 Single Score Intraclass Correlation

   Model: twoway 
   Type : agreement 

   Subjects = 10 
     Raters = 3 
   ICC(A,1) = 0.942

 F-Test, H0: r0 = 0 ; H1: r0 > 0 
  F(9,18.7) = 46.1 , p = 7.12e-11 

 95%-Confidence Interval for ICC Population Values:
  0.844 < ICC < 0.984

結果のみでよければ、ここまでです

ここからは、もう少し詳しく学習してみましょう

統計モデル

モデルは繰り返し測定のない二元配置分散分析です。Shrout(1979)の論文では、以下のようになっています（ただし記号は異なります）。

$y_{ij}=\mu +{\alpha }_i+{\beta }_j+(\alpha \beta {)}_{ij}+w_{ij}$

$i=1,\dots ,nj=1,\dots ,k$

このモデルを簡単に説明すると、「検査結果＝全体の平均＋患者の影響＋評価者の影響+ 患者と評価者の組み合わせで生じる影響＋誤差による影響」となります。しかし繰り返し測定がないので、$(\alpha \beta {)}_{ij}\mathrm{と}{w}_{ij}$ を個別に推定することは困難です（Shrout, 1979, p422）。ここでは以下のモデルで考えます（ICCの結果には問題ありません）。

$y_{ij}=\mu +{\alpha }_i+{\beta }_j+e_{ij}$

${\alpha }_i\sim N(0,{\sigma }_T^2){\beta }_j\sim N(0,{\sigma }_U^2)e_{ij}\sim N(0,{\sigma }_E^2)$

• ${\sigma }_T^2$：真の分散（真の検査結果のバラツキ）
• ${\sigma }_U^2$：評価者による分散（測定誤差）
• ${\sigma }_E^2$：誤差による分散

$ICC(2,1)={\displaystyle \frac{{\sigma }_T^2}{{\sigma }_T^2+{\sigma }_U^2+{\sigma }_E^2}}$

The Expected Mean Squares（平均平方の期待値）

評価者数$＝ｋ$、患者数$＝ｎ$

患者平均平方和 (between-targets mean square) の期待値

$BMS=k{\sigma }_T^2+{\sigma }_E^2$

評価者平均平方和 (between-judges mean square) の期待値

$JMS=n{\sigma }_U^2+{\sigma }_E^2$

誤差平均平方和 (Mean Square Error) の期待値

$MSE={\sigma }_E^2$

$ICC(2,1)={\displaystyle \frac{BMS-MSE}{BMS+(k-1)MSE+{\displaystyle \frac{k(JMS-MSE)}{n}}}}$

（期待値の求め方は以下のページをご参照ください）

級内相関係数, Case1, ICC(1, 1), ICC(1, k)

Rを使いながら級内相関係数,（ICC(1,1)）を解説します。分散比の推定からICCの95％信頼区間まで求めます。リハビリテーション統計学では、信頼性の評価などに適用されます。

y2pt.com

2023.02.19

分散分析表

fit <- lm(data ~ ID + PT, data=dat2)
anova(fit)

> anova(fit)
Analysis of Variance Table

Response: data
          Df  Sum Sq Mean Sq F value   Pr(>F)    
ID         9 22162.7 2462.52  46.051 1.33e-10 ***
PT         2    19.5    9.73   0.182   0.8351    
Residuals 18   962.5   53.47      

分散分析表から平均平方和のみを取り出してICCを求めてみます

BMS <- anova(fit)[1, 3]
JMS <- anova(fit)[2, 3]
MSE <- anova(fit)[3, 3]
(ICC2 <- (BMS - MSE)/(BMS + (3-1)*MSE + 3*(JMS-MSE)/10))

> (ICC2 <- (BMS - MSE)/(BMS + (3-1)*MSE + 3*(JMS-MSE)/10))
[1] 0.9423787

Rの関数iccで算出した解と同じですね

ICC(2,1) の95％信頼区間

$ICC(2,1)={\displaystyle \frac{BMS-MSE}{BMS+(k-1)MSE+{\displaystyle \frac{k(JMS-MSE)}{n}}}}$

$={\displaystyle \frac{n(BMS-MSE)}{nBMS+(kn-k-n)MSE+kJMS}}$

ここからは、かなり難解ですので・・・

Shrout(1979)の論文から引用させていただきます

$ICC=ICC(2,1)$

$F_j={\displaystyle \frac{JMS}{MSE}}$

$v={\displaystyle \frac{(k-1)(n-1)(k\ast ICC\ast F_j+n(1+(k-1)ICC)-k\ast ICC{)}^2}{(n-1)k^2\ast IC{C}^2\ast F_j^2+(n(1+(k-1)\ast ICC)-k\ast ICC{)}^2}}$

ICC(2,1)の95％信頼区間は以下のようになるそうです

$Fu=F_{0.975}(n-1,v)$

$Fl=F_{0.975}(v,n-1)$

下限：$Lo={\displaystyle \frac{n(BMS-F^{\ast }MSE)}{F^{\ast }(kJMS+(kn-k-n)EMS)+nBMS)}}$

上限：$Up={\displaystyle \frac{n(F_{\ast }BMS-MSE)}{kJMS+(kn-k-n)EMS+n{F}_{\ast }BMS}}$

そのまま計算してみます（正確には「Rに計算してもらいます」です）

k <- 3
n <- 10
ICC <- ICC2
f <- JMS/MSE
v1 <- ((k-1)*(n-1))*(k*ICC*f+n*(1+(k-1)*ICC)-k*ICC)^2
v2 <- (n-1)*k^2*ICC^2*f^2+(n*(1+(k-1)*ICC)-k*ICC)^2
v <- v1/v2
F1 <- qf(0.975, 9, v)
F2 <- qf(0.975, v, 9)
(Lo <- n*(BMS-F1*MSE)/(F1*(k*JMS+(k*n-k-n)*MSE)+n*BMS))
(Up <- n*(F2*BMS-MSE)/(k*JMS+(k*n-k-n)*MSE+n*F2*BMS))

> (Lo <- n*(BMS-F1*MSE)/(F1*(k*JMS+(k*n-k-n)*MSE)+n*BMS))
[1] 0.8439671
> (Up <- n*(F2*BMS-MSE)/(k*JMS+(k*n-k-n)*MSE+n*F2*BMS))
[1] 0.9839527

Rの関数で算出した値と同じですね

詳細は下記の論文をご参照ください
SATTERTHWAITE, Franklin E. An approximate distribution of estimates of variance components. Biometrics bulletin, 1946, 2.6: 110-114.
FLEISS, Joseph L.; SHROUT, Patrick E. Approximate interval estimation for a certain intraclass correlation coefficient. Psychometrika, 1978, 43: 259-262.