陳文松 劉 曼 劉玉秀,4△ 許敏怡 熊 殷 鞏浩雯 李維勤
【提 要】 目的 探討觀察性研究中用于混雜偏倚控制的傾向性評分匹配、馬氏距離匹配和遺傳匹配三種方法的性能。方法 針對連續型結局變量,設定混雜變量與處理分組變量之間具有不同復雜度的回歸模型結構,采用Monte-Carlo模擬方法比較三種匹配方法在處理組間效應估計和匹配前后自變量均衡的區別,進而對三種方法性能進行評估。結果 在給定的模擬情形下,相比于傾向性評分匹配和馬氏距離匹配,遺傳匹配法得出的效應估計偏差最小,匹配后兩組自變量均衡性最好。結論 遺傳匹配在三種匹配方法中表現出較好的統計性能,可考慮作為觀察性研究中控制混雜偏倚優先推薦的匹配方法。
在非隨機觀察性研究中,當需要考慮某目標處理變量的效應時,由于混雜偏倚的客觀存在,常常導致估計出的組間差別不能反映真實情況,此時混雜偏倚的控制就尤為關鍵。除已知的分層分析、多因素分析等傳統方法,近年來,利用距離度量(metric)對研究個體進行匹配的傾向性評分匹配法(propensity score matching,PSM)以及馬氏距離匹配法(Mahalanobis distance matching,MDM)、遺傳匹配法(genetic matching,GM)在醫學研究中的應用越來越廣泛[1-3],然而有關這幾種匹配方法在統計推斷時的性能比較研究尚未見詳細文獻報道。
本文將在介紹基于進化算法的遺傳匹配方法的基礎上,采用Monte-Carlo模擬方法對三種匹配方法的統計性能進行比較研究,并通過一個實例進行分析說明。
理論上,將同質性相近的個體在不同的組間進行匹配,組間的自變量的分布將趨于均衡,從而減少或抵消兩組混雜因素對治療效果的影響,這就是匹配的基本原理[4-5]。在醫學研究中,除規定匹配個體所有變量取值均相同的精確匹配外,以馬氏距離為基礎的多元變量匹配和將多維數據轉化為一維的傾向性評分匹配最為常見。有關傾向評分匹配已有較多文獻介紹,此不贅述。以下重點介紹基于馬氏距離發展的廣義馬氏距離,以及以廣義馬氏距離為度量的遺傳匹配。
1936年,P.C.Mahalanobis提出了馬氏距離(Mahalanobis distance,MD),用以表示數據的協方差距離[6]:
(1)
其中Xi和Xj分別為第i個個體和第j個個體的自變量組成的向量,S-1為樣本協方差矩陣的逆。個體間馬氏距離越小,表示個體間越相似。
2013年,Diamond及其同事提出了基于進化算法的遺傳匹配(genetic matching)來控制組間混雜因素的影響。遺傳匹配在馬氏距離的基礎上創立了更一般化的距離度量,即廣義馬氏距離(generalized Mahalanobis distance,GMD)[7]:
(2)

個體的自變量集合X除自變量以外,還可以包含傾向性評分。此時方程(1)中的X可以用Z代替,其中Z是一個由傾向性評分π(X)和自變量X組成的矩陣。如果將PS以外的其他變量的權重設為0,此時GMD等同于傾向性評分;如果將PS的權重設為0,其他自變量的各個權重設為1,此時GMD等價于馬氏距離。PSM和MDM都是遺傳匹配的極限情況。
常見的計算傾向性評分的方法有logistic回歸、Probit回歸,以及機器學習估計中的神經網絡、支持向量機、分類與回歸樹、隨機森林等方法[8-9]。其中利用logistic回歸計算傾向性評分最為常見,即利用組別為因變量,建立logistic模型,本文就采用這種方法計算傾向性評分。
遺傳匹配根據每個變量的相對重要性對其加權來獲取最佳的自變量平衡。遺傳匹配利用遺傳算法(genetic algorithms)搜索一系列距離度量,通過最小化損失函數,找到使得自變量均衡性最高的加權矩陣W。該算法允許使用者自定義損失函數,默認損失函數定義為匹配后兩組自變量的Kolmogorov-Smirnov檢驗的P值及配對t檢驗的P值最小值,用以描述兩組自變量的最大差異。
遺傳算法是一種通過模擬達爾文生物進化論的自然選擇和遺傳學機理的生物進化過程搜索最優解的方法[10-11]。圖1總結了遺傳算法的流程。該算法產生一批初代權重Ws,每個W對應一個不同的距離度量。在每一代中,計算每個匹配樣本對應的損失函數。選擇使得損失函數最小的權重作為父代,產生子代權重并不斷迭代這一流程,直至漸進收斂于最優解。每一代產生的權重數量稱為種群規模,在迭代過程中種群規模保持不變。增大種群規模通常會改善遺傳匹配的效果,但會增加算法尋找最優解所花費的時間。

