復(fù)旦大學(xué)公共衛(wèi)生學(xué)院生物統(tǒng)計(jì)學(xué)教研室和公共衛(wèi)生安全教育部重點(diǎn)實(shí)驗(yàn)室(200032)

孫 婷 秦國友 武振宇△ 趙耐青

?

不同混雜結(jié)構(gòu)條件下各傾向性評分方法的模擬比較研究*

復(fù)旦大學(xué)公共衛(wèi)生學(xué)院生物統(tǒng)計(jì)學(xué)教研室和公共衛(wèi)生安全教育部重點(diǎn)實(shí)驗(yàn)室(200032)

孫 婷 秦國友 武振宇△趙耐青

目的 通過構(gòu)建不同混雜結(jié)構(gòu)的處理因素模型和結(jié)局模型、不同相關(guān)性的協(xié)變量，比較多種傾向性評分方法在結(jié)局模型為線性回歸模型的情況下估計(jì)處理效應(yīng)的優(yōu)劣。方法 采用Monte Carlo模擬方法，通過構(gòu)建四種由簡單到復(fù)雜的不同結(jié)構(gòu)的混雜模型，生成相應(yīng)的數(shù)據(jù)集，再分別應(yīng)用傾向性評分匹配、回歸調(diào)整、加權(quán)以及分層的方法估計(jì)處理效應(yīng)并進(jìn)行比較。評價(jià)指標(biāo)包括點(diǎn)估計(jì)、標(biāo)準(zhǔn)誤、相對偏倚、均方誤差。結(jié)果 在結(jié)局模型為線性回歸模型情況下，傾向性評分回歸調(diào)整法估計(jì)的相對偏倚最小，穩(wěn)定性也最好。匹配法卡鉗值取0.02較卡鉗值取傾向性評分標(biāo)準(zhǔn)差的0.2倍估計(jì)的相對偏倚更小。當(dāng)處理因素模型中含有非線性效應(yīng)時(shí)，用逆概率加權(quán)法估計(jì)的偏倚較大，并且加權(quán)法估計(jì)的標(biāo)準(zhǔn)誤也最大。傾向性評分分層法在各種情況下估計(jì)的相對偏倚都較大。結(jié)論 傾向性評分回歸調(diào)整法能夠較好地估計(jì)處理效應(yīng)，并且在各種情況下估計(jì)都較為穩(wěn)健。建議當(dāng)協(xié)變量與處理因素和結(jié)局變量的關(guān)系無法確定時(shí)，這四種方法中可以考慮優(yōu)先使用回歸調(diào)整法。

傾向性評分 混雜因素 Monte Carlo模擬 偏倚

隨機(jī)對照試驗(yàn)(randomized controlled trial,RCT)被認(rèn)為是評價(jià)干預(yù)作用設(shè)計(jì)的金標(biāo)準(zhǔn)，但由于RCT其自身的限制，在很多情況下無法實(shí)施[1]。而觀察性研究或回顧性電子醫(yī)療病歷的數(shù)據(jù)相對較易獲得，因此越來越多地被應(yīng)用于不同干預(yù)效果的比較、食品和藥物安全研究、藥物不良反應(yīng)監(jiān)測等。在觀察性研究中，受試對象的治療分組往往不是隨機(jī)分配的，因此會存在某些重要混雜因素在不同組間分布不均衡，從而無法準(zhǔn)確評價(jià)干預(yù)因素對結(jié)局的作用[2]。在不同調(diào)整混雜的方法中，由Rosenbaum和Rubin提出的傾向性評分(propensity score,PS)方法逐漸受到重視[3]，傾向性評分是指在給定一組協(xié)變量的情況下，任意一個(gè)研究對象被分到處理組的條件概率，計(jì)算PS時(shí)不存在未觀察到的混雜因素的假設(shè)條件下，基于PS的匹配法、分層法、回歸調(diào)整法均能得到處理效應(yīng)的相合估計(jì)[3]。

