劉 妍 王紅宇 田嬌妮 劉桂芬△
無金標準條件下含協變量的ROC曲線估計方法*
劉 妍1王紅宇2田嬌妮1劉桂芬1△
目的闡明無金標準條件下,考慮協變量后估計ROC曲線的兩部貝葉斯模型。方法 介紹兩部貝葉斯模型,結合實例,篩選無金標準條件下ROC曲線的影響因素,考慮協變量影響后,估計ROC曲線。結果 兩部貝葉斯模型不僅可探討協變量對疾病狀態的影響,而且可探討協變量對診斷試驗結果的影響,同時可計算不同協變量取值條件下ROC曲線下面積。結論 兩部貝葉斯模型可有效地解決無金標準條件下,考慮協變量影響的ROC曲線估計問題。
兩部貝葉斯模型 ROC曲線 無金標準 協變量
*:國家自然科學基金項目(編號81172774);山西省自然科學基金項目(編號2009011005-2)
1.山西醫科大學衛生統計教研室(030001)
2.山西醫科大學第二醫院
△通訊作者:劉桂芬,E-mail:liugf66@gmail.com
診斷試驗是臨床研究的重要組成部分,它可在篩檢試驗的基礎上,進一步把患者、疑似病例和需鑒別的其他疾病區別開來,且病人的療效評價、預后估計等在一定程度上也都依賴于診斷試驗,所以科學地評價診斷試驗,可提供給患者關于疾病的可靠信息、影響醫生制定治療計劃,有效避免不必要的資源浪費。受試者工作特征曲線(receiver operating characteristic curve,ROC曲線)是目前評價診斷試驗常用的一種方法,有關含協變量的ROC曲線分析,一般方法都必須基于金標準存在的假設。然而在許多疾病狀態下,由于“金標準”不存在或十分昂貴或執行起來不符合實際,有時很難甚至不可能建立一個權威的“金標準”。這使在許多ROC曲線分析中,不自覺地采用不完善評價診斷試驗標準,從而導致估計的ROC曲線估計偏差加大。
本文擬闡明一種無金標準條件下考慮協變量影響后,估計ROC曲線的方法——兩部貝葉斯模型,它可在考慮協變量對疾病狀態影響的同時,考慮協變量對診斷試驗結果的影響,并可計算不同協變量取值條件下ROC曲線下面積(AUC),從而進行不同試驗準確度的比較。
兩部貝葉斯模型(two-part bayesian model)首先篩選影響疾病狀態的協變量,該過程可用logistic回歸模型來擬合;第二個過程即篩選影響試驗結果的協變量,該過程可采用線性模型來擬合。
隨機抽取含量為n的樣本,設第i個個體(i=1,…n),Di(Di=0或1)為其客觀但未知的二分類疾病狀態,Ti為目標診斷試驗的試驗結果,Ri為采用不完善標準(參照試驗)進行診斷的試驗結果(Ri=1:陽性,Ri=0:陰性)。假定兩試驗相互獨立,即在Di條件下Ti與Ri相互獨立。
對于第 i個個體,設有 K 個協變量 Xi,1,…,Xi,K,可能影響疾病狀態Di或目標診斷試驗結果Ti,或者二者都受影響。令 Xi=(1,Xi,1,…,Xi,K)',為篩選影響疾病狀態的協變量,建立logistic回歸模型,定義為:

為篩選在Di條件下,影響目標試驗結果Ti的協變量,建立Ti條件均值線性模型,記作:

式中:p:任意截斷點處的(1-特異度),α={E(Ti|Di=1)-E(Ti|Di=0)}/σD,β = σH/σD,由此計算得到的ROC曲線下面積記作:

對于參照試驗結果Ri,建立自然誤分類模型(naturemisclassification model),模型如下:

式中:θD和θH:該參照試驗的假陰性率和假陽性率,二者同樣也是未知的。對于第i個個體,Di、Ti與Ri的聯合似然函數表示為:


為每個參數設定先驗分布,假定 λD、λH、1/σ、1/σ、θD、θH和 φ 的先驗分布相互獨立。對于每個 λD,0,λD,1…λD,K、λH,0,λH,1…λH,K和 φ0,φ1…φK選用正態無信息先驗分布 N(0,1000),1/σ和 1/σ選用伽馬無信息先驗分布γ(0.001,0.001),θD和 θH選用貝塔無信息先驗分布β(0.5,0.5)。根據先驗分布和似然函數,采用MCMC方法得出后驗分布的參數估計值。本研究使用Gibbs抽樣構造馬爾科夫鏈來模擬參數的后驗分布,得到所有待估參數的完全條件分布。根據估計獲得的模型參數,繪制出協變量取值分層的多條ROC曲線,并計算相應的ROC曲線下面積,進而準確地評價該指標分層診斷的價值。整個過程應用MATLAB7.8實現。
在金標準D和參照標準R對應的θD、θH均已知的情況下進行模擬,隨機產生服從兩部貝葉斯模型要求的隨機變量T,條件均值結構為

