胡良平

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

1 概 述

1.1 基本概念

本期《基于SAS與R軟件的主成分分析》一文介紹了“單組設(shè)計(jì)多元定量資料”這種特殊的數(shù)據(jù)結(jié)構(gòu)，并指出：分析這種數(shù)據(jù)結(jié)構(gòu)的統(tǒng)計(jì)分析方法占全部多元統(tǒng)計(jì)分析方法中的絕大部分，主成分分析法是其中一個(gè)最基本的方法，可以被靈活地運(yùn)用于多重回歸分析之中。

1.2 何為主成分回歸分析

主成分回歸分析(the Principal Components Regression Analysis)是將作為自變量的多個(gè)定量因素轉(zhuǎn)換為全部互相獨(dú)立的綜合變量(即主成分變量)，構(gòu)建定量因變量依賴所求得的“全部主成分變量”變化而變化的回歸模型，這樣完成的多重線性回歸分析被稱為“主成分回歸分析”[1]。

1.3 為何要使用主成分回歸分析

由于在經(jīng)典統(tǒng)計(jì)學(xué)中，要求自變量互相獨(dú)立，此時(shí)，才可以構(gòu)建多重線性回歸模型。當(dāng)自變量之間存在多重共線性時(shí)，經(jīng)典統(tǒng)計(jì)學(xué)理論認(rèn)為：所建立的多重線性回歸模型的質(zhì)量就不高，甚至可能是不能解決實(shí)際問題或違反專業(yè)知識的回歸模型(指某些回歸系數(shù)的正負(fù)號不符合基本常識和專業(yè)知識要求)。消除自變量之間的多重共線性關(guān)系，有多種具體方法，本文擬采用“主成分回歸分析法”。

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

【例1】為推算成年人的收縮壓(SBP,mmHg)，研究者測量并收集了50名成年人的年齡(age)、身高(height,inches)、體質(zhì)量(weight,pounds)和體質(zhì)量指數(shù)(BMI)。結(jié)果見表1[1]。試以SBP為因變量、其他四個(gè)變量為自變量建立多重線性回歸模型。

表1 50名成年人的測量數(shù)據(jù)

續(xù)表1：

33397321027．7113534776213825．2415035397323030．3412536406917025．1012637446211521．039938276114026．4511439297322029．0313940786311019．4915041626520834．6111242227112517．4312743376417630．2112544387219526．4513645226514023．3010846796112523．6215647246214626．7010848326714122．0810549427019227．5512150426818528．13126

1.5 對數(shù)據(jù)結(jié)構(gòu)的分析

在表1中，單從數(shù)據(jù)角度看，基本滿足進(jìn)行多重線性回歸分析的要求；但嚴(yán)格地說，題目中并沒有交代清楚：這50例受試者是從一個(gè)什么樣的總體中以什么方式抽取的，更具體的疑問是：他們在收縮壓(SBP)這個(gè)指標(biāo)上是否具有同質(zhì)性？若有些人屬于“嚴(yán)重高血壓患者”、有些人屬于“中度高血壓患者”、有些人屬于“輕度高血壓患者”、有些人屬于“臨界高血壓狀態(tài)人群”，還有些人屬于“正常血壓人群”，而且，這50人組成的樣本不是從各群體中按相同比例隨機(jī)抽取的，那么，他們就不具有同質(zhì)性，此資料也就不值得進(jìn)行任何統(tǒng)計(jì)分析。

在本文中，假定資料具有同質(zhì)性，且確定樣本含量也是有一定依據(jù)的。再結(jié)合醫(yī)學(xué)專業(yè)知識，人的年齡、身高、體質(zhì)量和體質(zhì)量指數(shù)確實(shí)對收縮壓有一定的影響(嚴(yán)格地說，應(yīng)找全所有可能影響收縮壓的因素，應(yīng)確定一個(gè)同質(zhì)的研究總體，按不少于80%的效能從此總體中隨機(jī)或分層隨機(jī)抽取足夠大的樣本含量)，筆者僅從統(tǒng)計(jì)建模角度展開下面的論述。

2 探索性分析

2.1 構(gòu)建基本的多重線性回歸模型

基于本刊2018年第1期《基于經(jīng)典統(tǒng)計(jì)思想實(shí)現(xiàn)多重線性回歸分析》一文介紹的方法，利用下面的SAS程序構(gòu)建基本的多重線性回歸模型：

data a1;

input id age height weight bmi sbp;

cards;

(此處輸入表1中50行6列數(shù)據(jù))

；

run;

proc reg data=a1;

model sbp=age height weight bmi;

run;

因篇幅所限，輸出結(jié)果從略。從輸出的結(jié)果中可看到：截距項(xiàng)無統(tǒng)計(jì)學(xué)意義，需要將其刪除后重新擬合多重線性回歸模型。在前面的過程步程序的“model語句”的“；”之前增加一個(gè)選擇項(xiàng)“/NOINT”，得到的結(jié)果表明：weight和BMI兩個(gè)自變量都無統(tǒng)計(jì)學(xué)意義，若將它們都從回歸模型中刪除，則可被利用的自變量就很少了。于是，考慮這4個(gè)自變量之間是否存在某種線性關(guān)系，從而導(dǎo)致不夠理想的結(jié)果出現(xiàn)。