目前傾向性評分主要有傾向性評分匹配法、回歸調(diào)整法、分層法和加權(quán)法四種應(yīng)用方式。有文獻(xiàn)總結(jié)了2001-2009年高影響力的醫(yī)學(xué)和公共衛(wèi)生雜志中使用傾向性評分的文獻(xiàn)情況，其中匹配法是使用最多的方法，占所有查閱文獻(xiàn)的三分之一，其次是回歸調(diào)整法，約占24%，分層法和加權(quán)法分別占22%和18%[4]。評價(jià)各種傾向性評分方法估計(jì)效果的文獻(xiàn)有很多，我們可以看到在估計(jì)不同指標(biāo)的情況下各種估計(jì)方法的優(yōu)劣有所不同[5-9]。在線性回歸模型情況下，有文獻(xiàn)指出PS回歸調(diào)整法可以得出相合估計(jì)[10]。另外，當(dāng)結(jié)局變量為連續(xù)性變量時(shí)，基于不同組別的傾向性評分重疊情況，幾種不同的PS方法估計(jì)的優(yōu)劣情況也不同[4]，當(dāng)處理組傾向性評分分布被包含在對照組中時(shí)，匹配法估計(jì)偏倚較小，但如果處理組和對照組傾向性評分分布重疊很小時(shí)，回歸調(diào)整法表現(xiàn)更好。在醫(yī)學(xué)研究中，有很多情況下結(jié)局為連續(xù)性變量，并且協(xié)變量之間存在不同情況的相關(guān)性，當(dāng)協(xié)變量與處理因素和結(jié)局變量之間的混雜結(jié)構(gòu)不同時(shí)，例如協(xié)變量與處理因素或結(jié)局變量之間不僅含有線性關(guān)系，還有非線性關(guān)系時(shí)，這四種基于PS的方法估計(jì)處理效應(yīng)的偏倚大小和穩(wěn)定情況還沒有明確的結(jié)論。因此，本研究通過設(shè)置不同混雜結(jié)構(gòu)的模型，在連續(xù)性結(jié)局變量、不同相關(guān)程度自變量的情況下，比較各種傾向性評分方法估計(jì)真實(shí)處理效應(yīng)優(yōu)劣的情況，為不同條件下觀察性數(shù)據(jù)的分析處理提供依據(jù)。

原理與方法

傾向性評分是指在給定一組協(xié)變量(Xi)條件下，將任意一個(gè)研究對象i(i=1,2,…,n)分配到處理組(Zi=1)的條件概率。第i個(gè)研究對象被分配到處理組的條件概率可以表示為：

e(Xi)=P(Zi=1|Xi)

(1)

其中，e(Xi)被稱為傾向性評分。傾向性評分相同的兩個(gè)不同組別的研究對象，其擁有的多個(gè)協(xié)變量整體上分布是相同的[11]。因此，組間協(xié)變量的不均衡性對處理效應(yīng)估計(jì)的干擾被消除了。

傾向性評分的基本原理是用一個(gè)分值來替代多個(gè)協(xié)變量，均衡處理組和對照組間協(xié)變量的分布。對非隨機(jī)化研究中的混雜因素進(jìn)行類似隨機(jī)化的均衡處理，減少選擇偏倚。計(jì)算得出PS分值后，可采用匹配、回歸調(diào)整、加權(quán)、分層的方法來均衡各組間協(xié)變量的差異，最終估計(jì)處理效應(yīng)。

1.傾向性評分匹配法

傾向性評分匹配是從對照組中選出與處理組中某一個(gè)體傾向性評分值相同或相近的個(gè)體進(jìn)行配對，常用的匹配方法有最鄰近匹配、卡鉗匹配、全局最優(yōu)匹配等。本研究采用的是最鄰近卡鉗匹配法，即將處理組和對照組的研究對象分別進(jìn)行隨機(jī)排序，然后從處理組中依次選出一個(gè)研究對象，從對照組中選出與其最接近的傾向性評分值的研究對象進(jìn)行匹配。配對時(shí)設(shè)置兩組傾向性評分的差值在一定范圍內(nèi)，即卡鉗值(caliper)。根據(jù)以往研究建議[12]，本研究設(shè)置兩種卡鉗值，分別為0.02和兩組傾向性評分標(biāo)準(zhǔn)差的0.2倍。在匹配的過程中一個(gè)關(guān)鍵的問題是是否允許放回。有無放回是指配對后的對照組對象是否參加下一組的配對。允許放回在方差估計(jì)的過程中需要考慮某一研究對象被多次使用的事實(shí)[13]。本研究統(tǒng)一采用無放回匹配。