式中,X1,X2均服從均勻分布,X1,X2~ U(0,1),(T|D=1),(T|D=0) ~ N(0,1)。取 θD=0.1,θH=0.2 進行模擬研究,采用兩部貝葉斯模型進行分析,取無信息先驗分布(如前文),模擬50000次,退火算法(burnin)迭代次數取5000次,退火后迭代次數即 Monte Carlo樣本量為45000,參數估計結果見表1和表2。

表1 兩部貝葉斯模型模擬研究參數估計結果分析(1)

表2 兩部貝葉斯模型模擬研究參數估計結果分析(2)
由表模擬結果可見:(1)隨樣本含量的增加,參數估計的準確度越來越好,當樣本含量增加到150時,參數估計值與真值接近,但Chong Wang等人研究表明,為了將居住環境、生活習慣等混雜因素對診斷結果的影響降到最小,應在允許的范圍內盡量抽取不同居住地、不同生活習慣的樣本,也即盡可能增大樣本含量。(2)兩部貝葉斯估計的中位數和均數雖相差不大,但中位數更接近真實值。(3)在樣本含量低于100時,模型估計得到的值與設定的真實值的一致率不足85%;當樣本含量逐漸增大到150時,二者的一致率達90%;當樣本含量大于300時,參數估計結果幾乎與真值接近。由此可見,兩部貝葉斯模型對于無金標準條件下考慮協變量后的ROC曲線估計,方法可行,結果準確,樣本含量在150例以上,診斷試驗準確度更高。
國際公認只有冠狀動脈造影(CAG)才是冠心病診斷的金標準。但由于CA G術是利用導管對冠狀動脈解剖進行的放射影像學檢查,屬一種創傷性介入診斷技術,患者依從性較差,這為接受冠心病的診斷試驗帶來困難。欲在未進行CAG檢查條件下,正確地對冠心病做出診斷,本研究收集了有臨床癥狀自我感覺不適的疑似冠心病患者168例,檢測其收縮壓、甘油三酯含量,并以心電圖作為參照標準。以24h動態心電圖(Holter)中心率變異指標——竇性心搏RR間期標準差(SDNN),作為診斷區分冠心病的主要指標。收集分析數據與統計描述結果見表3。

表3 168例疑似冠心病患者部分檢測結果及統計描述表
采用MCMC法進行模型參數估計,先驗分布選取無信息先驗,如前文,模擬50000次,退火算法迭代次數取5000次,退火后迭代次數即Monte Carlo樣本量為45000,按所得金標準分組的統計描述及分析結果見表4和表5。
表4 168例疑似冠心病患者各指標統計描述(±s)

表4 168例疑似冠心病患者各指標統計描述(±s)
指標 冠心病患者(n=53) 非冠心病患者(n=115)31收縮壓(mmHg) 138.47±19.24 131.99±16.96甘油三酯(mmol/L) 2.14±1.55 1.23±0.87 SDNN(ms) 118.30±42.83 150.04±36.年齡(歲) 62.08±15.10 59.11±14.07

表5 SDNN診斷冠心病參數估計結果
由表5結果可見,參數φ2和φ3的95%可信區間未包括0,可以認為收縮壓和甘油三酯含量是影響冠心病患病與否的因素,兩指標95%可信區間均大于0,表明冠心病患者組的收縮壓和甘油三酯含量均高于對照組;λD,3和 λH,3的 95% 可信區間未包括 0,表明甘油三酯含量對SDNN診斷冠心病有影響,二者均小于0,表明甘油三酯含量越高,對應的SDNN值越小,心率變異越小。同時,尚不能認為協變量年齡和收縮壓是心率變異SDNN診斷冠心病時的有意義因素。
若以甘油三酯含量分層,繪制相應的ROC曲線(見圖1)可見,甘油三酯含量越高,相應的ROC曲線下面積越大;AUC分別為0.7965、0.8359和0.8792,變化較為明顯,表明甘油三酯含量在SDNN診斷冠心病方面,對診斷結果準確性有較大影響。甘油三酯含量偏離正常值越遠,SDNN用于冠心病診斷的準確性就越高。

