行政院國家科學委員會專題研究計畫成果報告 計畫中文名稱: 利用微晶片資料診斷疾病準確性測度量之研究
計畫英文名稱: On the accuracy of diagnosing disease for microarray data 計 畫 編 號 : NSC93-2118-M006-005
執 行 期 限 : 93 年 8 月 1 日至 94 年 7 月 31 日 主 持 人: 馬瀰嘉 國立成功大學統計學研究所 碩士兼任助理: 蔡旻熹
中文摘要
對於一種診斷測量值的準確性,通常利用 ROC 曲線下面積來評估,而目前的醫學診斷 中,利用多個診斷測量值同時評估某種疾病 之方法是常見的,例如: 生物微晶片
(microarray)上數個基因的表現值,同時被用 來診斷某種疾病。針對同時使用多種不同診 斷測量值去評估某種疾病,本文提出一多變 量 ROC 曲線下面積來評估準確性,並在隨機 效應模型下導出其估計量及其信賴區間。另 外一無母數估計方法也被提出並利用統計模 擬比較兩估計方法。
關鍵字:ROC 曲線,無母數的估計量,多個 診斷測量值
Abstract
A statistical tool which is popular for describing the accuracy of diagnostic is the receive operating characteristic (ROC) curve. In fact, it is possible that several diagnostic
markers are used to detect a disease, e.g.
microarray data, and how one can measure their simultaneous diagnostic accuracy. Based on the random effect model, we derive an estimate of the area under multivariate ROC curve. A model-based estimation and an approximated confidence interval are derived by combining the confidence intervals of two parameters using Bonferroni method. A nonparametric
estimation will also be proposed and compared to the model-based estimate by a simulation study.
Keywords: AUROC, nonparametric estimation, multiple markers.
1. 前言
假設 X 為一個診斷測量值,若 則判斷為
有病。敏感度(sensitivity)的定義為病人真有病 經由診斷也為有病的機率,明確度(specificity) 的定義為病人沒有病經由診斷亦為沒病的機 率,令給定
c X ≥
=1
Y 時隨機變數 X 為變數 1,給 定
X
=0
Y 時為變數 。令 和 為 和
的累積分配函數(c.d.f.), 和
X2 G(⋅) F(⋅) X1
X2 g(⋅) f(⋅)為
和 的機率密度函數(p.d.f.),c 為切斷點 (cutoff point),由這些假設,敏感度和明確度 可以被表示如下: 敏感度=
X1
X2
) ( 1 ) (
) 1
|
(X c Y P X1 c G c
P ≥ = = ≥ = −
明確度=P(X <c|Y =2)=P(X2 <c)=F(c). AUROC 的公式如下:
∫
∫
≥
=
≥
=
−
−
=
1
0 1 1 2
1 0
) (
) ( ) (
)) ( 1 ( )) ( 1 (
X X P dc c f c X P
c F d c θ G
(1) 本文第二節針對同時使用 種不同診斷 測量值評估某種疾病時,提出多變量 ROC 曲 線及其曲線下面積(AUMROC);第三節介紹 提出的兩種估計方法;第四節以模擬實驗來 比較兩種估計方法;第五節為結論與討論。
l
1
2. 多變量 ROC 曲線下面積 假設Xr1
為有病這群人所得的 個診斷值,也 就是
l
) , , , ( ) 1
|
( 11 21 1
1 ≡ X Y = = X X Xl ′
Xr r L
和Xr2
為健康這群人所得的 個診斷值,也就 是
l
) , , , ( ) 2
|
( 12 22 2
2 ≡ X Y = = X X Xl ′
Xr r L
而Xr =(X1,X2,L,Xl)。假設Xr1
的 c.d.f.為 )
(x
G v 和 p.d.f.為 g v ;(x) Xr2
的 c.d.f.為 F v 和(x) p.d.f.為 f v 。如果(x) 或 …或
,則被檢驗者會被判定為有病的,意
即如果 和 和…和 ,則
被檢驗者會被判定為健康的。此時多變量敏 感度和多變量明確度可表示如下:
1
1 c
X ≥ X2 ≥c2
l
l c
X ≥
1
1 c
X < X2 <c2 Xl <cl
多變量敏感度= P(Xv1>cv) =1−G(cr) 多變量明確度= P(Xv2≤cv) =F(cr)
多變量 ROC 曲線可定義為從F(cr)至G(cr)的 函數 R(.),意即G(cr =) R(F(cr)),其中
l, cl
c c
cr =( 1, 2,L, )∈(−∞,∞) )
( 1 )
(c F c
F r = − r ,G(cr)=1−G(cr)。如果
和 …和 ,則被檢驗者
被判定為有病的,可令
和 ,
。因此,多變量敏感度和多變量明 確度可以被表示成
1
1 c
X ≥ X2 ≥c2 Xl ≥cl
l
l X
V X V X
V1 =− 1, 2 =− 2,L, =− ci' =−ci
l i=1,L,
多變量敏感度
=P(V11 ≤c1' andV21 ≤c2' andLandVl1 ≤cl')
=G*(cr')
1-多變量明確度
=1−P(V12 >c1' orV22 >c2' orLorVl2 >cl')
=P(V12 ≤c1' andV22 ≤c2' andLandVl2 ≤cl')
=F*(cr')
其中G* 和F*分別是 和 的累積分配函數。
) , , (V11 LVl1 )
, , (V12 L Vl2
那麼 AUMROC 的公式能被表示如下:
) (
) (
) (
) (
) ( ) (
2 1 12
11
2 1 12
11
' 2 '
1 12 1
0
' 1 '
1 11 1 0
'
* '
*
l l
l l
l l
l l
X X and and
X X P
V V and and
V V P
c V and and
c V dP
c V and and
c V P
c dF c G
≥
≥
=
≤
≤
=
≤
≤
≤
≤
=
=
∫
∫
L L
L L r θ r
(2)
令 是有病這群人中第 個病人
第
1
Xij i(i=1,K,n1) )
, , 1 (j l
j = K 個診斷測量值, 是健康這
群人中第
'j2
Xi
) , , 1
( ' 2
' i n
i = K 個病人第
) , , 1 (j l
j = K 個診斷測量值,那麼 AUMROC 也能被表示為 θ = P[Xij1 > Xi'j2, j=1,K,l]。 3. AUMROC 的估計方法
本文提出兩種不同估計 AUMROC 的方法。第 一種方法是假設在隨機效果模型下估計 AUMROC,第二種方法是使用無母數統計方 法估計。介紹如下:
(1) 隨機效果模型下的估計
假設Xijk符合隨機效果模型,表示如下:
2 , 1
, , , 1 , , , 1 ,
=
=
= +
+ +
= k
l j
n i
e a t
Xijk μk ik jk ijk K k K 其中μ 是第 種健康狀況的母體平均,也就k 是
k
μ 為有病這群人的總平均;1 μ 為健康這群2
人的總平均; 是第 種健康狀況下第 個人
的隨機效應,意即 ; 是第 種
健康狀況下個體經由第
tik k i
) , 0 (
~ t2
iid
ik N
t σ ajk k
j 種診斷測量值的隨
2
機效應,意即 ; 為誤差項,
假設 ;利用上面的假設,可以
計算此二因子隨機效果模型的變異數,令 為病人的均方和,MSA為診斷測量值的
均方和,而 為誤差均方和,公式如下:
) , 0 (
~ a2 iid
jk N
a σ eijk
) , 0 (
~ e2 iid
ijk N
e σ
MST
MSE
∑ ∑ ∑
= = = ⋅⋅ 2
k)
⋅ −X + −
= 2
1 1 1
2
1 2)
/(
(
k n
i l
j k i
k
n n X
MST
∑ ∑ ∑
= = = ⋅ − ⋅⋅ −
= 2
1 1 1
2 /2( 1) )
(
k n
i l
j
k jk
k
l X
X MSA
∑∑∑= = = − ⋅ − ⋅ + ⋅⋅ − + −
=
= 2
1 1 1
2 1 2
2 ( ) /( 1)( 2)
ˆ
k n
i l
j
k jk k i ijk e
k
n n l X X X X σ MSE
在上述模型假設下,AUMROCθ 如下:
[
l]
Z Zr
E ( ( ξ))
θ = Φ + , 及其估計量
[
l]
Z
m E ( (Zrˆ ˆ))
ˆ ξ
θ = Φ + , (3) 其中 Z 呈標準常態分配, 為標準常態分配 的累積分配函數,而
Φ
2 2 2
ˆ ˆ ˆ ˆ
e a
r t
σ σ σ +
= ,
σe
μ ξ μ
2ˆ ˆ ˆ= ˆ1− 2 ,
l MSE MST
t
= − ˆ2
σ 和
2 / ) ˆ (
2 1 2
n n
MSE MSA
a +
= −
σ 。
由於θ 會隨 增加而遞減,故將l θ 開 次方,得 到一個修正的測度量
l ϕ,其
估計量為 m ml
1
ˆ θˆ
ϕ = 。 (4) 由 delta 方法(delta method)可推導出 的變異 數估計量如下:
θˆm
ˆ), ˆ, ( ˆ)) ˆ, ( ˆ))(
ˆ, ( ( 2
ˆ) ( ˆ))
ˆ, ( ( ˆ) ( ˆ)) ˆ, ( ( ˆ ) (
' '
2 '
2 '
r Cov r g r g
r Var r g Var
r g Var
r
r m
ξ ξ
ξ
ξ ξ
ξ θ
ξ ξ
∧
∧
∧
∧
+
+
=
{ }
∫∞ Φ + + −Φ − Φ
=
= (ˆ,ˆ) 0 [ (ˆ ˆ)] [1 (ˆ ˆ)] ( )
ˆm g ξ r rZ ξ l rZ ξ l d Z
θ
, ) ( } ˆ) ˆ ( ˆ)]
ˆ ( 1 [
ˆ) (ˆ ˆ)]
(ˆ [ {
)| ,ˆ (ˆ ˆ)
ˆ, (
1 0
1 ) , ( ˆ) ˆ, ( '
dZ Z Z Z r Z
r l
Z Z r Z
r l
r r r g
g
l l
r r r
φ ξ φ
ξ
ξ φ
ξ
ξ ξ ξ ξ
−
− Φ
−
−
+ +
Φ
=
∂
= ∂
−
∞ −
=
∫
, ) ( ˆ)}
(ˆ ˆ)]
(ˆ 1 [
ˆ) (ˆ ˆ)]
(ˆ [ {
)| ,ˆ (ˆ ˆ)
ˆ, (
1 0
1 ) , ( ˆ) ˆ, ( '
dZ Z Z
r Z
r l
Z r Z
r l
r r g
g
l l
r r
φ ξ φ
ξ
ξ φ
ξ ξ
ξ ξ ξ ξ
ξ
−
− Φ
− +
+ +
Φ
=
∂
= ∂
−
∞ −
=
∫
而φ 是標準常態分配的機率密度函數,
⎪⎭
⎪⎬
⎫
⎪⎩
⎪⎨
⎧
⎥⎦
⎢ ⎤
⎣
⎡ Γ
− Γ
− +
= −
∧ 2
2 2 1 2 2
2 2
) (
) ˆ (
) 2 1 ˆ 2( 2ˆ
ˆ ˆ)
( υ
δ υ
δ υ υ
υ ξ σ
e
Var k
) 4 ( ) 2 ˆ (
) ˆ 2
( 2 2 ˆ
) 1 (ˆ
2 2 2 1
2 1 2 2 2
) ( 2
2 1
2
1 − −
− +
⎟⎟
⎟
⎠
⎞
⎜⎜
⎜
⎝
⎛
≈ −
+ + +
∧
m m
m
m m r m
Var
n n l
n n
ϕ l
) ˆ (
ˆ , ), 2 )(
1 (
2 1 2 2 2
2 1
1 ˆ 1
ˆ 2 ˆ
ˆ ˆ 2
1
n n l l
k
e t
k a
n n l
+
=
=
− +
−
=
+ +
− σ
σ σ
μ
δ μ
υ
2 2
2 2
1 2
2
2 2 2 2 1
)]
( [ 2
)]
( [ ) ( ) 1 (
ˆ ˆ ) ˆ
2( ) 1 ,ˆ (ˆ
υ υ υ
υ
σ υ σ σ μ μ ξ
Γ Γ
− Γ
− Γ
+
−
=
−
∧ −
e a
r t
Cov
,
其中
) /(
)
ˆ1 ( 1 2 2 2 222( 1)2
2 1
2 2
1+ − + −
+
= d MST d MSA nd MSTn d MSAl
m ,
) 2 )(
1
( 1 2
2 = l− n +n −
m , 和
1 2
2 1 n
d = n+ ,d2 =l,d3 =l(n1 +n2)/2, 由 delta 方法可得到
ˆ ) ˆ (
ˆ )
( 2 1
1 2
m l m
m Var
Var∧ ϕ = θ l−l ∧ θ (5) 由於θ 的信賴區間(confidence interval)可由兩
部分組成,一部份是ξ的信賴區間,另一部份 是 r 的信賴區間。ξ之100(1−α)%信賴區間為
⎟⎟⎠
⎜⎜ ⎞
⎝
⎛ ˆ± ∧ (ˆ)
2
ξ ξ Zα Var
3
r 之100(1−α)%信賴區間如下:
⎥⎥
⎦
⎤
⎢⎢
⎣
⎡
+ +
− + +
+
− +
− (
ˆ ( ) ;
(
) 2 ˆ (
2 1
2 1 ,
, 2 1
1 2 1 ,
, 1 2
2 2 2 1
n n l
n n F
n n l
l n n
Fα m m α m m
ψ ψ
) ) 2l
其中 d MSE
MSA d MST d
3 2
ˆ 1 +
ψ = 。
利用 Bonferroni 方法計算聯合的信賴區間,
則(Zr+ξ)的100(1−α)%信賴區間 表 示如下:
) , (L U
ˆ) ˆ (
) (
) 2 ˆ (
4 2
4 1
2 1
2 1 ,
,
ξ ψ ξ
α α
− ∧
+ + +
− +
= Z Var
n n l
l n n Z F
L
m m
ˆ) ˆ (
) (
) 2 ˆ (
4 2
4 1 1 2
2 1 ,
, 1
ξ ψ ξ
α α
∧
−
+ + +
+
− +
= Z Var
n n l
l n n Z F
U
m m
則ϕ的100(1−α)%信賴區間如下:
[
{E[(Φ(L))l]}1/l,{E[(Φ(U))l]}1/l]
(6) (2)無母數方法估計無母數估計是指不需要任何分配的假設,本 文提出一無母數方法來估計θ (AUMROC),
估計量表示如下:
∑∑∏
= = =>
= 1 2
'
'
1 1 1
1 2 2
1
) 1 (
ˆ n
i n
i l
j
j ij i
s I X X
n
θ n (7)
估計 公式如下:
ˆ ) ( s Var θ
2 1
2 2
1 1)ˆ
ˆ ( ˆ ˆ )
( nn
n n
Var s s B s
θ θ =θ + − + −
∧
其中
∑∑∑∏
∑ ∑∑∏
≠ =
= ≠ =
>
>
+
>
>
=
1 1 2
'
' '
1 2
' 2
'
' '
)
&
1 (
)
&
1 ( ˆ
1 2 1 2
2 1 1
1
1 2 1 2
2 1 1
n
i n
m n
i
j mj i
j ij i
l
j n
i n
i n
m
j ij m
j ij i
l
j
X X
X X n I
n
X X X
X n I
B n
ˆ ) ( s
Var∧ θ 是動差的估計方法,因不具不偏性,
所以推導Var(θˆs)另一不偏估計量如下:
) 1 (
)ˆ 1 ˆ (
ˆ ˆ )
~(
2 1 2 1
2 2 1
− +
−
− +
−
= +
n n n n
n n
V s s B s
θ
θ θ (8) θ 的100(1−α)%信賴區間如下:
⎟⎠
⎜ ⎞
⎝⎛ ˆ − ~(ˆ ), ˆ + ~(ˆ )
2
2 s s s
s Z V θ θ Z V θ
θ α α
修正測度量ϕ =θ1l的無母數估計量及其變異 數估計公式如下所示:
l
s s
ˆ1
ˆ θ
ϕ = , ~(ˆ )
ˆ ) (1 ˆ )
( 1 2
1
s s
s V
Var∧ ϕ = lθl− θ (9) ϕ的100(1−α)%信賴區間如下:
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛⎜⎝⎛ −s Z V s ⎟⎠⎞l ⎜⎝⎛ +s Z V s ⎟⎠⎞l
1 1
ˆ )
~( , ˆ
ˆ )
~( ˆ
2 2
θ θ
θ
θ α α
4.模擬結果
因篇幅關係省略模擬過程,模擬結果顯示,
以ϕˆm估計出的準確值較接近真正的ϕ值,亦 即平均偏差是較小的,ϕˆm的變異數較ϕˆs的變 異數小,信賴區間的涵蓋機率也較接近 95%。
5.結論與建議
假設隨機效果模型下估計方法的缺點是必須 適用在有相同測度單位的診斷測量值上,而 藉由無母數方法估計則無此缺點,綜合而 言,藉由無母數方法估計準確性之方法ϕˆs雖 然在模擬上其偏誤和均方誤都較ϕˆm稍大,但 其不需要任何分配或模型的假設,又考量到 不同診斷測量值在同一個體上測量的相關 性。當同時使用多個診斷測量值去評估某種 疾病時,本文推薦利用ϕˆs估計其準確性。
參考文獻
DeLong, E., DeLong, D., and Clarke-Pearson, D., “Comparing the areas under two or more correlated receiver operation characteristic curves : A nonparametric approach”, Biometrics, 44, 837-845, 1988.
Ma, M. C. and Ou, J. C., “The eatimation of the area under multivariate ROC Curve”, submitted, 2003.
4