2.傾向性評分回歸調(diào)整法

在四種基于傾向性評分的方法中，回歸調(diào)整法是最直接、使用最為方便的一種方法，因此在臨床醫(yī)學(xué)研究中使用也較為廣泛[4]。傾向性評分回歸調(diào)整法是指將估計(jì)的傾向性評分值作為一個(gè)協(xié)變量與處理因素一起納入到估計(jì)處理效應(yīng)的回歸模型中。對于本研究中連續(xù)性的結(jié)局變量使用線性回歸方法，處理效應(yīng)的估計(jì)是調(diào)整的均值差。盡管在某些條件下使用PS方法可以得到相應(yīng)估計(jì)，但是當(dāng)結(jié)局變量和PS值之間的線性關(guān)系不成立時(shí)，用這種方法估計(jì)的處理效應(yīng)可能是有偏的[4]。

3.傾向性評分逆概率加權(quán)法(Inverse probability of weighting，IPW)

傾向性評分逆概率加權(quán)，是邊緣結(jié)構(gòu)模型這類因果推斷方法中的一種[14]，其基本原理與傳統(tǒng)的標(biāo)準(zhǔn)化法類似。根據(jù)傾向性評分值賦予每個(gè)研究對象一個(gè)相應(yīng)的權(quán)重，從而構(gòu)建出一個(gè)虛擬的人群，在這個(gè)虛擬人群中，協(xié)變量的組間分布沒有差異，因此消除了混雜因素的影響。在逆概率加權(quán)的方法中，權(quán)重被定義為研究對象實(shí)際分組情況的概率的倒數(shù)，計(jì)算如下：

(2)

計(jì)算權(quán)重后，再應(yīng)用加權(quán)回歸的方法估計(jì)處理效應(yīng)。

4.傾向性評分分層法(PS stratification)

分層法是非隨機(jī)化研究中控制偏倚的重要手段。傾向性評分分層法是指在估計(jì)出每個(gè)研究對象的傾向性評分值后，根據(jù)傾向性評分值將研究對象分為若干層。文獻(xiàn)指出，當(dāng)估計(jì)線性處理效應(yīng)的時(shí)候，將傾向性評分值分為五層可以消除組間近90%的混雜偏倚[15]。經(jīng)過分層后，每一層內(nèi)處理組與對照組的協(xié)變量分布應(yīng)該是均衡的。分析過程中，先在每一層內(nèi)估計(jì)處理效應(yīng)，最后將每層的效應(yīng)整合成總的處理效應(yīng)。

模擬實(shí)驗(yàn)

在之前研究的基礎(chǔ)上[16]，本研究構(gòu)建了四種由簡單到復(fù)雜，不同結(jié)構(gòu)的混雜模型，并且考慮協(xié)變量之間不同的相關(guān)性，分別使用傾向性評分回歸調(diào)整、匹配、加權(quán)以及分層的方法估計(jì)處理效應(yīng)并通過Monte Carlo模擬方法進(jìn)行比較。

1.數(shù)據(jù)集生成和參數(shù)設(shè)置

(1)生成自變量

本研究共模擬20個(gè)協(xié)變量，其中X1～X10設(shè)置為混雜因素，這10個(gè)協(xié)變量的產(chǎn)生我們設(shè)置了三種不同的情形：

情形I:X1～X5為服從N(0,1)標(biāo)準(zhǔn)正態(tài)分布的連續(xù)性變量，并且自變量之間相互獨(dú)立；X6～X10為服從Bern(0.5)分布的二分類變量，且相互獨(dú)立。

情形II:X1～X5為服從N(0,1)標(biāo)準(zhǔn)正態(tài)分布的連續(xù)性變量，X1和X2具有相關(guān)性，相關(guān)系數(shù)R12=0.2；X6～X10為服從Bern(0.5)分布的二分類變量，且相互獨(dú)立。

情形III:X1～X5為服從N(0,1)標(biāo)準(zhǔn)正態(tài)分布的連續(xù)性變量，X1和X2、X1和X3、X2和X4分別具有相關(guān)性，相關(guān)系數(shù)分別為R12=0.2，R13=0.3，R24=0.4；X6～X10為服從Bern(0.5)分布的二分類變量，且相互獨(dú)立。