圖1 不同甘油三酯下SDNN診斷冠心病的ROC曲線
模擬研究結果表明,隨樣本含量的增加,兩部貝葉斯模型對于無金標準條件下考慮協變量后的ROC曲線估計,參數估計的準確度越來越好,當樣本含量增加到150時,參數估計值與真值接近。在樣本含量低于100時,模型估計得到的值與設定的真實值的一致率不足85%;當樣本含量逐漸增大到150時,二者的一致率達90%;當樣本含量大于300時,參數估計結果與真值接近。
24 h動態心電圖(Holter)作為無創檢查很易被大家接受,但一直被臨床工作者認為,不能作為有價值的診斷指標。SDNN是24 h動態心電圖心率變異的一項基礎指標,容易檢測。實例分析結果可見,在考慮協變量對疾病狀態影響的同時,考慮了協變量對診斷試驗結果影響的兩部貝葉斯模型,可作為無金標準條件下,冠心病診斷準確性達到80%左右的無創診斷指標。本研究將影響診斷結果的甘油三酯含量納入模型,得到含量偏離正常值(<1.7mmol/L)越遠,診斷準確性越高的結論,解釋更接近臨床實際,并可計算出不同甘油三酯含量條件下,ROC曲線診斷的準確度。
通常的ROC分析,反應變量若不服從雙正態分布,需要對數據進行轉換,而兩部貝葉斯模型對于反應變量不服從雙正態分布的情況,可以通過改變第二個過程的擬合模型來解決,對于反應變量分布沒有特定要求。
兩部貝葉斯模型作為無金標準診斷試驗模型,不僅可考慮協變量的影響,同時可較好地解決無金標準時疾病的診斷問題。模擬研究表明,基于MCMC算法的貝葉斯估計,參數估計結果準確。SDNN是考慮甘油三酯含量影響后,冠心病無創診斷試驗穩健性更好的一項診斷指標。有關含有缺失值的兩部貝葉斯模型估計,建議采用隨機缺失數據多重填補后,再進行參數估計,討論見另文。基于MCMC算法的貝葉斯估計,退火迭代次數應保證迭代鏈達到收斂,結果可由模擬結果迭代圖來判斷。若不能達到收斂,就不可進行參數估計,需考慮修正模型或者選擇其他形式的先驗分布(見另文)。
1.Zhou XH.Statistical methods in diagnostic medicine.New York:John Wiley & Sons,2002.
2.O'Malley AJ,Zou KH,Fielding JR,et al.Bayesian regression methodology for estimating a receiver operating characteristic curve with two radiologic applications.Academic Radiology,2001,8(8):713-725.
3.Wang C,Turnball BW,Grohn YT,et al.Estimating receiver operating characteristic curves with covariates when there is no perfect reference test for diagnosis of Johne's disease.American Dairy Science Association,2006,89:3038-3046.
4.Goetqhebenr E,liinev J,Boelaert M,et al.Diagnostic test analyses in search of their gold standard:latent class analyses with random effects.Statistical Methods in Medical Research,2000,9(3):231-248.
5.Pepe MS.A regression modelling framework for receiver operating characteristic curves in medical diagnostic testing.Biometrika,1997,84:595-608.
6.宇傳華.ROC分析方法及其在醫學研究中的應用.陜西:第四醫科大學(博士論文),2000.
7.陳衛中,潘小平,倪宗瓚.Logistic回歸模型在ROC分析中的應用.中國衛生統計,2007,24(1):22-24.
Estimating Receiver Operating Characteristic Curves with Co- variates in the Absence of Gold Standard Test
LiuYan,WangHongyu,TianJiaoni,etal.TheDepartmentofHealthStatistics,ShanxiMedicalUniversity(030001),Taiyuan
ObjectiveTo introduce the method for estimation of receiver operating characteristic curves with covariates when there is no gold standard test.MethodsTo estimate the ROC curve after considered the effects of covariates,using two-part Bayesian model screening the impact factors for the ROC curve without gold standard test.ResultsTwopart Bayesian model can detect the impact of covariates not only on disease status but also on test scores.Moreover,the areas under the ROC curve with different values of covariates could be calculated.The method was also treated as the proof of statistical analysis for clinical diagnostic test.ConclusionTwo-part Bayesian model can be effectively used to solve the problem of estimating ROC curves impacted by covariates when there is no gold standard test.
Two-part Bayesian model;ROC curve;No golden standard test;Covariate;