胡良平
(1.軍事醫(yī)學(xué)科學(xué)院研究生院,北京 100850;2.世界中醫(yī)藥學(xué)會聯(lián)合會臨床科研統(tǒng)計學(xué)專業(yè)委員會,北京 100029
在進行回歸分析時,對“多值名義(也包括多值有序)自變量”進行“啞變量變換”已成為統(tǒng)計學(xué)界公認的“標準做法”。若采用某種篩選自變量的方法,就可能保留一部分啞變量而淘汰掉另一部分啞變量,這樣做的前提是“假定它們互相獨立”,但事實上,源自一個多值名義自變量的多個啞變量之間是有聯(lián)系的,即它們之間并不滿足“互相獨立”的要求;按理說,由一個多值名義自變量產(chǎn)生的全部啞變量,要么全部被保留在回歸模型中,要么全部被剔除到回歸模型之外。那么問題又出現(xiàn)了:當它們?nèi)砍霈F(xiàn)在回歸模型中時,有些啞變量可能就沒有統(tǒng)計學(xué)意義,這樣的回歸模型就不是最節(jié)儉的,因此,回歸模型的擬合優(yōu)度不高;若全部剔除這些啞變量,就相當于沒有發(fā)揮該多值名義自變量的作用??傊?,“啞變量變換”存在一些弊端,是一個值得深究的“統(tǒng)計疑難問題”。
1.2.1 何為“算術(shù)均值變換”
1.2.2用“算術(shù)均值變換”取代“啞變量變換”的合理性
將一個“多值名義自變量”變換成一個“定量自變量”,這似乎很不合理。若孤立地考察“多值名義自變量”各水平之間的關(guān)系,前述的變換的確很不合理。因為“多值名義自變量”的各水平之間只是“名稱或符號上的不同(例如ABO血型系統(tǒng)中的A型、B型、AB型和O型血型之間的關(guān)系)”,而沒有“數(shù)量或程度上的差別”。然而,既然叫做“自變量”,它一定要影響其他變量,其中,定量結(jié)果變量(也叫定量因變量)是進行定量資料差異性分析和回歸分析時最關(guān)注的“變量”。換言之,無論是進行t檢驗或方差分析,還是進行回歸分析或判別分析,永遠不可能孤立地去考察“自變量”各水平之間的關(guān)系,而不可避免地要建立起“自變量各水平”與因變量的具體取值之間的數(shù)量聯(lián)系。
在對定量資料進行差異性分析(如t檢驗或方差分析)時,統(tǒng)計學(xué)上也是這樣做的。例如:設(shè)“藥物種類(因素A)”有兩個水平:試驗藥(A1水平)和對照藥(A2水平)。兩個藥物組都有較多數(shù)量的受試對象參與試驗,假定評價藥物療效的指標為“收縮壓下降值”,測定并計算出兩組“收縮壓下降值”的“算術(shù)平均值”。顯然,這就是采用了各自的“算術(shù)平均值”代替“試驗藥(A1水平)”與“對照藥(A2水平)”的“療效”水平。顯而易見,“A1”與“A2”之間是沒有數(shù)量大小之分的,但它們各自組中療效的“算術(shù)平均值”是有數(shù)量大小之分的。在多因素定量資料的方差分析中,情況也是如此。
由此可知,“算術(shù)均值變換”是用“動態(tài)思維(利用了因果關(guān)系)”建立起自變量各水平與定量結(jié)果變量之間數(shù)量關(guān)系的一種處理方法;而“啞變量變換”是用“靜態(tài)思維(孤立地看待自變量)”處置自變量各水平之間關(guān)系的方法,此法割裂了“多值名義自變量”與“定量因變量”之間的數(shù)量聯(lián)系。因此,“算術(shù)均值變換”比“啞變量變換”更具有合理性。
沿用本期科研方法專題第一篇文章《提高回歸模型擬合優(yōu)度的策略(Ⅰ)——啞變量變換與其他變量變換》中的“實際問題與數(shù)據(jù)結(jié)構(gòu)”[1],此處不再贅述。
2.1.1對“燃油種類(fuel)”這個“6值名義自變量”進行“算術(shù)均值變換”
求出“燃油種類(fuel)”各水平下“氧化氮釋放量(nox)”的算術(shù)均值所需要的SAS程序如下:
/*下面的SAS程序計算出多值名義變量各水平下定量因變量的算術(shù)均值*/
data aa;
set sashelp.gas;
proc sort data=aa;
by fuel;
run;
proc univariate data=aa noprint;
var nox;
by fuel;
output out=aaa means=m_nox;
run;
proc print data=aaa;
run;
【SAS輸出結(jié)果】
ObsFuelm_nox182rongas3.52589294%Eth2.087883Ethanol1.957384Gasohol3.349385Indolene3.546596Methanol1.55967
以上就是6種燃油對應(yīng)的“氧化氮釋放量(nox)”的算術(shù)均值。
將“燃油種類(fuel)”的各水平變換成與定量因變量相應(yīng)的“算術(shù)均值”所需要的SAS程序如下:
/*下面的SAS程序?qū)⒍嘀得x變量各水平變換成與定量因變量相應(yīng)的算術(shù)均值*/
data a1;
set sashelp.gas;
if fuel=' 82rongas' then mfuel=3.52589;
else if fuel=' 94%Eth' then mfuel=2.08788;
else if fuel=' Gasohol' then mfuel=3.34938;
else if fuel=' Indolene' then mfuel=3.54659;
else if fuel=' Methanol' then mfuel=1.55967;
else if fuel=' Ethanol' then mfuel=1.95738;
run;
通過運行上面的SAS程序,在原數(shù)據(jù)集sashelp.gas中就增加了一個定量自變量mfuel,將用它取代多值名義自變量“燃油種類(fuel)”。
2.1.2 對定量因變量和自變量不進行任何變換
在進行回歸分析時,通常都不對定量因變量和自變量做任何變換。然而,由基本常識和統(tǒng)計學(xué)知識可知,這樣的建模結(jié)果往往不夠理想。
2.1.3 僅對定量自變量進行變換
善于思考問題的分析者會依據(jù)探索性分析結(jié)果,對某些定量自變量進行合適的變量變換,以獲得更好的回歸建模效果。包括單一變換,如進行“對數(shù)變換、平方根變換、指數(shù)變換等”之一的變換;或引入派生變量,即同時使用多項變量變換的結(jié)果,如引入某變量的“平方項、立方項、交叉乘積項等”。
2.1.4 僅對定量因變量進行變換
在進行簡單直線回歸分析和/或多重回歸分析之前,應(yīng)考察“誤差項應(yīng)服從正態(tài)分布”的前提條件是否成立。若不成立,就需要通過探索性分析了解定量因變量可能的分布規(guī)律,從而采取合適的變量變換方法,以使其基本滿足或近似滿足進行相應(yīng)回歸分析的前提條件。然而,原始數(shù)據(jù)中并沒有“誤差項”這個變量,只能先假定資料滿足回歸分析的前提條件,基于創(chuàng)建的回歸模型計算出因變量的預(yù)測值,再計算出各觀測點上的“殘差”,用此“殘差”取代統(tǒng)計學(xué)理論上所指的“誤差”,最后再檢驗“誤差”是否服從正態(tài)分布。
由此可知:在很多實際問題中,有必要對定量因變量進行某些可能的變量變換,以獲得更好的擬合效果。
2.1.5 同時對定量自變量和因變量進行變換
事實上,在對數(shù)據(jù)進行回歸分析之前,除了對定性自變量進行必要的變量變換以外,還應(yīng)進行一系列探索性分析,以了解定量因變量與定量自變量各自的分布情況,以及定量因變量與各定量自變量之間的相互關(guān)系和變化趨勢,以便對它們分別選擇不同的變量變換方法。目的是獲得在專業(yè)上和統(tǒng)計學(xué)上都成立的最佳回歸模型。
2.2.1對定量自變量進行多種變量變換,以便產(chǎn)生派生變量
所需要的SAS程序如下:
/*在數(shù)據(jù)集a1的基礎(chǔ)上增加定量自變量的各種派生變量18個,形成數(shù)據(jù)集a2*/
data a2;
set a1;
x1=log(cpration);x2=sqrt(cpration);x3=exp(cpration);
x4=cpratio**2;x5=x4*cpratio;
w1=log(eqration);w2=sqrt(eqration);w3=exp(eqration);
w4=eqratio**2;w5=w4*eqratio;
z1=log(mfuel);z2=sqrt(mfuel);z3=exp(mfuel);
z4=mfuel**2;z5=z4*mfuel;
m1=cpration*eqration;m2=cpration*mfuel;
m3=eqration*mfuel;
run;
運行以上SAS程序后,就創(chuàng)建了數(shù)據(jù)集a2,它在數(shù)據(jù)集a1基礎(chǔ)上增加了由三個定量自變量“cpratio”“eqratio”和“mfuel”派生出來的18個新自變量,它們分別是每個定量自變量的自然對數(shù)變換、平方根變換、指數(shù)變換、平方變換和立方變換的結(jié)果;還有三個定量自變量兩兩交叉乘積變換的結(jié)果。
2.2.2 對定量因變量進行5種變量變換
所需要的SAS程序如下:
/*在數(shù)據(jù)集a2的基礎(chǔ)上增加定量因變量的5種變量變換結(jié)果,形成數(shù)據(jù)集a3*/
data a3;
set a2;
y1=log(nox);y2=sqrt(nox);y3=exp(nox);
y4=1/nox;y5=exp(nox)/(1+exp(nox));
run;
運行以上SAS程序后,就創(chuàng)建了數(shù)據(jù)集a3,它在數(shù)據(jù)集a2基礎(chǔ)上增加了由定量因變量“nox”派生出來的5個新因變量,它們分別是自然對數(shù)變換(y1)、平方根變換(y2)、指數(shù)變換(y3)、倒數(shù)變換(y4)和Logistic變換(y5)的結(jié)果。
【說明】在以下的建模策略中,先對多值名義自變量進行“算術(shù)均值變換”;然后,在下面的每種情形中都將分別在“包含截距項”與“不含截距項”的條件下,分別采取“前進法”“后退法”和“逐步法”篩選自變量。
2.3.1 以“氧化氮釋放量(nox)”為定量因變量
以“cpratio”“eqratio”和“mfuel”為三個定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型1”與“模型2”。
以“cpratio”“eqratio”“mfuel”和18個派生變量為定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型3”與“模型4”。
2.3.2以“氧化氮釋放量的自然對數(shù)變換結(jié)果(y1)”為定量因變量
以“cpratio”“eqratio”和“mfuel”為三個定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型5”與“模型6”。
以“cpratio”“eqratio”“mfuel”和18個派生變量為定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型7”與“模型8”。
2.3.3以“氧化氮釋放量的平方根變換結(jié)果(y2)”為定量因變量
以“cpratio”“eqratio”和“mfuel”為三個定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型9”與“模型10”。
以“cpratio”“eqratio”“mfuel”和18個派生變量為定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型11”與“模型12”。
2.3.4以“氧化氮釋放量的指數(shù)變換結(jié)果(y3)”為定量因變量
以“cpratio”“eqratio”和“mfuel”為三個定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型13”與“模型14”。
以“cpratio”“eqratio”“mfuel”和18個派生變量為定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型15”與“模型16”。
2.3.5以“氧化氮釋放量的倒數(shù)變換結(jié)果(y4)”為定量因變量
以“cpratio”“eqratio”和“mfuel”為三個定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型17”與“模型18”。
以“cpratio”“eqratio”“mfuel”和18個派生變量為定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型19”與“模型20”。
2.3.6以“氧化氮釋放量的Logistic變換結(jié)果(y5)”為定量因變量
以“cpratio”“eqratio”和“mfuel”為三個定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型21”與“模型2”。
以“cpratio”“eqratio”“mfuel”和18個派生變量為定量自變量建模,在“包含截距項”與“不含截距項”的條件下,選取最佳模型編號分別為“模型23”與“模型24”。
以摘要形式呈現(xiàn)選出的24個擬合效果較好的回歸模型。見表1。
表1 反映24個多重回歸模型擬合優(yōu)度的計算結(jié)果
注:第1組模型對應(yīng)的因變量為“氧化氮釋放量(nox)”;第2組模型對應(yīng)的因變量為“氧化氮釋放量的自然對數(shù)變換結(jié)果(y1)”;第3組模型對應(yīng)的因變量為“氧化氮釋放量的平方根變換結(jié)果(y2)”;第4組模型對應(yīng)的因變量為“氧化氮釋放量的指數(shù)變換結(jié)果(y3)”;第5組模型對應(yīng)的因變量為“氧化氮釋放量的倒數(shù)變換結(jié)果(y4)”;第6組模型對應(yīng)的因變量為“氧化氮釋放量的Logistic變換結(jié)果(y5)”
3.2.1第1組模型的擬合效果評價
第1組模型對應(yīng)的因變量為“氧化氮釋放量”,模型1與模型2都是僅基于3個定量自變量進行變量篩選,其區(qū)別在于模型1假定包含截距項,而模型2假定不含截距項;模型3與模型4都是基于3個定量自變量及其18個派生變量進行變量篩選,其區(qū)別在于模型3假定包含截距項,而模型4假定不含截距項。由表1中前4行結(jié)果可知:模型2優(yōu)于模型1、模型4優(yōu)于模型3,即在相同情況下,假定不含截距項的擬合結(jié)果優(yōu)于假定包含截距項的擬合結(jié)果;進一步比較可知:模型4優(yōu)于模型2,即引入派生變量的擬合結(jié)果優(yōu)于不引入派生變量的擬合結(jié)果。
3.2.2 第2組模型的擬合效果評價
第2組模型對應(yīng)的因變量為“氧化氮釋放量的自然對數(shù)變換結(jié)果(y1)”,模型5與模型6都是僅基于3個定量自變量進行變量篩選,其區(qū)別在于模型5假定包含截距項,而模型6假定不含截距項;模型7與模型8都是基于3個定量自變量及其18個派生變量進行變量篩選,其區(qū)別在于模型7假定包含截距項,而模型8假定不含截距項。由表1中第5~8行結(jié)果可知:模型6優(yōu)于模型5、模型8優(yōu)于模型7,即在相同情況下,假定不含截距項的擬合結(jié)果優(yōu)于假定包含截距項的擬合結(jié)果;進一步比較可知:模型8優(yōu)于模型6,即引入派生變量的擬合結(jié)果優(yōu)于不引入派生變量的擬合結(jié)果。
3.2.3 第3組模型的擬合效果評價
第3組模型對應(yīng)的因變量為“氧化氮釋放量的平方根變換結(jié)果(y2)”,模型9與模型10都是僅基于3個定量自變量進行變量篩選,其區(qū)別在于模型9假定包含截距項,而模型10假定不含截距項;模型11與模型12都是基于3個定量自變量及其18個派生變量進行變量篩選,其區(qū)別在于模型11假定包含截距項,而模型12假定不含截距項。由表1中第9~12行結(jié)果可知:模型10優(yōu)于模型9、模型12優(yōu)于模型11,即在相同情況下,假定不含截距項的擬合結(jié)果優(yōu)于假定包含截距項的擬合結(jié)果;進一步比較可知:模型12優(yōu)于模型10,即引入派生變量的擬合結(jié)果優(yōu)于不引入派生變量的擬合結(jié)果。
3.2.4 第4組模型的擬合效果評價
第4組模型對應(yīng)的因變量為“氧化氮釋放量的指數(shù)變換結(jié)果(y3)”,模型13與模型14都是僅基于3個定量自變量進行變量篩選,其區(qū)別在于模型13假定包含截距項,而模型14假定不含截距項;模型15與模型16都是基于3個定量自變量及其18個派生變量進行變量篩選,其區(qū)別在于模型15假定包含截距項,而模型16假定不含截距項。由表1中第13~16行結(jié)果可知:模型14優(yōu)于模型13、模型16優(yōu)于模型15,即在相同情況下,假定不含截距項的擬合結(jié)果優(yōu)于假定包含截距項的擬合結(jié)果;進一步比較可知:模型16優(yōu)于模型14,即引入派生變量的擬合結(jié)果優(yōu)于不引入派生變量的擬合結(jié)果。
3.2.5 第5組模型的擬合效果評價
第5組模型對應(yīng)的因變量為“氧化氮釋放量的倒數(shù)變換結(jié)果(y4)”,模型17與模型18都是僅基于3個定量自變量進行變量篩選,其區(qū)別在于模型17假定包含截距項,而模型18假定不含截距項;模型19與模型20都是基于3個定量自變量及其18個派生變量進行變量篩選,其區(qū)別在于模型19假定包含截距項,而模型20假定不含截距項。由表1中第17~20行結(jié)果可知:模型18優(yōu)于模型17、模型20優(yōu)于模型19,即在相同情況下,假定不含截距項的擬合結(jié)果優(yōu)于假定包含截距項的擬合結(jié)果;進一步比較可知:模型20優(yōu)于模型18,即引入派生變量的擬合結(jié)果優(yōu)于不引入派生變量的擬合結(jié)果。
3.2.6 第6組模型的擬合效果評價
第6組模型對應(yīng)的因變量為“氧化氮釋放量的Logistic變換結(jié)果(y5)”,模型21與模型22都是僅基于3個定量自變量進行變量篩選,其區(qū)別在于模型21假定包含截距項,而模型22假定不含截距項;模型23與模型24都是基于3個定量自變量及其18個派生變量進行變量篩選,其區(qū)別在于模型23假定包含截距項,而模型24假定不含截距項。由表1中第21~24行結(jié)果可知:模型22優(yōu)于模型21、模型24優(yōu)于模型23,即在相同情況下,假定不含截距項的擬合結(jié)果優(yōu)于假定包含截距項的擬合結(jié)果;進一步比較可知:模型24優(yōu)于模型22,即引入派生變量的擬合結(jié)果優(yōu)于不引入派生變量的擬合結(jié)果。
從以上的“評價結(jié)果”可知:模型4、模型8、模型12、模型16、模型20和模型24分別是從6組模型中挑選出來的“最優(yōu)模型”,現(xiàn)將它們從表1中摘錄出來,以便直觀比較和判斷。見表2。
表2 各組挑選出來的6個“最優(yōu)”多重回歸模型擬合優(yōu)度的計算結(jié)果
由表2可知:模型24是6個“最優(yōu)”模型中“最佳”的。該模型的因變量為“氧化氮釋放量的Logistic變換結(jié)果(y5)”,從全部(3+18=21個)自變量中篩選出了12個具有統(tǒng)計學(xué)意義的自變量,模型中不含截距項。具體計算結(jié)果如下:
方差分析源自由度平方和均方FPr>F模型12126.0596110.5049717367.6<0.0001誤差1570.094960.00060486未校正合計169126.15458
變量參數(shù)估計值標準誤差I(lǐng)I 型 SSFPr>FCpRatio-0.305280.151010.002474.090.0449EqRatio-294.7077832.775120.0489080.85<0.0001x21.516950.692470.002904.800.0300x40.004580.002060.003004.970.0273w187.098569.411430.0518085.65<0.0001w2-349.8387438.958570.0487780.64<0.0001w3284.6897431.703120.0487780.64<0.0001w5-130.5704714.465230.0492881.48<0.0001z20.342150.078430.0115119.03<0.0001z50.002360.000988630.003465.710.0180m1-0.024140.002950.0403766.74<0.0001m3-0.132250.018390.0312951.72<0.0001
輸出以上結(jié)果的“SAS過程步程序”如下:
/*模型24:R2=0.9992,調(diào)整R2=0.9992,MSE=0.00060486,Cp=14.9351,niv=12,無截距項*/
proc reg data=a3;
model y5=cpratio eqratio mfuel x1-x5 w1-w5 z1-z5 m1-m3/noint selection=backward sls=0.05 r;
/*模型24*/
run;
在對定量因變量構(gòu)建多重回歸模型的過程中,摒棄了傳統(tǒng)統(tǒng)計思維下的理論和方法(對定量因變量和定量自變量保持一次方形式,即不做任何變量變換,也不產(chǎn)生派生變量;對多值名義自變量進行啞變量變換,構(gòu)建所謂的“多重線性回歸模型”),而引入了動態(tài)統(tǒng)計思維下的理論和方法(對定量因變量分別采取不做變量變換和進行對數(shù)變換、平方根變換、指數(shù)變換、倒數(shù)變換和Logistic變換);淘汰了對“定量自變量不做任何變換”和“永遠固定為一次方形式”的僵化思維,不僅對其做“對數(shù)變換、平方根變換和指數(shù)變換”,還引入了“平方項、立方項和交叉乘積項”;本文提出了一種新的變量變換方法,即對“多值名義自變量”進行“算術(shù)均值變換”,不僅將其變換成“定量自變量”,還產(chǎn)生出多項派生變量。
由表1和表2可知:基于傳統(tǒng)統(tǒng)計思維創(chuàng)建的回歸模型擬合效果非常差(模型編號分別為1、5、9、13、17、21),而基于動態(tài)統(tǒng)計思維創(chuàng)建的回歸模型擬合效果很好(表2中除模型16之外)。其中,“算術(shù)均值變換”“引入派生變量”“假定回歸模型中不含截距項”和“找到定量因變量合適的變量變換方法(就本文實例而言,除了‘指數(shù)變換’外,其他4種變量變換方法的擬合效果都相當好,其中,最佳的是Logistic變換)”是動態(tài)統(tǒng)計思維建模策略中的“核心”。
有興趣的讀者還可以在本文的基礎(chǔ)上,對定量自變量增加其他一些變量變換(例如倒數(shù)變換、Logistic變換等),并將它們作為“派生變量”引入回歸建模的過程中,有可能獲得擬合效果更好的回歸模型。
本文構(gòu)建的是“多重參數(shù)非線性回歸模型”,從擬合優(yōu)度上來看,可與現(xiàn)代的“機器學(xué)習(xí)”[2-4]建模效果媲美;眾所周知,“機器學(xué)習(xí)”回歸建模的結(jié)果,只有“誤差”可以達到足夠小的程度這個“唯一優(yōu)勢”,而幾乎無法寫出其回歸模型。即便能寫出回歸模型,其“參數(shù)的個數(shù)”可能接近“無窮大”。從模型必須“精簡實用”且“便于呈現(xiàn)”的角度來考量,“機器學(xué)習(xí)”回歸建模效果似乎要遜色于本文介紹方法的建模效果。