為了模擬實(shí)際情況，又生成另外10個(gè)相互獨(dú)立的噪聲協(xié)變量(X11～X20)，與處理因素和結(jié)局變量均無關(guān)。其中X11～X15為服從N(0,1)標(biāo)準(zhǔn)正態(tài)分布的連續(xù)性變量，X16～X20為服從Bern(0.5)分布的二分類變量。

(2)生成處理因素和結(jié)局變量

為了模擬不同復(fù)雜程度的結(jié)構(gòu)，我們設(shè)置了四種混雜模型來生成處理因素和結(jié)局變量。其中，簡單混雜結(jié)構(gòu)中只有線性關(guān)系，而復(fù)雜混雜結(jié)構(gòu)分別是處理因素模型中含有非線性關(guān)系、結(jié)局變量模型中含有非線性關(guān)系以及兩個(gè)模型中均含有非線性關(guān)系。

簡單混雜結(jié)構(gòu)：在這種簡單線性情況下，處理因素由公式(3)產(chǎn)生，結(jié)局變量由線性回歸模型公式(4)產(chǎn)生。其中α、β系數(shù)的值均從Unif(-1,1)中隨機(jī)產(chǎn)生，并且在之后的模擬中固定，αZ為處理因素效應(yīng)，設(shè)為0.5，i～N(0,1)代表隨機(jī)測量誤差。

logit{ei(Xi;θ1)}=β0+β1X1,i+…+β10X10,i

(3)

yi=αzZi+α1X1,i+…+α10X10,i+i,(i=1,…,n)

(4)

復(fù)雜混雜結(jié)構(gòu)I：保持公式(4)不變，但處理因素改為由公式(5)生成。在這種情況下，X2和處理因素不僅存在線性關(guān)系，還存在二次項(xiàng)關(guān)系。

(5)

復(fù)雜混雜結(jié)構(gòu)II：保持公式(3)不變，但結(jié)局變量由公式(6)生成。在這種情況下，X1與結(jié)局變量之間的關(guān)系是非線性的。

(6)

復(fù)雜混雜結(jié)構(gòu)III：處理因素由公式(5)生成，結(jié)局變量由公式(6)生成。在這種情況下，處理因素模型與結(jié)局變量模型中均含有非線性關(guān)系。

按上述組合，分別生成了樣本量為500和2000的數(shù)據(jù)集。因?yàn)閷?shí)際情況中我們很難準(zhǔn)確判斷協(xié)變量因素與處理因素之間的關(guān)系，所以在估計(jì)傾向性評分值的過程中，我們將按照通常分析此類數(shù)據(jù)的做法，線性納入所有20個(gè)協(xié)變量，再利用傾向性評分匹配法(兩種卡鉗值)、回歸調(diào)整法、逆概率加權(quán)法、分層法分別估計(jì)處理效應(yīng)。各種情況下分別重復(fù)模擬1000次。

2.評價(jià)指標(biāo)

本研究評價(jià)指標(biāo)包括處理效應(yīng)的點(diǎn)估計(jì)(Average(αZ))及其標(biāo)準(zhǔn)誤(SE(αZ))、相對偏倚(RB)、均方誤差(MSE)。相對偏倚(RB)是點(diǎn)估計(jì)與真實(shí)效應(yīng)之差的絕對值占真實(shí)效應(yīng)的百分比，均方誤差為偏倚的平方與標(biāo)準(zhǔn)誤平方之和。

3.軟件實(shí)現(xiàn)

本研究采用R 3.2.3軟件進(jìn)行模擬實(shí)驗(yàn)。其中傾向性評分匹配法使用R軟件中的MatchIt包實(shí)現(xiàn)。下面這段程序可以分別實(shí)現(xiàn)本研究中的幾種傾向性評分方法：

#估計(jì)傾向性評分值

>psfit=glm(treat～x1+x2+x3+x4+x5+x6+x7+x8+x9+x10+x11+x12+x13+x14+x15+x16+x17+x18+x19+x20,data,family=binomial())

>ps=fitted(psfit)

#傾向性評分回歸調(diào)整法

