合理進(jìn)行多重線性回歸分析
——結(jié)合傾向性評分分析

胡純嚴(yán) ，胡良平 ，2*

（1.軍事科學(xué)院研究生院，北京 100850；2.世界中醫(yī)藥學(xué)會聯(lián)合會臨床科研統(tǒng)計學(xué)專業(yè)委員會，北京 100029*通信作者：胡良平，E-mail：lphu927@163.com）

在開展臨床研究和流行病學(xué)研究中，研究者收集到的資料通常都是多因素多指標(biāo)的資料。當(dāng)結(jié)果變量為計量變量時，研究者常選擇多重線性回歸模型來擬合資料［1-2］。當(dāng)資料中存在處理變量時，研究者最關(guān)注的是處理變量對結(jié)果變量所產(chǎn)生的因果效應(yīng)［3］。若資料為基于前瞻性隨機(jī)對照研究設(shè)計收集到的大樣本臨床資料，協(xié)變量在處理組與對照組之間的可比性較好，通常可直接進(jìn)行多重線性回歸分析。然而，若資料為基于回顧性研究設(shè)計或前瞻性調(diào)查研究設(shè)計收集到的資料，即使樣本含量非常大，也很難保證大多數(shù)協(xié)變量在處理組與對照組之間具有很好的可比性。直接基于可比性差的資料進(jìn)行任何統(tǒng)計分析，其結(jié)果的可信度都會大打折扣。本文將介紹如何結(jié)合傾向性評分分析［3-4］，合理進(jìn)行多重線性回歸分析。

1 基本概念

1.1 處理變量

在新藥的臨床試驗中，常規(guī)的做法是將某疾病患者隨機(jī)均分入試驗組與對照組，試驗組患者接受新藥治療，對照組患者接受對照藥物（通常為安慰劑）治療。這里，“藥物種類”就是研究者關(guān)注并考查的重要影響因素，它在統(tǒng)計學(xué)上常被稱為“處理變量”或“試驗因素”。在流行病學(xué)的調(diào)查研究中，也常涉及處理變量。例如，在研究吸煙與患肺癌關(guān)系的調(diào)查研究中，“是否吸煙”就是一個處理變量；又例如，在戒煙與體重改變量之間關(guān)系的調(diào)查研究中，“是否戒煙”就是一個處理變量。

1.2 處理變量兩水平組之間協(xié)變量的可比性

來自受試者的各種非試驗因素常被稱為協(xié)變量，例如年齡、性別、種族、血型、職業(yè)、生活方式、行為習(xí)慣、身高、體重、體重指數(shù)、血壓和血脂等。處理組與對照組中的受試者若在諸多協(xié)變量上的取值不同，對結(jié)果變量的影響也必然不同。因此，即使研究者觀察到兩組受試者在某些結(jié)果變量上存在明顯差別，也不完全是由處理水平與對照水平之間的效應(yīng)之差所致。

由統(tǒng)計學(xué)的理論可知，在樣本含量足夠大的隨機(jī)對照臨床試驗中，處理變量T的兩個水平組（T=1代表處理組、T=0代表對照組）的協(xié)變量之間通常是具有可比性的；而在觀察性研究中，處理變量兩個水平組之間協(xié)變量的可比性通常較差。一旦兩組的協(xié)變量之間出現(xiàn)嚴(yán)重的不均衡性，若直接基于這樣的資料進(jìn)行統(tǒng)計分析，其結(jié)果的可信度就會大打折扣。

1.3 傾向性評分分析法

消除非隨機(jī)對照試驗資料或觀察性研究資料中處理變量兩個水平組協(xié)變量之間不均衡性的一種有效方法為傾向性評分分析法。該方法的基本思想是使試驗組與對照組中的受試對象形成配對組，處在同一個配對組中的受試對象在所有協(xié)變量上的取值盡可能接近。形成配對組的依據(jù)是個體進(jìn)入相應(yīng)組的概率相等或差距在事先確定的一個很小的區(qū)間內(nèi)。

