江蘇省計(jì)劃生育科學(xué)技術(shù)研究所(210036)

施雯慧 許豪勤△ 巴 磊

【提 要】 目的 通過(guò)蒙特卡羅模擬比較基于錯(cuò)誤發(fā)現(xiàn)率(FDR)控制的信號(hào)檢測(cè)方法與傳統(tǒng)方法的檢測(cè)效能,驗(yàn)證基于錯(cuò)誤發(fā)現(xiàn)率控制的方法是否能解決傳統(tǒng)方法假陽(yáng)性信號(hào)比例高的問(wèn)題。方法 基于不良反應(yīng)監(jiān)測(cè)數(shù)據(jù)的特點(diǎn),建立自發(fā)呈報(bào)系統(tǒng)模型,合理設(shè)置參數(shù),產(chǎn)生模擬數(shù)據(jù),使用R軟件進(jìn)行數(shù)據(jù)分析,進(jìn)行傳統(tǒng)方法和基于錯(cuò)誤發(fā)現(xiàn)率控制的方法的比較,并考察判定閾值對(duì)信號(hào)檢測(cè)結(jié)果的影響。結(jié)果 傳統(tǒng)頻數(shù)法(PRR和ROR法)的靈敏度高于貝葉斯法(BCPNN和GPS法),特異度稍低;基于FDR控制的方法較傳統(tǒng)方法靈敏度提高,特異度降低,且FDR控制對(duì)貝葉斯法的影響較頻數(shù)法更大。結(jié)論 若檢測(cè)要求較高的靈敏度和較大的ROC曲線下面積,可優(yōu)先選擇基于FDR控制的GPS法,若要求較高的陽(yáng)性預(yù)測(cè)值,則傳統(tǒng)BCPNN法為首選。

藥品不良反應(yīng)/醫(yī)療器械不良事件自發(fā)呈報(bào)系統(tǒng)目前是國(guó)際上藥品/醫(yī)療器械上市后監(jiān)測(cè)的最主要手段[1-2],自20世紀(jì)60年代以來(lái),隨著監(jiān)測(cè)體系和相關(guān)法規(guī)的不斷完善,各國(guó)的自發(fā)呈報(bào)系統(tǒng)收集到的不良反應(yīng)報(bào)告數(shù)大幅增長(zhǎng),以我國(guó)為例,截至2016年國(guó)家藥品不良反應(yīng)監(jiān)測(cè)中心已累計(jì)收到藥品不良反應(yīng)/事件報(bào)告近1075萬(wàn)份[3]、可疑醫(yī)療器械不良事件報(bào)告近168萬(wàn)份[4]。如何從海量報(bào)告數(shù)據(jù)中及時(shí)、準(zhǔn)確地發(fā)現(xiàn)不良反應(yīng)信號(hào),是當(dāng)前藥物警戒研究領(lǐng)域的熱點(diǎn)問(wèn)題。目前常用的信號(hào)檢測(cè)方法主要是基于不相稱測(cè)定原理,可分為：(1)頻數(shù)法,如比例報(bào)告比值比法(proportional reporting ratio,PRR)、報(bào)告比值比法(reporting odds ratio,ROR);(2)貝葉斯法,如貝葉斯置信傳播神經(jīng)網(wǎng)絡(luò)(Bayesian confidence propagation neural network,BCPNN)、經(jīng)驗(yàn)貝葉斯伽馬泊松縮減(empirical Bayes gamma Poisson shrinker,GPS)。這些方法應(yīng)用已經(jīng)較為成熟,但存在多重假設(shè)檢驗(yàn)問(wèn)題[5]。信號(hào)檢測(cè)過(guò)程一般是根據(jù)目標(biāo)藥品/醫(yī)療器械和目標(biāo)不良反應(yīng)/事件,計(jì)算各方法下的不相稱測(cè)定指標(biāo)(如PRR、ROR、IC等),再依據(jù)判定規(guī)則判斷該目標(biāo)組合是否為可疑信號(hào)。假設(shè)數(shù)據(jù)庫(kù)中有i種藥品/醫(yī)療器械和j種不良反應(yīng)/事件,則相應(yīng)的檢驗(yàn)進(jìn)行了i×j次,如將每次比較的錯(cuò)誤率控制在0.05的水準(zhǔn),則會(huì)增大總Ⅰ型錯(cuò)誤率,產(chǎn)生較多的假陽(yáng)性信號(hào)。