>yfit=lm(outcome～treat+ps,data)

>cov.est=coef(summary(yfit))[2,1]

#傾向性評分逆概率加權(quán)法

>dataf$wgt=treat/ps+(1-treat)/(1-ps)

>iptw.fit=lm(outcome～treat,data,weights=wgt)

>iptw.est=coef(summary(iptw.fit))[2,1]

#傾向性評分分層法

>quintiles=quantile(data$ps,prob=seq(from=0,to=1,by=0.2),na.rm=T)

>data$pstrata=cut(data$ps,breaks=quintiles,labels=1:5,include.lowest=T)

>stratified=dlply(.data=data,.variables=“pstrata”,.fun=function(DF){lm(outcome～treat,data=DF)})

>sub.est=mean(sapply(stratified,function(mod)mod$coef[“treat”]))

#傾向性評分匹配法

>psmatch=matchit(treat～x1+x2+x3+x4+x5+x6+x7+x8+x9+x10+x11+x12+x13+x14+x15+x16+x17+x18+x19+x20,distance=“l(fā)ogit”,caliper=0.2,method=“nearest”,data)

>psmatch.data=match.data(psmatch)

>match.fit=lm(outcome～treat,data=psmatch.data)

>match.est=coef(summary(match.fit))[2,1]

4.模擬實(shí)驗(yàn)結(jié)果

從表1可見，當(dāng)樣本量為500且自變量之間相互獨(dú)立時(shí)，傾向性評分回歸調(diào)整法在四種不同的方法中估計(jì)效果最好，相對偏倚范圍在0.17%～0.50%，并且隨著模型復(fù)雜度的增加，即模型中含有二次項(xiàng)時(shí)，這種優(yōu)勢保持不變。匹配法設(shè)置了兩種卡鉗值，其中匹配法-1代表卡鉗值為兩組傾向性評分標(biāo)準(zhǔn)差的0.2倍，匹配法-2代表卡鉗值為0.02。由結(jié)果可見，當(dāng)卡鉗值設(shè)為0.02時(shí)，匹配法估計(jì)的較為準(zhǔn)確，相對偏倚在2.86%～7.55%，在混雜結(jié)構(gòu)最為復(fù)雜的情況下，即處理模型和結(jié)局模型均含有二次項(xiàng)時(shí)，相對偏倚達(dá)到最大(7.55%)。而匹配法-1估計(jì)的準(zhǔn)確性則較差，相對偏倚均在10%以上，并且隨著混雜結(jié)構(gòu)復(fù)雜性的增加，相對偏倚也逐漸增加。加權(quán)法在僅處理模型中含有非線性關(guān)系的情況下，估計(jì)情況較差，相對偏倚達(dá)到14.73%，在另外三種情況下偏倚相對較小，均保持在10%以內(nèi)。分層法在各種情況下估計(jì)的偏倚較大，其中在復(fù)雜混雜結(jié)構(gòu)II情況下偏倚達(dá)到了19.78%。在估計(jì)的穩(wěn)定性上，回歸調(diào)整法也較其他三種方法好，均方誤差最小，而加權(quán)法標(biāo)準(zhǔn)誤和均方誤差最大。

當(dāng)增加了X1與X2之間的相關(guān)性，結(jié)果見表2，當(dāng)處理因素模型和結(jié)局變量模型中均含有非線性關(guān)系時(shí)，回歸調(diào)整法的估計(jì)稍有偏倚，相對偏倚為3.58%，其余情況下，回歸調(diào)整方法仍估計(jì)準(zhǔn)確，相對偏倚范圍為0.31%～0.46%，穩(wěn)定性也較其他方法好，均方誤差最小。匹配法-2在簡單混雜結(jié)構(gòu)下，偏倚相對較小，而只要有一個(gè)模型中含有非線性關(guān)系時(shí)，相對偏倚增大，在復(fù)雜混雜結(jié)構(gòu)III情況下，相對偏倚達(dá)到11.39%。與表1結(jié)果相似，匹配法-1估計(jì)的偏倚依然較大。僅處理因素模型中含有非線性關(guān)系時(shí)，逆概率加權(quán)法估計(jì)的相對偏倚達(dá)到13.6%，在其他三種情況下，相對偏倚都較小，但標(biāo)準(zhǔn)誤和均方誤差在四種方法中為最大。分層法估計(jì)的準(zhǔn)確性仍舊較差，相對偏倚大部分都在10%以上，標(biāo)準(zhǔn)誤相對較小。繼續(xù)增加自變量之間的相關(guān)性，結(jié)果見表3，回歸調(diào)整法在準(zhǔn)確性和穩(wěn)定性方面，表現(xiàn)均最好。匹配法-1估計(jì)的相對偏倚隨著自變量之間相關(guān)性的增加而增大，匹配法-2在簡單混雜結(jié)構(gòu)情況下，估計(jì)也較前兩種情況準(zhǔn)確性差?？傮w來看，隨著協(xié)變量之間相關(guān)性的增加，這四種方法的估計(jì)的準(zhǔn)確性和穩(wěn)定性優(yōu)劣模式差異不大。

