组内相关系数在脑功能磁共振(fMRI)的应用

组内相关系数（The Intraclass Correlation Coefficient）是测量实验可重复性的指标。例如某项实验有n个测量目标值，重复测量了k次。
举个例子有一个nk的矩阵。BMS表示测量目标(row)之间的组间均方差（mean square)。WMS表示测量目标的组内均方差。JMS表示重复测量（column）之间的组间均方差。EMS表示残差的均方差。
以下是常用的三种ICC的计算方式：

假设有一个64的矩阵，

根据公式计算可得，

因此，有：

下面贴出ICC的计算matlab脚本。

function [ICC] = ICC(X, model, rating)
% ---------------------------------------------------------------
% Intraclass Correlation Coefficient
%
% [ICC] = ICC(X, model, rating)
%         返回矩阵 X的ICC值
%
% The input X is an n x k matrix, where n is the number of tar-
% gets and k the number of judges. ICCs are divided into three 
% models for two types of ratings.
%
% Model Types:
% '1' = random effects model => a different group of k judges, 
%       who have been selected from a larger population of jud-
%       ges, rate each target
% '2' = mixed effects model => the same group of k judges, who 
%       have been selected from a larger population of judges, 
%       rate each target
% '3' = 固定效应模型 => the same group of k judges, who 
%       are the only judges of interest, rate each targets.
%
% Rating Types:
% '1' = individual rating
% 'k' = mean of several ratings
% ---------------------------------------------------------------
% Script: 
   [n,k] = size(X);
   table = ANOVA_for_ICC(X, 'off');

   BMS = table{6,4}; WMS = table{7,4};
   JMS = table{3,4}; EMS = table{8,4};
   
   if isequal(rating,'1')
        if isequal(model,'1')
             ICC = (BMS - WMS) / ...
                   (BMS + (k-1)*WMS);
        elseif isequal(model,'2')
             ICC = (BMS - EMS) / ...
                   (BMS + (k-1)*EMS + k*(JMS-EMS)/n);
        elseif isequal(model,'3')
             ICC = (BMS - EMS) / ...
                   (BMS + (k-1)*EMS);
        else
             error('Model type does not exist.');
        end     
   elseif isequal(rating,'k')
        if isequal(model,'1')
             ICC = (BMS - WMS) / ...
                    BMS;
        elseif isequal(model,'2')
             ICC = (BMS - EMS) / ...
                   (BMS + (JMS-EMS)/n);
        elseif isequal(model,'3')
             ICC = (BMS - EMS) / ...
                    BMS;
        else
             error('Model type does not exist.');
        end
   else
        error('Rating type does not exist.');
   end
% ---------------------------------------------------------------

调用的函数为：

function [table] = ANOVA_for_ICC(X, displayopt)
% ---------------------------------------------------------------
% Analysis of Variance for Ratings (for Intraclass Correlations)
%
% [table] = ANOVA_for_ICC(X, displayopt)
%           returns the ANOVA summary table as a cell array
%
% The input X is an n x k matrix, where n is the number of tar-
% gets and k is he number of judges. The option of displaying
% an ANOVA summary table is handled by displayopt, which can be
% set to 'on' (default) or 'off.'
% ---------------------------------------------------------------
% Script:
   [m,n] = size(X);
   total = m*n;

   A  = sum(sum(X.^2));
   Bc = sum(sum(X,1).^2)/m;
   Br = sum(sum(X,2).^2)/n;  
   D  = sum(sum(X))^2/total;

   SSbc = Bc-D;
   SSwc = A-Bc;   
   SSbr = Br-D;
   SSwr = A-Br;
   SSe  = A-Br-Bc+D;

   DFbc = n-1;
   DFwc = n*(m-1);
   DFbr = m-1;   
   DFwr = m*(n-1);
   DFe  = DFbc*DFbr;

   MSbc = SSbc/DFbc;
   MSwc = SSwc/DFwc;
   MSbr = SSbr/DFbr;
   MSwr = SSwr/DFwr;
   MSe  = SSe/DFe;

   Fbc = MSbc/MSe; Pbc = 1-fcdf(Fbc, DFbc, DFe);
   Fbr = MSbr/MSe; Pbr = 1-fcdf(Fbr, DFbr, DFe);
       
   table = { 'Source of Variance' 'SS' 'df' 'MS' 'F' 'Prob>F'
             'Columns'             []   []   []   []  []
             '     betweeen'       SSbc DFbc MSbc Fbc Pbc
             '     within'         SSwc DFwc MSwc []  []
             'Rows'                []   []   []   []  []
             '     between'        SSbr DFbr MSbr Fbr Pbr 
             '     within'         SSwr DFwr MSwr []  []
             'Residual'            SSe  DFe  MSe  []  [] };

   if isequal(nargin, 1) || isequal(displayopt, 'on')
        statdisptable(table, 'ANOVA for ICC', ...
                             'ANOVA Summary Table', '');
   end
% ---------------------------------------------------------------

运行结果，如图所示：

FMRI工具包afni中有命令3dICC和3dlme以及3dlmer可计算ICC的值。