傳統(tǒng)的多重比較方法(如Bonferroni法、sidak法等)控制的是總Ⅰ型錯(cuò)誤率[6],如將其定為0.05,那么每次檢驗(yàn)的水準(zhǔn)就極低,則會(huì)導(dǎo)致只能發(fā)現(xiàn)少數(shù)強(qiáng)關(guān)聯(lián)組合,檢驗(yàn)效能降低。FDR(false discovery rate)稱為錯(cuò)誤發(fā)現(xiàn)率或陽(yáng)性結(jié)果錯(cuò)誤率,表示陽(yáng)性檢驗(yàn)結(jié)果中判斷錯(cuò)誤的比例,由Benjaminni和Hochberg首先提出[7]。相比傳統(tǒng)假設(shè)檢驗(yàn)的檢驗(yàn)水準(zhǔn)取值固定,FDR能靈活取值,作為假設(shè)檢驗(yàn)錯(cuò)誤率的控制指標(biāo);此外,相比總Ⅰ型錯(cuò)誤率主要用于控制I類錯(cuò)誤,FDR的意義更為明確,可用于評(píng)價(jià)篩選出來(lái)的差異變量,因而常用于高微陣列數(shù)據(jù)分析的多重比較[8]。

FDR在藥品不良反應(yīng)信號(hào)檢測(cè)領(lǐng)域的應(yīng)用始于法國(guó)的Ahmed[9-11],但在醫(yī)療器械不良事件信號(hào)檢測(cè)方面的效果尚未可知?？紤]到信號(hào)檢測(cè)方法沒(méi)有金標(biāo)準(zhǔn),本文采用蒙特卡羅模擬的方法來(lái)比較四種常用信號(hào)檢測(cè)方法(PRR、ROR、BCPNN、GPS)及基于FDR控制后的四種方法在宮內(nèi)節(jié)育器不良事件數(shù)據(jù)中的檢驗(yàn)效能。

模型構(gòu)建

Emmanuel Roux等的研究表明[12],自發(fā)報(bào)告系統(tǒng)的藥物暴露頻數(shù)服從Poisson分布,報(bào)告數(shù)nij在一定時(shí)間Δt內(nèi)服從δij的Poisson分布:

nij～Poi(δij)=Poi(RRijIjEiPij)

其中RRij為藥物i與不良事件j組合的相對(duì)危險(xiǎn)度,Ij為不良事件j的背景發(fā)生率,Ei為藥物i的使用人數(shù),Pij為藥物i與不良事件j組合的報(bào)告概率。

國(guó)內(nèi)已有的幾項(xiàng)藥品不良反應(yīng)模擬數(shù)據(jù)研究中均沿用了上述模型,假定報(bào)告影響因素包括用藥人數(shù)、藥品上市時(shí)間、不良事件背景發(fā)生率、不良事件嚴(yán)重程度、報(bào)告概率[13-14]。而國(guó)家計(jì)劃生育藥具不良反應(yīng)監(jiān)測(cè)中心既往幾年的監(jiān)測(cè)工作情況表明,宮內(nèi)節(jié)育器不良事件報(bào)告并不完全符合上述特點(diǎn)。首先,嚴(yán)重的不良事件報(bào)告概率未必高于一般不良事件,因?yàn)榇蟛糠直O(jiān)測(cè)點(diǎn)設(shè)立于縣、鄉(xiāng)兩級(jí)的婦幼保健計(jì)劃生育服務(wù)機(jī)構(gòu),能收集到的嚴(yán)重傷害事件不多(由于醫(yī)療水平的限制,大多數(shù)嚴(yán)重傷害事件會(huì)被轉(zhuǎn)診至三級(jí)醫(yī)療機(jī)構(gòu));其次,由于宮內(nèi)節(jié)育器產(chǎn)品的特殊性(使用時(shí)間長(zhǎng),研發(fā)周期長(zhǎng),在市場(chǎng)上流通種類遠(yuǎn)低于藥品),產(chǎn)品上市時(shí)間與報(bào)告概率的關(guān)系也不是很大;因此,參考相關(guān)文獻(xiàn)和專家意見(jiàn),假定宮內(nèi)節(jié)育器不良事件自發(fā)呈報(bào)系統(tǒng)的報(bào)告數(shù)服從Poisson分布,影響因素包括使用人數(shù)、不良事件背景發(fā)生率和報(bào)告概率。