2.2 對全部自變量進(jìn)行共線性診斷

采用下面的SAS過程步程序?qū)θ孔宰兞窟M(jìn)行共線性診斷：

proc reg data=a1;

model sbp=age height weight bmi/noint collin;

run;

【SAS主要輸出結(jié)果】

共線性診斷個(gè)數(shù)特征值條件指數(shù)偏差比例ageheightweightbmi13．827911．000000．008220．002410．000861790．0005833420．139685．235040．721180．010990．014740．0031430．0269511．918190．040000．933340．106010．0249840．0054726．457450．230600．053260．878380．97130

由上面的最后一行后4列的計(jì)算結(jié)果可知：weight和BMI的“偏差比例”都大于0.5，說明它們之間確實(shí)存在較嚴(yán)重的共線性關(guān)系。

3 主成分回歸分析

3.1 對全部自變量進(jìn)行主成分分析[2]

所需要的SAS過程步程序如下：

proc princomp data=a1 out=b1 prefix=z;

var age height weight bmi;

run;

此程序輸出結(jié)果為：產(chǎn)生4個(gè)主成分變量z1-z4及其取值，這是中間結(jié)果(暫不呈現(xiàn)了)存儲(chǔ)在名為b1的數(shù)據(jù)集中。

3.2 基于z1-z4作主成分回歸分析

所需要的SAS過程步程序如下：

proc reg data=b1;

model sbp=z1-z4;

run;

輸出結(jié)果(此處從略)中包含無統(tǒng)計(jì)學(xué)意義的變量，需要采取變量篩選方法淘汰掉無統(tǒng)計(jì)學(xué)意義的主成分變量。所需要的SAS過程步程序如下：

proc reg data=b1;

model sbp=z1-z4/selection=backward sls=0.05;

run;

輸出結(jié)果(此處從略)表明：此模型的復(fù)相關(guān)系數(shù)R平方=0.3262，殘差方差=200.41879。

3.3 僅對有共線性的自變量進(jìn)行主成分分析

所需要的SAS過程步程序如下：

proc princomp data=a1 out=b2 prefix=z;

var weight bmi;

run;

此程序輸出結(jié)果為：產(chǎn)生2個(gè)主成分變量z1-z2及其取值，這是中間結(jié)果(暫未呈現(xiàn))存儲(chǔ)在名為b2的數(shù)據(jù)集中。

3.4 基于獨(dú)立自變量age、height和主成分變量z1-z2作主成分回歸分析

所需要的SAS過程步程序如下：

proc reg data=b2;

model sbp=age height z1-z2/selection=backward sls=0.05;

run;

輸出結(jié)果(此處從略)表明：此模型的復(fù)相關(guān)系數(shù)R平方=0.3470，殘差方差=194.24348。

3.5 改變篩選變量的策略

前面篩選自變量時(shí)假定需要保留“截距項(xiàng)”，下面假定不需要保留“截距項(xiàng)”，所需要的SAS過程步程序如下：

proc reg data=b2;

model sbp=age height z1-z2/noint selection=backward sls=0.05;

run;

【SAS主要輸出結(jié)果】

方差分析源自由度平方和均方F值Pr>F模型37408772469591285．20<0．0001誤差479031．34023192．15618未校正合計(jì)50749908

變量參數(shù)估計(jì)值標(biāo)準(zhǔn)誤差I(lǐng)I型SSF值Pr>Fage0．522340．114394006．5347820．85<0．0001height1．516260．0695091470476．02<0．0001z13．209081．44250951．010674．950．0309

此模型的復(fù)相關(guān)系數(shù)R平方=0.9880，殘差方差=192.15618。

3.6 基于兩種途徑得到的主成分回歸模型的比較

前面獲得兩個(gè)二重回歸模型擬合效果接近，第2個(gè)比第1個(gè)稍好，因R平方更大、殘差方差更小，但它們之間的差別很微?。坏?個(gè)模型為三重線性回歸模型，殘差方差192.15618比第2個(gè)模型的194.24348略小，然而，R平方=0.9880比第2個(gè)模型的0.3470更大。

4 慎用主成分回歸分析

4.1 通過改變篩選自變量的策略來提升回歸模型的擬合效果

一般來說，不要急于采取主成分回歸分析方法，而應(yīng)首先考慮改變篩選自變量策略來提升模型擬合效果。就本例而言，直接采取前進(jìn)法、后退法和逐步法篩選全部自變量，并分別假定模型中包含“截距項(xiàng)”與不包含“截距項(xiàng)”。

經(jīng)嘗試，包含截距項(xiàng)時(shí)，三種篩選方法得出的結(jié)果一致，結(jié)果如下：模型對資料的擬合效果不太好，R平方=0.3723，殘差方差=186.73058。

