張?zhí)灬?/p>
系統(tǒng)評(píng)價(jià)和Meta分析是循證醫(yī)學(xué)重要的證據(jù)信息來源[1,2],通過Meta分析對(duì)多個(gè)研究結(jié)果的統(tǒng)計(jì)綜合被廣泛用于評(píng)價(jià)不同干預(yù)的相對(duì)效應(yīng)[3]。在牛津循證醫(yī)學(xué)中心2011年發(fā)布的證據(jù)水平分級(jí)標(biāo)準(zhǔn)中,將以RCT為基礎(chǔ)的系統(tǒng)評(píng)價(jià)作為評(píng)價(jià)某種干預(yù)措施受益或危害的1級(jí)水平[4]。雖然RCT是公認(rèn)的最佳治療性研究設(shè)計(jì)方案,但來自RCT的證據(jù)常缺乏或不足,難以回答所有臨床問題,或有時(shí)RCT不需要、不可行或不恰當(dāng),如新藥1期臨床試驗(yàn),則可以采用單臂試驗(yàn)設(shè)計(jì)方法[5,6],國、內(nèi)外均有依據(jù)單臂試驗(yàn)批準(zhǔn)新藥上市的案例[7]。所謂單臂試驗(yàn)是指將具有某種靶健康狀況的個(gè)體樣本給予試驗(yàn)性治療并觀察某種治療效應(yīng)的設(shè)計(jì)[6],所有受試對(duì)象都在同一觀察組,不設(shè)立其他實(shí)驗(yàn)組或?qū)φ战M[5]。近年來,國內(nèi)外循證醫(yī)學(xué)工作者開始探索單臂試驗(yàn)數(shù)據(jù)的Meta分析方法。單臂試驗(yàn)的觀察指標(biāo)按資料性質(zhì)可分為定量與定性資料。對(duì)于定性資料,多以比例、率、比數(shù)和發(fā)病密度等為效應(yīng)量[5],目前有眾多研究文獻(xiàn)介紹此類數(shù)據(jù)在貝葉斯和頻率學(xué)框架下的Meta分析方法,以及基于R、Stata、SAS、RevMan等軟件實(shí)現(xiàn)[5, 8-13];而對(duì)連續(xù)型數(shù)據(jù)的Meta分析方法則研究或介紹不多,有少數(shù)介紹使用Stata軟件[5]和Meta-Analyst軟件[14]基于頻率學(xué)框架下采用倒方差法進(jìn)行Meta分析的文獻(xiàn),難以滿足實(shí)踐需要。本文介紹基于經(jīng)典的正態(tài)-正態(tài)層次模型(NNHM)框架下單臂試驗(yàn)連續(xù)型數(shù)據(jù)的貝葉斯Meta分析方法,并通過實(shí)例來介紹R軟件實(shí)現(xiàn)過程。
1.1 研究數(shù)據(jù) 數(shù)據(jù)來源于2個(gè)Meta分析。數(shù)據(jù)1:來自Dai等的單臂試驗(yàn)的Meta分析文獻(xiàn)[15],主要觀察來氟米特治療銀屑病性關(guān)節(jié)炎的有效性和安全性,本文選取的觀察指標(biāo)為銀屑病皮損面積和嚴(yán)重度指數(shù)(PASI)評(píng)分,為治療后評(píng)分與基線評(píng)分的變化值,具體數(shù)據(jù)信息為每個(gè)試驗(yàn)的樣本量、評(píng)分均數(shù)及標(biāo)準(zhǔn)差,在表1中分別以“n”、“mean”和“sd”表示;數(shù)據(jù)2:來自Nagpal等的Meta分析文獻(xiàn)[16],主要觀察造血干細(xì)胞療法對(duì)卒中早期臨床干預(yù)的有效性和安全性,該文獻(xiàn)納入了單臂試驗(yàn)和對(duì)照試驗(yàn)的數(shù)據(jù),本文選取單臂試驗(yàn)的相關(guān)數(shù)據(jù),觀察指標(biāo)為美國國立衛(wèi)生研究院卒中量表(NIHSS)評(píng)分,為治療后評(píng)分與基線評(píng)分的變化值,具體數(shù)據(jù)信息為評(píng)分的均數(shù)及其95%CI下限和上限,在表1中分別以“mean”、“l(fā)l”、“ul”表示。數(shù)據(jù)見表1。
表1 分析數(shù)據(jù)
1.2 模型與方法 經(jīng)典的Meta分析隨機(jī)效應(yīng)模型基于NNHM框架[17]:假設(shè)納入Meta分析有i=1,…,N個(gè)獨(dú)立研究,每個(gè)研究的真實(shí)值為θi,其估計(jì)值及估計(jì)值標(biāo)準(zhǔn)誤分別為yi和Si,研究間真實(shí)效應(yīng)變異為τ2,總的合并效應(yīng)為μ,則模型有兩層。
yi~N(θi,Si2)為抽樣模型,假定yi服從未知均值θi和已知標(biāo)準(zhǔn)誤Si的正態(tài)分布;θi~ N(μ,τ2)為參數(shù)模型,假定θi服從均數(shù)為μ、方差為τ2(標(biāo)準(zhǔn)差為τ)的正態(tài)分布,其中未知參數(shù)μ為主要感興趣的參數(shù),而異質(zhì)性參數(shù)τ則一般被認(rèn)為是討厭參數(shù)。在頻率學(xué)框架下,該模型可以使用Stata、R、SAS、RevMan等軟件來擬合和實(shí)現(xiàn),常用的參數(shù)估計(jì)方法有最大似然法(ML)和約束最大似然法(REML)等。在貝葉斯框架下,通常需要根據(jù)實(shí)際情況為參數(shù)μ和τ指定合理的先驗(yàn)。一般情況下,對(duì)于參數(shù)μ,指定為無信息先驗(yàn)分布;而對(duì)τ,可以指定為無信息先驗(yàn)分布如均勻分布,或弱信息先驗(yàn)分布如半正態(tài)分布等,該模型可由貝葉斯軟件WinBugs等軟件擬合。
1.3 模型擬合 在貝葉斯框架下,NNHM可由R軟件(3.4.4版)的bayesmeta包的bayesmeta()函數(shù)來擬合[18]。在本研究中,參數(shù)μ設(shè)定為無信息先驗(yàn)μ~N(0,104);參數(shù)τ指定為弱信息先驗(yàn),分別為半正態(tài)(HN)先驗(yàn)和半柯西(HC)先驗(yàn),尺度均設(shè)為1。以forestplot()函數(shù)繪制森林圖。同時(shí),為了比較與頻率學(xué)框架下Meta分析結(jié)果,以metafor包的rma.uni()函數(shù)擬合該模型,選用REML法估計(jì)參數(shù)τ。具體代碼及簡(jiǎn)要解釋見表2。
表2 單臂試驗(yàn)連續(xù)型數(shù)據(jù)貝葉斯Meta分析代碼(R軟件)及簡(jiǎn)要解釋
續(xù)表2 tau.prior=function(x){(dhalfcauchy(x, scale=1))})print(bma22)##繪制森林圖forestplot(bma22)##擬合模型3:頻率學(xué)框架,隨機(jī)效應(yīng)模型,采用REML法估計(jì)參數(shù)τ。meta23<-rma.uni(yi=score,sei=se,method="REML",data=stemcell)summary(meta23)
2個(gè)研究數(shù)據(jù)貝葉斯和頻率學(xué)Meta分析方法對(duì)效應(yīng)參數(shù)μ和研究間異質(zhì)性τ的估計(jì)結(jié)果如表3所示,2種方法獲得參數(shù)μ點(diǎn)估計(jì)值有差異,95%CI寬度有明顯的不同;參數(shù)τ估計(jì)值差異較大,基于HN(1)先驗(yàn)的貝葉斯方法似乎可以獲得更為精確的估計(jì)值。需要指出的是,bayesmeta包提供了更多的貝葉斯分析結(jié)果信息,如在其“marginal posterior summary”結(jié)果顯示部分,提供2個(gè)主要參數(shù)μ和τ后驗(yàn)眾數(shù)、后驗(yàn)均數(shù)、后驗(yàn)中位數(shù),及相應(yīng)95%CI的下限和上限,為了與森林圖結(jié)果一致,表1中所示選取的是后驗(yàn)中位數(shù)及95% CI的下限和上限結(jié)果。
數(shù)據(jù)1的森林圖如圖1和2所示。森林圖中上半部分顯示了納入Meta分析的每個(gè)研究的效應(yīng)量(即模型中的yi)的點(diǎn)估計(jì)及基于相應(yīng)標(biāo)準(zhǔn)誤獲得的95%CI,以黑色表示,即森林圖中的“quoted estimate”;以及每個(gè)研究真實(shí)效應(yīng)參數(shù)(即模型中的θ_i)的點(diǎn)估計(jì)及95%CI,以灰色表示,即森林圖中的“shrinkage estimate”;森林圖底部顯示合并效應(yīng)量(平均值)點(diǎn)估計(jì)(后驗(yàn)中位數(shù))和95%CI(與圖中“mean”相應(yīng)),以及預(yù)測(cè)值及區(qū)間(與圖中“prediction”相應(yīng))。
表3 基于NNHM框架下貝葉斯和頻率學(xué)Meta分析結(jié)果
圖1 數(shù)據(jù)1森林圖[τ為HN(1)先驗(yàn)]圖2 數(shù)據(jù)1森林圖[τ為HC(1)先驗(yàn)]
單臂試驗(yàn)在醫(yī)學(xué)研究中較為常見,如果其測(cè)量結(jié)局為連續(xù)型數(shù)據(jù),則具有一個(gè)共同特征[19]:針對(duì)某一測(cè)量結(jié)局指標(biāo),除了干預(yù)實(shí)施后進(jìn)行測(cè)量外,可能會(huì)對(duì)每一個(gè)參與者在干預(yù)實(shí)施前進(jìn)行基線測(cè)量,因此可以獲得最終測(cè)量值與基線測(cè)量值的差值,稱為基線改變值或改變?cè)u(píng)分,所以系統(tǒng)評(píng)價(jià)者在Meta分析提取數(shù)據(jù)時(shí)可能會(huì)考慮同時(shí)提取基線改變值和最終測(cè)量值。對(duì)于連續(xù)型數(shù)據(jù),原始研究對(duì)測(cè)量結(jié)果可能存在多種報(bào)告形式,常見的有:①測(cè)量指標(biāo)均數(shù)及相應(yīng)標(biāo)準(zhǔn)誤(或方差);②受試者樣本量、測(cè)量指標(biāo)均數(shù)及相應(yīng)標(biāo)準(zhǔn)差;③測(cè)量指標(biāo)的均數(shù)及95% CI;④測(cè)量指標(biāo)的中位數(shù)及上下4分位數(shù)間距等;⑤測(cè)量指標(biāo)的中位數(shù)、最大值及最小值等。無論何種形式的數(shù)據(jù),只要最終能獲得每個(gè)研究的效應(yīng)量及其方差或標(biāo)準(zhǔn)誤,即可進(jìn)行Meta分析,一般可以通過2種途徑:一是向原始文獻(xiàn)的通訊作者發(fā)郵件索要原始數(shù)據(jù);二是通過數(shù)據(jù)轉(zhuǎn)換,具體方法可以參看國內(nèi)外相關(guān)文獻(xiàn)[19-23]。本文數(shù)據(jù)1已知樣本量、測(cè)量指標(biāo)均數(shù)及相應(yīng)標(biāo)準(zhǔn)差,可以通過公式“均數(shù)標(biāo)準(zhǔn)誤=標(biāo)準(zhǔn)差/樣本量的平方根”獲得均數(shù)標(biāo)準(zhǔn)誤[24];數(shù)據(jù)2(因文獻(xiàn)[16]未報(bào)告每個(gè)研究的樣本量,本研究聚焦在貝葉斯Meta分析方法介紹,簡(jiǎn)單處理為假定每個(gè)研究樣本量很大,95%CI來自正態(tài)分布)已知測(cè)量指標(biāo)95%CI,則可以通過公式“均數(shù)標(biāo)準(zhǔn)誤=(95%CI上限-95%CI下限)/(2×1.96)”來獲得均數(shù)標(biāo)準(zhǔn)誤[19,22],如果研究樣本量小,或者95%CI由t分布計(jì)算所得,則上述公式中的1.96則需要一個(gè)更大一點(diǎn)t分布值代替[19,22],該值可以用自由度(研究樣本量-1)查t分布表或用Stata、R、Microsof Excel等軟件計(jì)算獲得。
一旦獲得每個(gè)研究觀察指標(biāo)測(cè)量值及其相應(yīng)標(biāo)準(zhǔn)誤,則可以在NNHM框架下采用隨機(jī)效應(yīng)模型來進(jìn)行Meta分析,主要有兩個(gè)參數(shù)真實(shí)效應(yīng)(real-valued effect)μ和異質(zhì)性(heterogeneity)τ,如果異質(zhì)性為0,則退變?yōu)楣潭ㄐ?yīng)模型。在NNHM框架下的分析方法有頻率學(xué)和貝葉斯兩大策略,目前更為常用的頻率學(xué)策略一般分2步,第1步估計(jì)異質(zhì)性,第2步基于異質(zhì)性估計(jì)派生出效應(yīng)估計(jì),但目前有很多估計(jì)異質(zhì)性算法,共同的問題是對(duì)于異質(zhì)性方差分量估計(jì)易趨向于0,導(dǎo)致選擇固定效應(yīng)模型。貝葉斯策略近年來也已逐步應(yīng)用于Meta分析中,綜合了未知參數(shù)的部體信息、樣本信息及先驗(yàn)信息,特別是可以允許(實(shí)際上也需要)指定先驗(yàn)信息用于數(shù)據(jù)計(jì)算。在貝葉斯隨機(jī)效應(yīng)模型中,根據(jù)需要很容易對(duì)效應(yīng)參數(shù)μ指定為無信息先驗(yàn)或弱信息先驗(yàn);對(duì)研究間異質(zhì)性參數(shù)τ,則常需要根據(jù)效應(yīng)量、研究數(shù)量等因素綜合考慮指定為信息先驗(yàn),特別是僅有幾個(gè)研究納入Meta分析時(shí)。有學(xué)者認(rèn)為,除非在研究數(shù)量少或期望信息先驗(yàn)等情況下,一般均宜將τ指定為均勻分布;研究數(shù)量少時(shí),推薦τ指定為重尾分布,如半學(xué)生t(half-Student-t)分布、HC分布、HN分布等弱信息先驗(yàn)[25,26]。本研究中所用數(shù)據(jù)納入Meta分析的研究數(shù)量均較少,因此將τ的先驗(yàn)指定為HN(1),另指定HC(1)作為敏感性分析,以比較不同先驗(yàn)對(duì)參數(shù)估計(jì)的影響。本研究發(fā)現(xiàn),采用貝葉斯分析方法,納入Meta分析的研究數(shù)量較少(如3個(gè)研究時(shí)),給定HN(1)和HC(1)先驗(yàn),未知參數(shù)μ和τ的后驗(yàn)分布差異較大;而研究數(shù)量為中等(如7個(gè)時(shí)),兩個(gè)未知參數(shù)后驗(yàn)分布則較為接近,提示,Meta分析納入研究數(shù)量少時(shí),對(duì)于研究間異質(zhì)性估計(jì)要謹(jǐn)慎,值得進(jìn)一步研究和探索。
貝葉斯方法通常需要高維積分方法等,在計(jì)算方面較難實(shí)現(xiàn),但隨著現(xiàn)代計(jì)算機(jī)技術(shù)的出現(xiàn)和發(fā)展,特別是馬爾可夫鏈蒙特卡洛(MCMC)方法的建立,使得計(jì)算方法制約問題得以解決,但使用MCMC方法需要一定量的專業(yè)知識(shí)如收斂性診斷和實(shí)踐經(jīng)驗(yàn)等。因在NNHM框架貝葉斯隨機(jī)效應(yīng)模型中,只有μ和τ兩個(gè)未知參數(shù)需要估計(jì)和推斷,可以利用數(shù)值積分方法或重要性抽樣使計(jì)算簡(jiǎn)單化[27],研發(fā)bayesmeta包的目的就是為了簡(jiǎn)單易行地實(shí)現(xiàn)基于NNHM框架下完全貝葉斯分析策略,采用稱為散度限制性條件鑲嵌(DIRECT)策略[28]的新算法獲得未知參數(shù)的后驗(yàn)分布,使用者無需擔(dān)心MCMC方法相關(guān)設(shè)置、診斷和后加工處理等貝葉斯方法專業(yè)問題。使用bayesmeta包可以快速獲得精確的后驗(yàn)概括性結(jié)果,并能提供類似于經(jīng)典頻率學(xué)Meta分析森林圖的圖示化結(jié)果;后驗(yàn)分布可以準(zhǔn)解析形式或預(yù)測(cè)估計(jì)和收縮估計(jì)(shrinkage estimation)等高級(jí)方法如存取。bayesmeta包還提供了其他結(jié)果信息,如預(yù)測(cè)未來一個(gè)新的研究真實(shí)效應(yīng)參數(shù)θk+1值,即森林圖中“prediction”相應(yīng)部分,也可以通過summary()函數(shù)獲得具體數(shù)字化結(jié)果。
最后,因?yàn)锽ayesmeta包具有以下優(yōu)點(diǎn):①計(jì)算速度快,精確度和可重復(fù)性強(qiáng);②操作簡(jiǎn)單,不需要深度掌握MCMC算法等專業(yè)知識(shí);③允許快速敏感性檢驗(yàn)、便捷地大規(guī)模模擬等,所以值得推廣和使用。該軟件包除了用于單臂試驗(yàn)連續(xù)型數(shù)據(jù)Meta分析外,還可以用于RCT等雙臂設(shè)計(jì)的二分類數(shù)據(jù)、連續(xù)型數(shù)據(jù)、生存數(shù)據(jù)等貝葉斯Meta分析。