當(dāng)樣本量增加(n=2000)，回歸調(diào)整法估計(jì)較為準(zhǔn)確，匹配法-2表現(xiàn)較匹配法-1好，估計(jì)更為準(zhǔn)確。加權(quán)法在處理因素模型中不含有二次項(xiàng)的情況下估計(jì)效果較好，相對偏倚都在1%以下。當(dāng)處理因素與協(xié)變量之間含有非線性關(guān)系時(shí)，即在復(fù)雜混雜結(jié)構(gòu)I和III兩種情況下，加權(quán)法估計(jì)偏倚較大，而且后者相對偏倚較前者小。分層法在樣本量增加后估計(jì)結(jié)果仍不準(zhǔn)確。在穩(wěn)定性上，隨著樣本量增加，各種方法的標(biāo)準(zhǔn)誤和均方誤差均減小，回歸調(diào)整法的標(biāo)準(zhǔn)誤最小，均方誤差也最小，加權(quán)法估計(jì)的標(biāo)準(zhǔn)誤最大。

表1 不同模型情況下各種傾向性評分方法的估計(jì)結(jié)果(自變量之間獨(dú)立)

討 論

傾向性評分方法作為控制混雜的分析方法在醫(yī)學(xué)研究中應(yīng)用越來越廣泛。本研究的主要目的是在不同混雜結(jié)構(gòu)情況下，比較傾向性評分匹配法、回歸調(diào)整法、逆概率加權(quán)法以及分層法對處理效應(yīng)的估計(jì)優(yōu)劣。Austin在多篇文獻(xiàn)中也比較過這幾種傾向性評分方法[7-11]，但本研究內(nèi)容和他之前的文獻(xiàn)有所不同。在本研究中我們設(shè)置了不同的生成處理因素的模型和結(jié)局變量模型，通過加入?yún)f(xié)變量與處理因素和結(jié)局變量的非線性關(guān)系逐步增加混雜結(jié)構(gòu)的復(fù)雜性，根據(jù)文獻(xiàn)檢索，還沒有這方面深入的探討。另外，還模擬了協(xié)變量之間存在相關(guān)性的各種不同情況，更加貼近了實(shí)際。

表2 不同模型情況下各種傾向性評分方法的估計(jì)結(jié)果(R12=0.2)

表3 不同模型情況下各種傾向性評分方法的估計(jì)結(jié)果(R12=0.2，R13=0.3，R24=0.4)

傾向性評分回歸調(diào)整法在不同的情況下估計(jì)的相對偏倚較分層法和加權(quán)法都小，并且MSE最小，這與Austin得出的結(jié)論也相同[7]。由模擬結(jié)果可見，通過設(shè)置不同的卡鉗值，估計(jì)的結(jié)果相差較大。卡鉗值設(shè)為0.02比卡鉗值設(shè)為傾向性評分標(biāo)準(zhǔn)差的0.2倍準(zhǔn)確性好。傾向性評分匹配法是醫(yī)學(xué)文獻(xiàn)中應(yīng)用較多的調(diào)整混雜的方法，而卡鉗值的選擇對結(jié)果的估計(jì)有很大的影響。本研究采用以往文獻(xiàn)中建議的兩種卡鉗取值分別進(jìn)行匹配估計(jì)，得出的結(jié)果相差很大，因此在應(yīng)用卡鉗匹配方法的過程中，卡鉗值應(yīng)當(dāng)根據(jù)實(shí)際匹配后組間均衡性來考慮，對于不同的數(shù)據(jù)集情況下，以往研究建議的參數(shù)需要慎重選擇。