具體做法：將處理變量視為結(jié)果變量，將所有需要納入分析的協(xié)變量視為自變量，構(gòu)建關(guān)于“T=1”的二值結(jié)果變量的多重Logistic回歸模型，計算數(shù)據(jù)集中每個個體的發(fā)生概率Pi（i=1，2，…，n）。當(dāng)個體屬于處理組時，保持Pi不變；當(dāng)個體屬于對照組時，取其概率為1-Pi。有鑒于此，可理解Rosenbaum和Rubin對傾向性評分的定義，即根據(jù)一組觀察到的基線協(xié)變量分配處理的概率［3］。

2 匹配方法

2.1 匹配的必要性

在隨機(jī)對照試驗研究中，人們普遍接受這樣的假定：即任何一個個體都以同等的機(jī)會進(jìn)入試驗組或?qū)φ战M。由基本常識可知，當(dāng)樣本含量足夠大時，被隨機(jī)分配到處理組或?qū)φ战M中的個體在所有協(xié)變量上的取值是基本均衡的、可比的。這就意味著具有相同協(xié)變量向量的個體被分配進(jìn)入處理組或?qū)φ战M的概率是基本相同的，由此可知，將處理組與對照組中具有相同概率或概率值近似相等的個體配成對子，就意味著這些個體在協(xié)變量上的取值接近?；谶@一思路形成的匹配后的新數(shù)據(jù)集，就類似于前瞻性隨機(jī)對照試驗所得到的數(shù)據(jù)集。值得一提的是，匹配會損失掉一些信息，導(dǎo)致匹配后新數(shù)據(jù)集的樣本含量明顯小于原始數(shù)據(jù)集。

SAS/STAT的proc psmatch過程中介紹了三種匹配方法，即貪婪的最近鄰匹配、替換匹配和最優(yōu)匹配。當(dāng)用戶指定match語句時，proc psmatch過程使用下文中詳細(xì)描述的三種匹配方法之一，將對照組與處理組的觀察個體進(jìn)行匹配［3，5］。

2.2 三種匹配方法

2.2.1 貪婪的最近鄰匹配

貪婪的最近鄰匹配由method=greedy選項請求，選擇傾向性評分與每個處理單元（即處理組中的個體）的傾向性評分最匹配的控制單元（即對照組中的個體）。貪婪的最近鄰匹配是按順序進(jìn)行匹配，不進(jìn)行替換。以下條件可用于貪婪最近鄰匹配：①與每個處理單元匹配的控制單元數(shù)量（用戶可以在K=子選項中指定此數(shù)量）；②處理單元傾向性評分的順序，可以是升序、降序或隨機(jī)（用戶可以在order=子選項中指定順序）。

2.2.2 替換匹配

替換匹配由method=replace選項請求，通過替換選擇傾向性評分最接近每個處理單元傾向性評分的控制單元。用戶可以在K=子選項中指定要與每個處理單元匹配的控制單元數(shù)。

2.2.3 最優(yōu)匹配

最優(yōu)匹配同時選擇所有匹配項，不進(jìn)行替換，以最小化所有匹配項中傾向性評分的絕對總差異。用戶可以請求以下最優(yōu)匹配方法：①固定比率匹配，由method=optimal選項請求，將固定數(shù)量的控制單元與每個處理單元相匹配；②可變比率匹配，由method=varratio選項請求，將一個或多個控制單元與每個處理單元相匹配；③完全匹配，由method=full選項請求，將每個處理單元與一個或多個控制單元匹配，或?qū)⒚總€控制單元與一個或多個處理單元匹配，通過另外指定kmean=、ncontrol=或pctcontrol=subsoptions，用戶可以請求約束完全匹配，其中匹配的控制單元數(shù)小于可用的控制單元的總數(shù)。

通常，用戶按傾向性評分值進(jìn)行匹配，但用戶還可以依據(jù)基于傾向性評分算得的logit值匹配，或者依據(jù)馬氏距離來匹配。

3 實例與SAS實現(xiàn)

3.1 問題與數(shù)據(jù)結(jié)構(gòu)

