A Macro to Calculate Kappa Statistics for Categorizations by Multiple Raters

原文載點：http://www2.sas.com/proceedings/sugi30/155-30.pdf

在類別資料分析中，我們經常會看到這樣的設計：訪問一組人兩次，然後檢定兩次的答案有沒有一致性。這樣的統計分析，可以用 Kappa 統計量來檢定。可是，當重複訪問的次數超過兩次，那 Kappa 統計量還可以用嗎？Bin Chen，Dennis Zaebst 和 Lynn Seel 等人在 2005 年的 SUGI 30 發表了利用 Kappa 統計量做多重一致性比較的分析，並分享一個新的 macro 在文章裡面。


KAPPA STATISTICS

Kappa 統計量首先是在 1960 年由 Cohen 所發明出來，之後一直被廣泛應用到現在。原始的理論是架構在 2X2 table 上，現在就先稍微介紹 nXn table (n≧2) 時的理論基礎。

假設有個評量（以下用 rater 來表示）裡面有 j=1~k 個選項，則在 n 個人重複測了 m 次之後，第 j 個選項被勾選的總比例是：



其中分子的 Xij 表示第 i 個人勾選第 j 個選項的次數。則第 j 個選項的 Kappa 統計量可以用下式表示：



因此，總 Kappa 統計量就是上式的加權平均：



其標準誤是：



MKAPPA MACRO

作者們開發了一個視窗化的 macro 程式，名為 MKAPPA。使用者可以直接在視窗內輸入一些設定即可完成多重一致性的檢定，而不需要寫上任何一個程式碼。視窗畫面如下所示：



APPLICATION EXAMPLE

作者簡單地提供了一個小小的範例來示範這個程式的用法。假設有下面這樣的資料：



每個人一份評量做三次，所以總共有三個 rater，每個 rater 內有五個選項。經過簡單的統計之後，每個 rater 內各選項被選到的次數分配表如下所示：



經過 MKAPPA 的計算，結果如下：



由於 p-value 小於顯著水準 0.05，拒絕虛無假設，因此可以下一個結論，那就是受測群三次做出來的結果並沒有一致性存在。

這個 MKAPPA 還是有一些使用上的限制，主要是處理 missing data 上。作者有提到一旦出現 missing data 的話程式便會強迫終止。但是在範例中的 rater1 有一個 missing data，程式卻仍舊可以執行。初步的研判應該是作者將 missing data 獨立成另一個選項來避免程式執行失敗。

最後，MKAPPA 的 macro 如下所示：

%window setmacro rows = 18
#2 "Data location : " dloc 80
#4 "Dataset name : " dset 40
#6 "Number of raters : " Nrater 8
#8 "Number of rating categories : " Ncat 8
#10 "Output file name : " outfile
;
%let dloc = c:\documents and settings\bin chen\my documents\projects\KappaMacro\ ;
%let dset = testno1 ;
%let Nrater = 3 ;
%let Ncat = 5 ;
%let outfile = output.rtf ;
%display setmacro ;
%macro mkappa;
libname kap "&dloc";
data temp; set kap.&dset;
subj=_n_;
data temp; set temp;

%do i=1 %to &nrater;
rate=rater&i; rater=&i; output; %end;
run;
proc freq data=temp noprint; tables subj*rate / out=xij ;
run;
data subj; set xij; keep subj;
proc sort data=subj nodupkey out=subj; by subj; run;
proc sql noprint; select count(*) into :tsub from subj;
data temp2; do subj=1 to &tsub; do rate=1 to &ncat; output; end; end; run;
data xij; set xij; drop percent;
proc sort data=xij; by subj rate;
proc sort data=temp2; by subj rate;
data xij; merge xij temp2(in=a); by subj rate; if a;
data xij; set xij; if count=. then count=0; run ;
data xij; set xij; xmx=count*(&nrater-count); run;
proc sort data=xij; by rate;
data xij; set xij; by rate;
retain x xmx2;
if first.rate then do; x=0; xmx2=0; end;
x=x+count; xmx2=xmx2+xmx;
if last.rate then output;
run;
data xij; set xij; keep rate x xmx2 ; run;
data xij; set xij;
p=x/(&tsub*&nrater);
pq=p*(1-P);
pqp=p*(1-p)*(1-2*p);
kj=abs(1 - xmx2/(&tsub*&nrater*(&nrater-1)*pq));
numj=kj*pq;
run;
data xij; set xij; if p=0 then do;
kj=0; numj=0; end; run;
data final; set xij end=end; retain num den pqpsum;
if _n_=1 then do; num=0; den=0; pqpsum=0; end;
num=num+numj;
den=den+pq;
pqpsum=pqpsum+pqp;
if end then output;
run;
data kap.koutput; set final; keep Kappa SE pvalue;
Kappa=num/den;
SE=sqrt(2*(den*den-pqpsum))/(den*sqrt(&tsub*&nrater*(&nrater-1)));
pvalue=1-probnorm(kappa/se);
run;
ods rtf file="&dloc&outfile" ;
options nonumber ;
goptions dev=png target=png htitle=1.5 htext=1.2 ftext=swissb ftitle=swissb
xmax=6 ymax=5 ;
proc print data=kap.koutput noobs label;
label kappa="Kappa" SE="Standard Error" pvalue="Prob>Z";
title "Kappa statistics for &Ncat-category ratings by &Nrater raters";
run;
ods rtf close;
%mend mkappa;

CONTACT INFORMATION
Your comments and questions are valued and encouraged. Please contact the authors at:
Bin Chen
WESTAT Inc.
5555 Ridge Ave
MS-R44
Cincinnati, OH 45213
Phone: 513-841-4316
Fax: 513-841-4470
Email: bchen2@cdc.gov