參數(shù)設(shè)置

本次模擬假定有40種宮內(nèi)節(jié)育器和30種不良事件,共計(jì)有1200種節(jié)育器與不良事件組合,假定每種宮內(nèi)節(jié)育器有4種與其有關(guān)聯(lián)的不良事件(相對(duì)危險(xiǎn)度RR從低到高設(shè)置4檔,分別為2、3、5和10),則共計(jì)有40×4=160個(gè)組合為真實(shí)信號(hào),剩余組合為虛假信號(hào),相對(duì)危險(xiǎn)度RR為1。

40種宮內(nèi)節(jié)育器中,按使用人數(shù)分為三類:5種為常用,使用人數(shù)估計(jì)為500萬(wàn);15種為一般,使用人數(shù)估計(jì)為80萬(wàn);20種為較少使用,使用人數(shù)估計(jì)為10萬(wàn)。30種不良事件中,按背景發(fā)生率分為三類,其中10種為1/100,10種為1/1000,10種為1/20000;報(bào)告概率設(shè)定為4個(gè)水平,分別為0.005、0.025、0.05和0.1。

步 驟

1.生成模擬數(shù)據(jù)庫(kù)

使用SAS9.3軟件隨機(jī)分配各節(jié)育器對(duì)應(yīng)的4種有關(guān)聯(lián)的不良事件,按照構(gòu)建的模型和參數(shù)設(shè)置節(jié)育器與不良事件組合發(fā)生的頻數(shù)。模擬產(chǎn)生100個(gè)對(duì)應(yīng)數(shù)據(jù)集,信號(hào)檢測(cè)結(jié)果取100個(gè)數(shù)據(jù)集結(jié)果的均值。

2.數(shù)據(jù)提取及信號(hào)檢測(cè)

提取信號(hào)檢測(cè)計(jì)算時(shí)所需要的四格表數(shù)據(jù),使用R軟件及PhViD包進(jìn)行信號(hào)檢測(cè)。PhViD包由Ismail Ahmed開(kāi)發(fā),使用其內(nèi)置的主要函數(shù)PRR()、ROR()、BCPNN()、GPS(),可根據(jù)給定數(shù)據(jù)以PRR、ROR、BCPNN、GPS四種方法快速計(jì)算出信號(hào)。這些函數(shù)中比較重要的參數(shù)包括MIN.n11、DECISION、DECISION.THRES等,其中MIN.n11指定目標(biāo)藥物與目標(biāo)不良事件組合的最小頻數(shù),DECISION指定生成信號(hào)的判定規(guī)則(基于傳統(tǒng)規(guī)則還是基于FDR控制),DECISION.THRES指定相應(yīng)規(guī)則的閾值。本軟件包中使用的FDR估計(jì)方法是基于LBE的估計(jì)方法[15]。

根據(jù)PhViD包的使用說(shuō)明,基于傳統(tǒng)判定規(guī)則生成信號(hào)時(shí),對(duì)應(yīng)的函數(shù)參數(shù)設(shè)置如下:

(1)PRR(DATABASE,RR0 = 1,DECISION = 3,DECISION.THRES = 0,RANKSTAT = 2)

(2)ROR(DATABASE,RR0 = 1,DECISION = 3,DECISION.THRES = 0,RANKSTAT = 2)

(3)BCPNN(DATABASE,RR0 = 1,DECISION = 3,DECISION.THRES = 0,RANKSTAT = 2,MC = TRUE,NB.MC = 10000)

(4)GPS(DATABASE,RR0 = 1,DECISION = 3,DECISION.THRES =1,RANKSTAT = 2,TRONC = FALSE,TRONC.THRES = 1,PRIOR.INIT = c(alpha1 = 0.2,beta1 = 0.06,alpha2 = 1.4,beta2 = 1.8,w = 0.1),PRIOR.PARAM = NULL)