【例1】研究者分析戒煙對個體體重的影響。本例數(shù)據(jù)是Hernán和Robins的NHANES I流行病學(xué)隨訪研究（NHEFS）數(shù)據(jù)的子集。在這項研究中，醫(yī)療和行為信息是在最初的身體檢查期間收集的，大約十年后再次在后續(xù)訪談中收集［3］。數(shù)據(jù)集SmokingWeight（樣本含量n=1 746）中部分?jǐn)?shù)據(jù)顯示在下面的數(shù)據(jù)步程序中，因篇幅所限，詳細(xì)數(shù)據(jù)見文獻(xiàn)［3］。

數(shù)據(jù)中的變量含義如下，Activity：日?；顒铀剑≈捣謩e為 0、1、2；Age：1971年時的年齡；BaseWeight：1971年時的體重（kg）；Change：隨訪和基線訪談中的體重差異（kg）；Education：教育水平，取值分別為0、1、2、3、4；Exercise：定期休閑運(yùn)動量，取值分別為0、1、2；PerDay：1971年每天吸煙數(shù)量；Quit：如果個體在初次訪談和后續(xù)訪談之間戒煙，取值為1，否則，取值為0；Race：白人取值為0，非白人取值為1；Sex：男性取值為0，女性取值為1；Weight：隨訪時的體重（kg）；YearsSmoke：吸煙的年數(shù)。

試構(gòu)建因變量Change關(guān)于全部自變量的多重線性回歸模型，著重考查處理變量Quit對因變量的影響程度和方向。

3.2 原始數(shù)據(jù)集中處理組與對照組協(xié)變量之間可比性比較

3.2.1 將10個協(xié)變量按性質(zhì)劃分成兩組

已知變量Quit為處理變量，Quit=1與Quit=0分別代表處理組與對照組；變量Change為結(jié)果變量，其他10個變量（Sex、Age、Race、Education、Exercise、BaseWeight、Weight、Activity、YearsSmoke、PerDay）均為協(xié)變量。由于Change=Weight-BaseWeight，故應(yīng)該將Weight和BaseWeight兩個變量從協(xié)變量的集合中刪除。為了顯示8個協(xié)變量在處理組與對照組中的可比性，需要將它們分成兩組：第一組為定性的協(xié)變量，包括Sex、Race、Education、Exercise、Activity；第二組為定量的協(xié)變量，包括Age、YearsSmoke、PerDay。

3.2.2 比較5個定性協(xié)變量在兩組之間的差別

設(shè)所需要的SAS程序如下：

【SAS主要輸出結(jié)果及解釋】兩組性別（χ2=8.810，P=0.003）、種族（χ2=6.820，P=0.009）、教育水平（χ2=12.678，P=0.013）構(gòu)成之間的差別均有統(tǒng)計學(xué)意義；兩組運(yùn)動量（χ2=4.490，P=0.106）、日常活動水平（χ2=2.681，P=0.262）構(gòu)成之間的差別均無統(tǒng)計學(xué)意義。由此可知：有3個定性協(xié)變量在兩組之間的差異具有統(tǒng)計學(xué)意義。

3.2.3 比較3個定量協(xié)變量在兩組之間的差別

設(shè)所需要的SAS程序如下：

【SAS主要輸出結(jié)果及解釋】兩組年齡（t=5.530，P＜0.01）、吸煙的年數(shù)（t=3.360，P=0.001）、1971年每天吸煙數(shù)量（t=3.990，P＜0.01）均值之間的差異均有統(tǒng)計學(xué)意義。由此可知：3個定量協(xié)變量在兩組之間的差異均有有統(tǒng)計學(xué)意義。

3.2.4 兩組協(xié)變量可比性的總結(jié)

由“第3.2.2節(jié)”和“第3.2.3節(jié)”的分析結(jié)果可知，本例資料中處理組（Quit=1）與對照組（Quit=0）之間有6個協(xié)變量不具有可比性。因此，基于這樣的資料構(gòu)建多重線性回歸模型，所得結(jié)果的可信度不高。

3.3 對原始數(shù)據(jù)集采用傳統(tǒng)的多重線性回歸分析

定量的結(jié)果變量為Change，著重考查的協(xié)變量Quit為處理變量，其他8個變量（Sex、Race、Education、Exercise、Activity、Age、YearsSmoke、PerDay）為一般協(xié)變量，構(gòu)建Change關(guān)于9個變量的多重線性回歸模型。設(shè)所需要的SAS程序如下：