當(dāng)處理因素模型中不含有二次項(xiàng)時(shí)，即在簡單混雜結(jié)構(gòu)和復(fù)雜混雜結(jié)構(gòu)II這兩種情況下，傾向性評分加權(quán)法估計(jì)效果較好，并且隨著樣本量的增加，估計(jì)結(jié)果更為準(zhǔn)確，均方誤差也減小。Austin也在文獻(xiàn)中提到，應(yīng)用傾向性評分逆概率加權(quán)法可以得到風(fēng)險(xiǎn)差(risk difference，RD)的相合估計(jì)[9]。但在復(fù)雜混雜結(jié)構(gòu)I、III情況下加權(quán)法估計(jì)的準(zhǔn)確性又較差。在這兩種情況下，真實(shí)的處理因素模型含有非線性關(guān)系，而估計(jì)傾向性評分的過程中只納入了線性關(guān)系，因此計(jì)算的權(quán)重有誤，從而導(dǎo)致估計(jì)的處理效應(yīng)偏差較大。Linder等在文獻(xiàn)中得出類似結(jié)論[17]，即使結(jié)局模型指定錯(cuò)誤，基于PS的逆概率加權(quán)方法也可以得出相合估計(jì)，而當(dāng)PS模型指定錯(cuò)誤時(shí)，基于PS的逆概率加權(quán)方法估計(jì)結(jié)果較差。模擬結(jié)果中，分層法估計(jì)的相對偏倚在各種情況下均較大，而這也是分層法的一個(gè)缺陷，有文獻(xiàn)顯示其估計(jì)的偏倚較其他傾向性評分法更大[18]。

本研究模擬了協(xié)變量之間不同相關(guān)程度的情況，但在不同的相關(guān)性情況下，各種方法的估計(jì)優(yōu)劣情況類似。傾向性評分是將多個(gè)協(xié)變量的影響綜合成一個(gè)變量，起到一個(gè)“降維”的作用，這也是傾向性評分方法優(yōu)于傳統(tǒng)多元回歸方法的一個(gè)優(yōu)點(diǎn)，傳統(tǒng)多元回歸方法中可能需要考慮協(xié)變量之間的相關(guān)性、共線性問題，但我們的模擬結(jié)果顯示在不同協(xié)變量相關(guān)性情況下使用傾向性評分方法估計(jì)的結(jié)果似乎差異不大。

總體來說，從估計(jì)的相對偏倚大小和穩(wěn)定性方面考慮，我們認(rèn)為傾向性評分回歸調(diào)整法能夠更好地估計(jì)處理效應(yīng)，并且在我們設(shè)置的各種混雜情況下的估計(jì)都較為穩(wěn)健，因此當(dāng)協(xié)變量數(shù)目較多并且協(xié)變量與處理因素和結(jié)局變量的關(guān)系無法確定時(shí)，這四種方法中可以考慮優(yōu)先使用回歸調(diào)整法。本研究存在的局限性是未考慮結(jié)局變量是二分類、計(jì)數(shù)資料或生存資料的情況下，各種傾向性評分方法估計(jì)效果的優(yōu)劣。此外，值得注意的是，傾向性評分方法一個(gè)前提假定就是所有的混雜因素都已經(jīng)觀察到，對于未觀察到的重要混雜因素并不能進(jìn)行均衡，因此在后續(xù)研究中可以進(jìn)一步探索未包含重要混雜協(xié)變量時(shí)傾向性評分估計(jì)的準(zhǔn)確性。

[1]Johnston SC,Rootenberg JD,Katrak S,et al.Effect of a US National Institutes of Health programme of clinical trials on public health and costs.Lancet,2006,367(9519):1319-1327.

[2]Sturmer T,Joshi M,Glynn RJ,et al.A review of the application of propensity score methods yielded increasing use,advantages in specific settings,but not substantially different estimates compared with conventional multivariable methods.J Clin Epidemiol,2006,59(5):437-447.

[3]Rosenbaum PR,Rubin DB.The central role of the propensity score in observational studies for causal effects.Biometrika,1983,70(1):41-55.