圖1 遺傳算法流程圖
當損失函數滿足迭代標準時停止迭代,此時的權重即為最優權重。利用最優權重計算廣義馬氏距離,根據此距離度量進行匹配,得到匹配后樣本。與傾向性評分一樣,此時可以使用不同的匹配策略,例如有無替換、有無卡鉗的最近鄰匹配或最優匹配等。
值得注意的是,遺傳匹配與其他匹配方法的使用前提一樣,即要事先確定匹配的具體變量、選定測量自變量平衡性的指標,以及指定合適的匹配方式。
本模擬實驗基于模擬生成的數據集,采用Monte-Carlo方法比較遺傳匹配與經典的馬氏距離匹配、傾向評分匹配在傾向性評分模型不同復雜程度下的效應估計精度和匹配前后自變量均衡程度,評價其統計性能,為方法學選擇提供理論依據。
Setoguchi及其同事在探究機器學習技術構建傾向性評分模型的研究中提出了一種模擬數據生成策略[12],該策略產生的模擬數據包含了實際中可能出現的多數情景,具有較好的代表性。本研究在此策略的基礎上進行了略微的修改。利用SAS軟件總共生成1000個模擬數據集,其中每個數據集的樣本量為1000。在每個數據集中,模擬生成10個自變量,其中4個混雜變量,3個暴露預測變量和3個結果預測變量。4個自變量(X2,X4,X7,X10)服從標準正態分布,其余6個二元變量(X1,X3,X5,X6,X8,X9)則由標準正態分布轉換得到,具體的模擬數據變量結構及自變量間相關系數見圖2。設定group為組別變量,并使其在自變量均值上約等于0.5以保證兩組人數相同。Y作為連續型結局變量,由group和其他自變量的線性組合生成。在文后附錄中提供了計算傾向性評分、結局生成及治療分配的詳細過程。

圖2 模擬數據各變量關系的基本結構
基于所設定的模擬數據結構,本研究指定了7種復雜程度各不相同的情形來指代真正的傾向性模型,以便評估三種匹配方法在傾向性模型不同復雜程度下的表現情況:
A:線性可加(僅主效應)
B:輕微非線性(包含1個二次項)
C:中等非線性(包含3個二次項)
D:輕微非可加性(包含3個兩因素交互作用項)
E:輕微非可加性和非線性(包含3個兩因素交互作用項和1個二次項)
F:中等非可加性(包含10個兩因素交互作用項)
G:中等非可加性和非線性(包含10個兩因素交互作用項和3個二次項)
為評價模擬生成的數據集的統計特征,本研究比較了兩組間自變量的均衡性。圖3展示了傾向性評分模型為情形A時兩組間傾向性評分分布與自變量均衡情況,其他情形下的模擬數據表現與其類似,這里不予羅列。從圖3可知,試驗組中傾向性評分大于0.5,對照組中傾向性評分小于0.5,兩組傾向性評分分布互有重疊;兩組間混雜變量與暴露預測變量的標準化差值大于0.1,混雜變量與暴露預測變量在兩組間分布不均衡。