【SAS主要輸出結(jié)果及解釋】由對整個回歸模型的假設(shè)檢驗結(jié)果可知，回歸模型有統(tǒng)計學(xué)意義（F=37.060，P＜0.01）；復(fù)相關(guān)系數(shù)平方的數(shù)值比較小，即R2=0.087?；貧w模型中各參數(shù)的假設(shè)檢驗結(jié)果見表1。

表1 回歸模型中各參數(shù)的假設(shè)檢驗結(jié)果Table 1 Hypothesis test results of the parameters in the regression model

以上是回歸模型中各參數(shù)的假設(shè)檢驗結(jié)果，截距項和4個協(xié)變量與0之間的差異均有統(tǒng)計學(xué)意義。值得一提的是：Quit的回歸系數(shù)為3.193＞0，表明戒煙使受試者的體重平均增加3.193 kg。

3.4 采用proc psmach過程進(jìn)行傾向性評分分析

3.4.1 調(diào)用傾向性評分過程

傾向性評分分析的目的是使所有的協(xié)變量在處理組與對照組之間的差別盡可能縮小，程序輸出結(jié)果主要是匹配后的數(shù)據(jù)集。設(shè)所需要的SAS程序如下：

3.4.2 基于傾向性評分分析的結(jié)果進(jìn)行后續(xù)分析

傾向性評分分析的結(jié)果儲存在輸出數(shù)據(jù)集OuttEx4中，主要內(nèi)容是在原始數(shù)據(jù)集中增加了一個新變量“MatchID”，即匹配變量。比較匹配后的數(shù)據(jù)集中協(xié)變量在處理組與對照組之間的差別。設(shè)所需要的SAS程序如下：

【SAS主要輸出結(jié)果及解釋】兩組性別（χ2=0.000，P=1.000）、種族（χ2=0.014，P=0.905）、教育水平（χ2=0.521，P=0.971）、運(yùn)動量（χ2=0.057，P=0.972）、日?；顒铀剑é?=0.327，P=0.849）構(gòu)成之間的差異均無統(tǒng)計學(xué)意義。兩組年齡（t=-0.150，P=0.879）、吸煙的年數(shù)（t=0.280，P=0.783）、1971年每天吸煙數(shù)量（t=0.230，P=0.821）均值之間的差異均無統(tǒng)計學(xué)意義。

由此可知：8個協(xié)變量在兩組之間的差異均無統(tǒng)計學(xué)意義。值得一提的是，本例原始數(shù)據(jù)集的8個協(xié)變量中有6個協(xié)變量在兩組之間的差異有統(tǒng)計學(xué)意義（見前文第3.2.4節(jié)）。這說明采用傾向性評分分析方法進(jìn)行匹配后獲得的數(shù)據(jù)集，所考查的協(xié)變量在處理組與對照組之間均達(dá)到了均衡性要求。

3.5 基于匹配后的數(shù)據(jù)集再擬合多重線性回歸模型

3.5.1 直接基于匹配后的數(shù)據(jù)集建模

設(shè)所需要的SAS程序如下：

【SAS主要輸出結(jié)果及解釋】由對整個回歸模型的假設(shè)檢驗結(jié)果可知，回歸模型有統(tǒng)計學(xué)意義（F=33.450，P＜0.01）；復(fù)相關(guān)系數(shù)平方的數(shù)值R2=0.110較之前的R2=0.087略有增大。回歸模型中各參數(shù)的假設(shè)檢驗結(jié)果見表2。

表2 回歸模型中各參數(shù)的假設(shè)檢驗結(jié)果Table 2 Hypothesis test results of the parameters in the regression model

以上是回歸模型中各參數(shù)的假設(shè)檢驗結(jié)果，截距項和3個協(xié)變量與0之間的差異均有統(tǒng)計學(xué)意義。值得一提的是：Quit的回歸系數(shù)為3.346＞0，表明戒煙使受試者的體重平均增加3.346 kg。

3.5.2 在匹配后的數(shù)據(jù)集中引入派生自變量后再建模