[4]Hade E M,Lu B.Bias associated with using the estimated propensity score as a regression covariate.Stat Med,2014,33(1):74-87.

[5]Austin PC.The performance of different propensity score methods for estimating marginal odds ratios.Stat Med,2007,26(16):3078-3094.

[6]Austin PC,Grootendorst P,Anderson GM.A comparison of the ability of different propensity score models to balance measured variables between treated and untreated subjects:a Monte Carlo study.Stat Med,2007,26(4):734-753.

[7]Austin PC.The performance of different propensity-score methods for estimating relative risks.J Clin Epidemiol,2008,61(6):537-545.

[8]Austin PC,Grootendorst P,Normand SL,et al.Conditioning on the propensity score can result in biased estimation of common measures of treatment effect:a Monte Carlo study.Stat Med,2007,26(4):754-768.

[9]Austin PC.The performance of different propensity-score methods for estimating differences in proportions(risk differences or absolute risk reductions)in observational studies.Stat Med,2010,29(20):2137-2148.

[10]Pfeiffer RM,Riedl R.On the use and misuse of scalar scores of confounders in design and analysis of observational studies.Stat Med,2015,34(18):2618-2635.

[11]Austin PC.An Introduction to Propensity Score Methods for Reducing the Effects of Confounding in Observational Studies.Multivariate Behav Res,2011,46(3):399-424.

[12]Austin PC.Some methods of propensity-score matching had superior performance to others:results of an empirical investigation and Monte Carlo simulations.Biom J,2009,51(1):171-184.

[13]Hill J,Reiter JP.Interval estimation for treatment effects using propensity score matching.Stat Med,2006,25(13):2230-2256.

[14]Thoemmes FJ,Kim ES.A Systematic Review of Propensity Score Methods in the Social Sciences.Multivariate Behav Res,2011,46(1):90-118.

[15]Rosenbaum PR,Rubin DB.Reducing bias in observational studies using subclassification on the propensity score.Journal of the American Statistical Association,1984,79(387):516-524.

[16]Zou B,Zou F,Shuster JJ,et al.On variance estimate for covariate adjustment by propensity score analysis.Stat Med,2016,35(20):3537-3548.

[17]Linden A,Uysal SD,Ryan A,et al.Estimating causal effects for multivalued treatments:a comparison of approaches.Stat Med,2016,35(4):534-552.

[18]Austin PC.The performance of different propensity score methods for estimating marginal hazard ratios.Stat Med,2013,32(16):2837-2849.

(責(zé)任編輯：張 悅)

Comparison of Propensity Score Methods Under Different Confounding Structures:A Simulation Study

Sun Ting,Qin Guoyou,Wu Zhenyu,et al

(DepartmentofBiostatistics,SchoolofPublicHealthandKeyLaboratoryofPublicHealthSafety,MinistryofEducation,FudanUniversity(200032),Shanghai)

Objective The performance of propensity score(PS)methods were compared through constructing different confounding structures and generating covariates with different correlations when the outcome model was linear.Methods Monte Carlo method was used to simulate the datasets by constructing four confounding structures from simple to complex.Then four PS-based methods including PS matching,covariate adjustment,inverse probability of weighting(IPW)and stratification were applied to estimate the treatment effect.The results were compared from different aspects including the point estimate,standard error,relative bias and mean square error.Results When the outcome model was linear,covariate adjustment showed the least biased and stable estimates among the four methods.PS matching with caliper 0.02 performed better than the other matching methods when the caliper is 0.2 of the standard deviation of the PS value.When there were nonlinear relationship in the treatment model,IPW showed biased results and largest standard error.PS stratification resulted in biased estimates in all settings.Conclusion Covariate adjustment by PS is robust to complex confounding structure and achieved the least biased estimates.We propose that when the relationships between confounding factors and treatment or outcome variable cannot be confirmed,using PS covariate adjustment seems a better choice.

Propensity score;Confounding factor;Monte Carlo simulation;Bias

國家自然科學(xué)基金(編號：11371100)；上海市科研計(jì)劃項(xiàng)目(編號：13411950406)

△通信作者：武振宇，E-mail:zyw@fudan.edu.cn