基于FDR控制生成信號(hào)時(shí),對(duì)應(yīng)的函數(shù)參數(shù)設(shè)置如下:

(1)PRR(DATABASE,RR0 = 1,DECISION.THRES=0.05)

(2)ROR(DATABASE,RR0 = 1,DECISION.THRES=0.05)

(3)BCPNN(DATABASE,RR0 = 1,DECISION.THRES=0.05)

(4)GPS(DATABASE,RR0 = 1,DECISION.THRES=0.05)

3.結(jié)果比較

將各種方法的檢測(cè)結(jié)果與模擬設(shè)置的真實(shí)信號(hào)(RR不為1的組合)進(jìn)行比較,計(jì)算靈敏度、特異度、陽(yáng)性預(yù)測(cè)值、陰性預(yù)測(cè)值、約登指數(shù)等指標(biāo),并比較不同判定閾值A(chǔ)、不同RR對(duì)各方法檢測(cè)結(jié)果的影響。

結(jié) 果

1.模擬數(shù)據(jù)描述

100個(gè)模擬數(shù)據(jù)集共有30792825條記錄,宮內(nèi)節(jié)育器與不良事件組合96312個(gè),平均每個(gè)數(shù)據(jù)集有30余萬(wàn)條記錄,涉及900余個(gè)組合。模擬數(shù)據(jù)中RR的分布如表1所示,各組合記錄數(shù)的分布如表2所示。

表1 模擬數(shù)據(jù)庫(kù)中宮內(nèi)節(jié)育器與不良事件組合RR的分布

表2 模擬數(shù)據(jù)庫(kù)中宮內(nèi)節(jié)育器-不良事件組合記錄數(shù)的分布

2.RR對(duì)各方法靈敏度的影響

考慮到在模擬數(shù)據(jù)庫(kù)參數(shù)設(shè)置中,將RR分為四級(jí),可針對(duì)不同RR分析各檢測(cè)方法的靈敏度,結(jié)果如表3所示。隨RR增大,各方法的靈敏度均明顯增加。在同一RR下,PRR和ROR法的靈敏度非常接近,而B(niǎo)CPNN和GPS法的靈敏度稍低;基于FDR控制的四種方法相較傳統(tǒng)方法靈敏度均有所提高,但BCPNN和GPS法的提高幅度大于PRR和ROR法。

表3 各檢測(cè)方法在不同RR下的靈敏度

3.判定閾值A(chǔ)對(duì)各方法檢測(cè)能力的影響

模擬數(shù)據(jù)庫(kù)中宮內(nèi)節(jié)育器與不良事件組合記錄數(shù)的最小值1,最大值50671,四分位數(shù)分別為4、25、108,選取A分別為1、2、3、4、10、25時(shí),比較各方法的檢測(cè)能力。

對(duì)模擬數(shù)據(jù)庫(kù)中的100個(gè)模擬數(shù)據(jù)集分別計(jì)算,求得不同A值下各方法檢測(cè)出信號(hào)數(shù)的均值?？傮w上,各方法生成的信號(hào)數(shù)隨A值增加而減少(BCPNN、GPS法中A≤3時(shí)例外);基于FDR控制的4種方法比各自對(duì)應(yīng)的傳統(tǒng)方法生成的信號(hào)數(shù)均有所增加,且BCPNN、GPS法增加的幅度更大;相同A值下,PRR、ROR法生成的信號(hào)數(shù)近似,基于FDR控制的PRR、ROR法生成信號(hào)數(shù)也近似,但基于FDR控制的BCPNN和GPS法之間信號(hào)數(shù)的差距較傳統(tǒng)BCPNN和GPS法間的差距縮小。

進(jìn)一步比較判定閾值A(chǔ)對(duì)各方法檢測(cè)能力的影響,結(jié)果表明,隨A增大,各方法的靈敏度降低、特異度升高、曲線下面積減少;相同A值下,傳統(tǒng)頻數(shù)法的靈敏度高于傳統(tǒng)貝葉斯法,基于FDR控制的頻數(shù)法靈敏度低于FDR控制的貝葉斯法;四種傳統(tǒng)方法經(jīng)FDR控制后,靈敏度均升高,特異度均降低,且FDR控制對(duì)貝葉斯方法的影響較頻數(shù)法更大,見(jiàn)表4。