引入3個派生自變量（age1、age2、age3）后重新構(gòu)建多重線性回歸模型［6-7］，設(shè)所需要的SAS程序如下：

【SAS主要輸出結(jié)果及解釋】由對整個回歸模型的假設(shè)檢驗結(jié)果可知，回歸模型有統(tǒng)計學(xué)意義（F=55.260，P＜0.01）；復(fù)相關(guān)系數(shù)平方的數(shù)值R2=0.214較之前的R2=0.110有較大增加?；貧w模型中各參數(shù)的假設(shè)檢驗結(jié)果見表3。

表3 回歸模型中各參數(shù)的假設(shè)檢驗結(jié)果Table 3 Hypothesis test results of the parameters in the regression model

以上是回歸模型中各參數(shù)的假設(shè)檢驗結(jié)果，無截距項，4個協(xié)變量與0之間的差異均有統(tǒng)計學(xué)意義。值得一提的是：Quit的回歸系數(shù)為3.269＞0，表明戒煙使受試者的體重平均增加了3.269 kg。

3.6 結(jié)論

由“第3.2節(jié)”的比較結(jié)果可知，原始數(shù)據(jù)集中處理組與對照組協(xié)變量之間的可比性較差。因此，“第3.3節(jié)”基于原始數(shù)據(jù)集構(gòu)建多重線性回歸模型所得結(jié)果和結(jié)論的可信度不高。由“第3.4節(jié)”的比較結(jié)果可知，匹配后的數(shù)據(jù)集中處理組與對照組協(xié)變量之間的可比性很好。因此，“第3.5節(jié)”基于匹配后的數(shù)據(jù)集構(gòu)建多重線性回歸模型所得結(jié)果和結(jié)論的可信度有了較大提高（模型的復(fù)相關(guān)系數(shù)的平方由R2=0.087提高到R2=0.214）。

4 討論與小結(jié)

4.1 討論

傾向性評分分析法是對非隨機(jī)試驗資料或調(diào)查研究資料進(jìn)行預(yù)處理的有效方法，它特別適用于具有一個二水平處理變量的多因素資料?；谄ヅ浜蟮臄?shù)據(jù)集，可以采用多重線性回歸分析（要求結(jié)果變量為計量變量）或多重Poisson回歸分析（要求結(jié)果變量為計數(shù)變量）或多重Logistic回歸分析（要求結(jié)果變量為定性變量）或多重Cox比例風(fēng)險回歸分析（要求結(jié)果變量為生存時間變量）。

匹配后的數(shù)據(jù)集樣本含量的大小取決于所采用的匹配方法，一般來說，匹配方法越嚴(yán)格，匹配后的數(shù)據(jù)集樣本含量就越小。就本例而言，原始數(shù)據(jù)集的樣本含量n=1 746，采用貪婪的最近鄰匹配法得到的匹配后數(shù)據(jù)集樣本含量n=854，采用替換匹配法得到的匹配后數(shù)據(jù)集樣本含量n=741，采用最優(yōu)匹配法得到的匹配后數(shù)據(jù)集樣本含量n=1 498。值得一提的是，基于不同的匹配后數(shù)據(jù)集進(jìn)行多重線性回歸分析，其結(jié)果可能會不同。

4.2 小結(jié)

本文介紹了與傾向性評分分析有關(guān)的基本概念、三種匹配方法以及基于一個實例的SAS實現(xiàn)?；靖拍畹膬?nèi)容涉及處理變量、處理變量兩水平組之間協(xié)變量的可比性和傾向性評分分析法；三種匹配方法分別是貪婪的最近鄰匹配、替換匹配和最優(yōu)匹配；運(yùn)用SAS軟件對一個實例進(jìn)行了全面分析，包括對原始數(shù)據(jù)集基本情況的描述、基于傾向性評分分析產(chǎn)生匹配后的數(shù)據(jù)集、分別對原始數(shù)據(jù)集和匹配后的數(shù)據(jù)集進(jìn)行多重線性回歸分析，并對基于不同數(shù)據(jù)集所得到的分析結(jié)果進(jìn)行比較和總結(jié)。