經(jīng)嘗試，不包含截距項(xiàng)時(shí)，三種篩選方法得出的結(jié)果一致，結(jié)果如下：模型對資料的擬合效果較好，R平方=0.9883，殘差方差=186.89456。這個(gè)結(jié)果比上面基于主成分回歸分析得到的最好結(jié)果更好。

4.2 通過引入派生變量和改變篩選自變量的策略來提升回歸模型的擬合效果

可以引入4個(gè)自變量的平方項(xiàng)和交叉乘積項(xiàng)(共10項(xiàng))，再加上原先的4個(gè)自變量，總共14個(gè)自變量參與自變量的篩選。篩選時(shí)，仍采取前面提及的策略，即基于包含“截距項(xiàng)”與不包含“截距項(xiàng)”且分別采用前進(jìn)法、后退法和逐步法篩選，得到最好的結(jié)果如下：

方差分析源自由度平方和均方F值Pr>F模型874477093096761．06<0．0001誤差425137．61575122．32418未校正合計(jì)50749908

變量參數(shù)估計(jì)值標(biāo)準(zhǔn)誤差I(lǐng)I型SSF值Pr>Fage1．821820．492941670．8335513．660．0006weight-88．0080126．676361331．3848910．880．0020x3-0．009710．00342986．535298．060．0069x60．645690．193051368．3934811．190．0017x74．324561．309171334．7560410．910．0020x8-0．058350．018841173．866879．600．0035x90．785300．258361130．125969．240．0041x10-2．624580．877151095．188378．950．0046

以上結(jié)果表明：模型對資料的擬合效果很好，R平方=0.9931，殘差方差=122.32418。

前面的“x3、x6～x10”分別代表：x3=age*weight、x6=height*weight、x7=height*bmi、x8=weight*weight、x9=weight*bmi、x10=bmi*bmi。

這個(gè)模型對資料的擬合效果要好于前面未引入派生變量的最好模型的擬合效果。

值得注意的是：weight的回歸系數(shù)為“-88.00801”，這個(gè)“負(fù)值”表明：體重越重的人收縮壓(SBP)越低，這似乎不符合臨床專業(yè)知識。但應(yīng)當(dāng)注意：模型中還包含了“x6=height*weight”和“x9=weight*bmi”這兩項(xiàng)，它們的系數(shù)都為正，且乘積的結(jié)果是很大的數(shù)值。所以，在整體上，此模型是不違反臨床專業(yè)知識的。

4.3 小結(jié)

顯然，上面的最后一個(gè)多重線性回歸模型中的很多“項(xiàng)”之間存在嚴(yán)重的多重共線性關(guān)系，若對它們采取“主成分分析”提取主成分變量，再進(jìn)行主成分回歸分析，最終結(jié)果如下：

(1)基于全部8個(gè)變量(age,weight,x3,x6-x10)構(gòu)建主成分回歸模型且假定包含截距項(xiàng)，R平方=0.5931，殘差方差=129.30182。

(2)基于全部8個(gè)變量(age,weight,x3,x6-x10)構(gòu)建主成分回歸模型且假定不包含截距項(xiàng)， 無法構(gòu)建多重線性回歸模型。

(3)基于6個(gè)有共線性關(guān)系的變量(weight,x6-x10)產(chǎn)生的6個(gè)主成分變量再加上2個(gè)獨(dú)立變量(age和x3)構(gòu)建主成分回歸模型且假定包含截距項(xiàng)，R平方=0.6095，殘差方差=1560.03944。

(4)基于6個(gè)有共線性關(guān)系的變量(weight,x6-x10)產(chǎn)生的6個(gè)主成分變量再加上2個(gè)獨(dú)立變量(age和x3)構(gòu)建主成分回歸模型且假定不包含截距項(xiàng)，R平方=0.8892，殘差方差=1731.12908。

以上的結(jié)果都不如僅引入派生變量但不采取主成分變量構(gòu)建的多重線性回歸模型的擬合效果好。故筆者建議：一般應(yīng)慎用主成分回歸分析。

類似文獻(xiàn)[3]的資料都可能用到多重線性回歸分析，當(dāng)然，也就有可能會(huì)用到主成分回歸分析。經(jīng)驗(yàn)表明：應(yīng)慎用主成分回歸分析，而在引入派生變量的前提下，采取不同的篩選自變量的策略，有可能獲得比較理想的多重線性回歸模型。

參考文獻(xiàn)

[1] 胡良平. 面向問題的統(tǒng)計(jì)學(xué)——(2)多因素設(shè)計(jì)與線性模型分析[M]. 北京: 人民衛(wèi)生出版社, 2012: 215-228, 264-272.

[2] 胡良平. 面向問題的統(tǒng)計(jì)學(xué)——(3)試驗(yàn)設(shè)計(jì)與多元統(tǒng)計(jì)分析[M]. 北京: 人民衛(wèi)生出版社, 2012: 19-39.

[3] 趙巍峰, 彭敏, 謝博, 等. 健康教育對精神分裂癥患者病恥感影響的持續(xù)性[J]. 四川精神衛(wèi)生, 2017, 30(6): 519-523.