討 論

從模擬結(jié)果可以看出,頻數(shù)法(PRR和ROR法)之間的一致性較高,其靈敏度高于貝葉斯法(BCPNN和GPS法),特異度較之稍低,這與既往研究結(jié)論類似。而基于FDR控制的四種方法,檢測(cè)出的信號(hào)數(shù)較傳統(tǒng)方法增加,靈敏度提高,特異度降低,且基于FDR控制的貝葉斯法靈敏度高于基于FDR控制的頻數(shù)法,說(shuō)明FDR控制對(duì)貝葉斯法的影響較頻數(shù)法更大。這一點(diǎn)與采用FDR控制后檢測(cè)信號(hào)數(shù)應(yīng)適當(dāng)減少的預(yù)期不太一致。國(guó)內(nèi)郭曉晶在相關(guān)預(yù)實(shí)驗(yàn)中也發(fā)現(xiàn)[16],將基于LBE的FDR估計(jì)方法與BCPNN法聯(lián)用,產(chǎn)生的信號(hào)數(shù)量多于單獨(dú)使用BCPNN法,具體原因有待進(jìn)一步研究。

表4 各檢測(cè)方法在不同A值下的檢測(cè)能力

而在不同A值下比較各方法檢測(cè)的靈敏度、特異度、曲線下面積等,可發(fā)現(xiàn)A=1時(shí)靈敏度高,曲線下面積較大,這一結(jié)果與既往的藥品不良反應(yīng)數(shù)據(jù)模擬研究并不相同。李嬋娟等[13]的研究表明,選擇A≥3的組合進(jìn)行信號(hào)檢測(cè),靈敏度要高于考慮任何A值的檢測(cè)結(jié)果,并由此得出兩份以下的報(bào)告對(duì)檢測(cè)意義不大的結(jié)論。分析原因,可能是節(jié)育器不良事件數(shù)據(jù)和藥品不良反應(yīng)數(shù)據(jù)的特點(diǎn)不同造成的。藥品不良反應(yīng)報(bào)告數(shù)據(jù)庫(kù)中,往往是報(bào)告例數(shù)小的組合占比較大,如我國(guó)2010-2011年藥品不良反應(yīng)報(bào)告數(shù)據(jù)庫(kù)中報(bào)告數(shù)小于3例的組合占64.35%[14],廣東省2002-2007年的藥品不良反應(yīng)報(bào)告數(shù)據(jù)庫(kù)中報(bào)告數(shù)小于3例的組合占76.24%[17],而2011-2015年宮內(nèi)節(jié)育器不良事件數(shù)據(jù)中報(bào)告數(shù)小于3例的組合僅占25.03%。

綜上所述,信號(hào)檢測(cè)方法應(yīng)用于宮內(nèi)節(jié)育器不良事件數(shù)據(jù)中時(shí),沒(méi)有必要限定于A≥3的組合。若檢測(cè)目的是要求較高的靈敏度和ROC曲線下面積,可優(yōu)先選擇基于FDR控制的GPS法,若重點(diǎn)關(guān)注檢測(cè)出的信號(hào)的可靠性,即要求較高的陽(yáng)性預(yù)測(cè)值,則傳統(tǒng)BCPNN法為首選。

模擬研究存在其自身局限性,真實(shí)世界中不良反應(yīng)/事件報(bào)告的影響因素眾多,不良事件背景發(fā)生率、報(bào)告概率等都未可知,很難找到完全理想的模型來(lái)構(gòu)建。本次模擬的數(shù)據(jù)與真實(shí)的宮內(nèi)節(jié)育器不良事件報(bào)告數(shù)據(jù)相比,節(jié)育器與不良事件的組合數(shù)稍多,各組合涉及報(bào)告數(shù)的分布較為接近,結(jié)果可供相關(guān)研究人員參考。后續(xù)可更詳盡地考慮其他影響因素,設(shè)置合理的模型和參數(shù),開(kāi)展更為充分的研究。