圖3 情形A時模擬數據集兩組自變量均衡性
本研究使用MDM、PSM和GM來估計每個情形中的group的效應。匹配方案均選擇有放回的最近鄰匹配(nearest neighbor matching),匹配比例為1:1,同時PSM的卡鉗值設為0.02,GM中每代的種群規模(pop.size)設置為500。由于實際數據的復雜性,本研究很難確定哪些變量應納入模型以及模型中變量間的具體關系,所以在PSM和GM計算傾向性評分時,采用兩種策略建立傾向性評分模型:在第一種策略中,默認所有自變量均以線性可加的形式包含在傾向性評分模型中;而在第二種策略中,仍舊將所有自變量納入傾向性評分模型,但對于連續型變量,使用限制性立方樣條(restricted cubic spline,RCS)來考慮連續型變量的非線性關系。MDM使用所有自變量計算馬氏距離。三種匹配方法均應用R軟件中的Matching包實現。此外,本研究還在遺傳匹配中去除傾向性評分,只使用所有自變量,探究其對結果的影響。其中,使用“PSM”表示使用線性可加模型進行的傾向性評分匹配,“PSM+”表示使用非線性可加模型進行的傾向性評分匹配,“MDM”表示馬氏距離匹配,“GM”表示使用線性可加模型得到的ps進行的遺傳匹配,“GM-”表示不使用ps進行的遺傳匹配,“GM+”表示使用非線性可加模型得到的ps進行的遺傳匹配。
考慮采用2個度量值作為效應估計性能的評價指標,分別為:①相對偏倚(average percent relative bias,APRB):估計值與真實值的相對差異百分比;其中采用配對t檢驗對匹配后樣本計算估計值。②均方根誤差(root mean squared error,RMSE):估計值與真實值偏差的平方與估計次數比值的平方根。另外,采用2個指標來評價匹配前后自變量均衡性,分別為:①標準化差值(standardized mean difference,SMD):用于定量評價均衡性。一般認為標準化差值小于0.1時具有良好的均衡性。②配對t檢驗或配對卡方檢驗的P值。對于連續性變量使用配對t檢驗,對于分類變量使用配對卡方檢驗。考慮到篇幅所限,本文中自變量均衡性結果僅提供標準化差值。
在不同的傾向性評分模型下,使用不同方法對模擬數據集分析,比較不同方法的效應估計性能評價指標和匹配前后自變量均衡性。
(1)不同匹配方法后的處理效應估計的比較
從表1可見,各方法在各情形下的總體偏倚從小到大依次為:GM≈GM+ 表1 不同匹配方法效應估計相對偏倚的比較 (2)不同匹配方法后組間自變量的均衡性比較 考慮結果的呈現性,本研究在此挑選了7種情形中的具有代表性的4種情形,并選擇了各類匹配方法中表現最優的3種方法。從圖4可見,經遺傳匹配處理后,兩組自變量的標準化差值最小,在各種情形中的自變量均衡情況表現均最優。在指定了正確的傾向性評分模型后(即情形A),PSM+與GM在整體上表現相近,但其SMD的分布略高于GM。隨著真實傾向性評分模型復雜度的增加,PSM+的兩組自變量均衡程度開始下降。MDM在控制分類變量的均衡性上表現較優,但在連續型變量上表現最差,其中在匹配前最不均衡的變量X7上控制表現最差。 圖4 三種匹配方法匹配各自變量均衡性情況:SMD (3)傾向性評分匹配中不同建模策略對統計性能的影響 圖5給出了在不同的真實傾向性評分模型下,使用兩種不同建模策略進行傾向性評分匹配的自變量均衡情況。從圖5可見,使用非線性可加模型進行傾向性評分匹配的組間自變量均衡程度更優,在任何情形下各變量SMD的上四分位數均小于0.1。當真實傾向性模型僅為簡單的線性可加模型時(情形A),PSM與PSM+的表現一致;當真實傾向性模型包含非線性項時(情形C、情形E、情形G),在二元變量上,PSM與PSM+的表現一致。但在連續型變量上,PSM+的表現優于PSM。從表1可見,在任何情形下,PSM+的療效相對偏倚小于PSM。 (4)遺傳匹配中使用傾向性評分與否對統計性能的影響 圖6給出了在不同的真實傾向性評分模型下,使用簡單線性相加模型得到的PS、使用非線性相加模型得到的PS與不使用PS進行遺傳匹配的自變量均衡情況,可見,不使用PS作為自變量進行遺傳匹配的組間自變量均衡程度更差。當計算PS的模型與真實傾向性評分模型更相近時,使用PS的遺傳匹配效果更優。GM與GM+的表現大體相似,在復雜情形下GM+的表現略優于GM。從表1可見,使用PS進行遺傳匹配的偏倚更小,從小到大依次為GM≈GM+ 圖5 不同傾向性評分建模策略匹配各自變量均衡情況:SMD 圖6 不同策略遺傳匹配后各自變量均衡情況:SMD Dehejia等人創建的數據集“lalonde”基于一個全國性的職業培訓研究,探討接受職業培訓是否會增加學員收入,該數據集內置于“Match”包中。該研究的研究對象共445例(試驗組185例,對照組260例),基線變量包含年齡、受教育時間、人種等8個變量。結局變量為研究對象在1978年的收入。本研究使用MDM、PSM與GM對此進行分析,匹配方案與先前模擬試驗中一致。匹配前后自變量均衡情況見表2及圖7。匹配前,有4個自變量在兩組間的SMD大于0.1(有無畢業文憑、族裔、受教育時間、年齡),其中“有無畢業文憑”在兩組間的P值小于0.001。使用三種方法匹配后自變量均衡情況均有所改善,但GM匹配后兩組自變量均衡性最高,匹配后兩組SMD最小。 表2 匹配前后自變量特征 圖7 匹配前后兩組標準化差異 本文在簡要介紹利用進化算法自動迭代地構建匹配度量的遺傳匹配的基礎上,借助Monte-Carlo模擬方法,通過效應估計的準確性及匹配后組間自變量的均衡性兩個方面評價了馬氏距離匹配、傾向性評分匹配和遺傳匹配在數據不同復雜度下的統計性能。模擬結果顯示,當指定的傾向性評分模型與真實的傾向性評分模型相近時,傾向性評分匹配表現出了較好的統計學性能。但當模型復雜度上升后,即指定的傾向性評分模型與真實模型相差較大時,傾向性評分匹配的偏倚開始增大。而遺傳匹配不會受真實模型的影響,在不同模型復雜度下得出的效應估計偏差最小,匹配后兩組自變量均衡性最高。實例分析結果也驗證了,遺傳匹配在保持較多匹配例數下仍能達到最好的組間均衡。 近年來,傾向性評分匹配在醫學研究中得到了廣泛的應用,但在實際使用中仍存在著如何確定合理傾向性評分模型的問題。Granger等人統計了在2014年到2016年發表在高影響因子的醫學期刊上使用傾向性評分進行分析的研究,其中20.9%(187/894)的研究并未闡述如何評估所建立的傾向性評分模型的合理性[13]。通常情況下,當研究人員根據數據特征并結合自身專業知識盡可能地納入混雜因素,確定傾向性評分模型后,可以通過評價匹配后自變量平衡情況來判斷是否指定正確并進行調整。Rosenbaum和Rubin曾建議應根據匹配后自變量均衡情況迭代修改傾向性評分模型已達到最優均衡,但在實際情況中許多研究者缺乏這種建立合理的傾向性評分模型的意識或技能,這無疑會導致傾向性評分匹配在醫學研究中的錯用,得到并不真實的研究結論。 遺傳匹配將損失函數定義為匹配后組間均衡性指標,利用遺傳算法自動地找到最優的匹配組合,從而實現自變量的均衡。與傾向性評分匹配相比,遺傳匹配不需要研究者嚴格指定傾向性評分模型,不需要研究者根據匹配后自變量平衡情況手動調整傾向性評分模型,這極大地降低了研究者的工作難度。目前尚無評價匹配后自變量均衡性的最優指標,除最常見的組間檢驗P值及標準化差值外,還有方差比、Hosmer-Lemeshow檢驗統計量、KS檢驗統計量及C統計量等[13]。在本文使用的R軟件包“Matching”中則可自定義損失函數,使用任何均衡性指標進行遺傳匹配。遺傳匹配獲得最優的匹配度量后,匹配的方案與其他類似。傳統的最近鄰匹配、卡鉗匹配、1:n匹配等在遺傳匹配中仍可使用。 當然,遺傳匹配與其他匹配方法一樣,只能處理已知的混雜因素,對未知混雜無法控制,也不適用于時依自變量所致混雜因素的處理[8-9];當其使用傾向性評分作為自變量進行匹配時,也面臨傾向性評分模型自變量的選取問題。 本研究評價了使用不同傾向性評分建模策略對傾向性評分匹配和遺傳匹配的影響程度,研究結果表明,在建立傾向性評分模型時,考慮自變量的非線性關系會提高匹配的效果。另外,與不納入PS的遺傳匹配相比,納入PS進行遺傳匹配具有更好的統計性能,且當使用的傾向性評分模型越正確,遺傳匹配的效果越好,這恰好支持了Diamond等人的建議。這一結論也符合雙重穩健估計的理論[14],即納入PS的遺傳匹配提升了效應估計的穩健程度。 在真實世界研究中,為去除混雜因素的影響,得到無偏的效應估計,眾多學者提出了許多方法,例如傳統的分層分析、多因素回歸分析、傾向性得分匹配、逆概率加權匹配、重疊加權等方法。近年來一些更為復雜的方法也被提出,例如G估計、雙重穩健估計等方法。本文因篇幅所限,僅比較了三種匹配方法的統計表現,遺傳匹配與其他方法的比較有待進一步研究。 本研究只比較了在設定的模擬情形下三種匹配方法的結果,此假設下模擬數據盡管能滿足組間傾向性評分的重疊假定(overlap assumption),但并未對方法的匹配度進行探討比較,也未探討兩組傾向性評分的重疊程度對匹配效能的影響;另外,本研究的模擬試驗只考慮了1∶1有放回的匹配策略,未進行其他比例的匹配,這些問題均有待進一步研究。



案例分析


討 論