——優(yōu)化計分變換與其他變量變換"/>
胡良平
(1.軍事醫(yī)學(xué)科學(xué)院研究生院,北京 100850;2.世界中醫(yī)藥學(xué)會聯(lián)合會臨床科研統(tǒng)計學(xué)專業(yè)委員會,北京 100029
1.1.1優(yōu)化計分變換
優(yōu)化計分變換就是采用Fisher[2]提出的給多值名義變量或多值有序變量各水平賦值的方法。Opscore(x)代表對變量x進行優(yōu)化計分變換,它可以被用于“多值名義變量”或“多值有序變量”。
1.1.2 單調(diào)變換
單調(diào)變換就是采用Kruskal提出的給多值有序變量各水平賦值的方法(被稱為次要最小平方單調(diào)變換)。Monotone(x)代表對變量x進行單調(diào)變換,它只能被用于“多值有序變量”。
1.1.3 “優(yōu)化計分變換”與“單調(diào)變換”的具體做法
兩者具體的變換方法是相同的,每一個不同的非缺失值形成一個不同的類,例如:“1,1,1,2,2,3”形成三類,即3個“1”算作一類、兩個“2”算作一類、一個“3”算作一類。具體賦值方法參見下面的實例:
設(shè)變量x的具體取值為( . . .A .A .B 1 1 1 2 2 3 3 3 4)' ;設(shè)變量y的具體取值為( 5 6 2 4 2 1 2 3 4 6 4 5 6 7)' 。
【說明】以上“x”與“y”都有14個取值,用“向量”形式表示,即這兩個向量都有14個分量。其中,變量x的前兩個分量都是缺失值“.”;第3和第4個分量都是“.A”;第5個分量是“.B”;第6~14個分量都是數(shù)字。而變量y的全部14個分量都是數(shù)字。
于是,以“y”為“因變量”,對“變量x”進行“優(yōu)化計分變換”和“單調(diào)變換”。用關(guān)鍵詞表示為:Opscore(x)和Monotone(x)。它們只有8個分量,即:
Opscore(x)=( 5 6 3 2 2 5 5 7)'
Monotone(x)=( 5 6 3 2 2 5 5 7)'
【說明】在“x”的向量中,最前面的兩個缺失值各形成一類,即有兩類;“.A”出現(xiàn)了兩次,形成一類;“.B”形成一類 ;三個“1”形成一類;兩個“2”形成一類;三個“3”形成一類;一個“4”形成一類,故總共形成8類。
在同一類中,求出“對應(yīng)于因變量y全部數(shù)值的算術(shù)平均值”,賦值給變換后的“變量x”,作為其“優(yōu)化計分變換”或“單調(diào)變換”的結(jié)果,例如:兩個“.A”對應(yīng)著變量y的兩個數(shù)值為2和4,其算術(shù)平均值為3;同理,三個“1”對應(yīng)著變量y的三個數(shù)值為“1、2、3”,其算術(shù)平均值為2;兩個“2”對應(yīng)著變量y的兩個數(shù)值為“4與6”,其算術(shù)平均值為5;三個“3”對應(yīng)著變量y的三個數(shù)值為“4、5、6”,其算術(shù)平均值為5;一個“4”對應(yīng)著變量y的一個數(shù)值為“7”,其算術(shù)平均值為7。
1.2.1 變量擴展
1.2.2 非優(yōu)化變換
非優(yōu)化變換就是用一個變換后的“新變量”取代“原變量”。這個“新變量”在后續(xù)的迭代算法中不再被重新變換了(對具有缺失值估計所做的可能的線性變換除外)。常用的“非優(yōu)化變換”有如下幾種。
ARSIN(x):此處的“x”為“百分率”數(shù)據(jù),通常其取值區(qū)間為“0 EXP(x):指數(shù)變換,通常為“ex”變換;若希望底數(shù)為任意的正實數(shù)a,即做“ax”變換,在使用SAS的“TRANSREG過程”中,需要另外再使用“Parameter=”來指定,即在“=”后寫出“a”的具體數(shù)值。 LOG(x):對數(shù)變換,通常取以“e”為底數(shù)的對數(shù)(e=2.718281828),即取自然對數(shù)變換;若希望底數(shù)為任意的正實數(shù)a,即做“l(fā)ogax”變換,在使用SAS的“TRANSREG過程”中,需要另外再使用“Parameter=”來指定,即在“=”后寫出“a”的具體數(shù)值。 LOGIT(x):logit變換,即做“l(fā)oga[x/(1-x)]”變換,這里的“a”通常為“e=2.718281828”。注意:x為數(shù)值型變量,其定義域為(0.0,1.0)。 POWER(x):冪變換,即做“xa”變換。例如:希望做“x1.5”變換,關(guān)鍵詞寫為POWER(x/PARAMETER=1.5)。顯然,當(dāng)“PARAMETER=”后數(shù)值為“2”時,就是平方變換;為“-1”時,就是倒數(shù)變換;為“0.5”時,就是開平方根變換。 RANK(x):秩變換。將變量x由小到大排序,再分別賦值為“1、2、3、4、5、……”。 1.2.3 非線性擬合變換 BOXCOX(x):這是Box等[5]提出的變量變換方法,實際應(yīng)用參見文獻[6]。此變換只適用于回歸分析中的“定量因變量”。 PBSPLINE(x):這是非迭代懲罰B樣條變換,使用時的特殊要求參見文獻[1]。此變換只適用于回歸分析中的“定量自變量”。 SMOOTH(x):這是一種非迭代光滑樣條變換,由Reinsch[7]提出的變量變換方法,使用時的特殊要求參見文獻[1]。此變換只適用于回歸分析中的“定量自變量”。 1.2.4 優(yōu)化變換 優(yōu)化變換就是采用迭代計算導(dǎo)出的變量變換方法。在SAS的“TRANSREG過程”中,涉及到以下六種優(yōu)化變換方法,即:LINEAR(x)、MONOTONE(x)、MSPLINE(x)、OPSCORE(x)、SPLINE(x)和UNTIE(x)。以上六種優(yōu)化變換在使用時都有具體要求,詳見文獻[1],此處從略。 1.3.1 恒等變換[IDENTITY(x)] 恒等變換,就是沒有改變變量的取值。通常情況下,就是“不做變換”。應(yīng)注意:在使用SAS的“TRANSREG過程”時,寫在“model語句”中的所有“因變量”和“自變量”之前必需要有“變量變換”的關(guān)鍵詞。若不想對“x4-x10”進行任何變量變換,必須按如下的方式寫:IDENTITY(x4-x10)。 使用此種變換,可以產(chǎn)生“交互項”,具體使用方法參見文獻[1] ,此處從略。 1.3.2 迭代光滑樣條變換[SSPLINE(x)] 這是對變量x進行迭代光滑樣條變換,此變量變換方法通常不會使誤差平方達到最小值。具體使用方法參見文獻[1] ,此處從略。 沿用本期科研方法專題第一篇文章《提高回歸模型擬合優(yōu)度的策略(Ⅰ)——啞變量變換與其他變量變換》中的“實際問題與數(shù)據(jù)結(jié)構(gòu)”[1],此處從略。 3.1.1所需要的SAS程序 對“燃油種類(fuel)” 這個“6值名義自變量”進行“優(yōu)化計分變換”;對“壓縮比(cpratio)” 這個“5值有序自變量”進行“單調(diào)變換”。所需要的SAS程序如下: ods graphics on; title' Gasoline Example'; title2' Iteratively Estimate NOx, CpRatio, EqRatio, and Fuel'; * Fit the Nonparametric Model; proc transreg data=sashelp.Gas solve test nomiss plots=all; ods exclude where=(_path_ ? ' MV' ); model mspline(NOx / nknots=9) = spline(EqRatio / nknots=9) monotone(CpRatio) opscore(Fuel); output out=aaa tdprefix=td tiprefix=ti; run; proc freq data=aaa; tables cpratio ticpratio fuel tifuel; run; 3.1.2 與所需變量變換有關(guān)的SAS輸出結(jié)果 Compression RatioCpRatio頻數(shù)百分比累積頻數(shù)累積百分比7.59354.399354.399179.9411064.33122414.0413478.36152011.7015490.0618179.94171100.00 以上是“壓縮比(cpratio)”的原始取值及其各水平的頻數(shù)分布。 Compression Ratio TransformationtiCpRatio頻數(shù)百分比累積頻數(shù)累積百分比7.15982275949354.399354.3910.953020684179.9411064.3312.9445525662414.0413478.3614.201537262011.7015490.0617.433541388179.94171100.00 以上是“壓縮比(cpratio)”經(jīng)過“單調(diào)變換”后的取值及其各水平的頻數(shù)分布。 由以上兩部分輸出結(jié)果可知:“單調(diào)變換”是將“多值有序變量”的水平做一一對應(yīng)的變換,各水平出現(xiàn)的“頻數(shù)”不會發(fā)生變化。 Fuel頻數(shù)百分比累積頻數(shù)累積百分比095.2695.2612514.623419.8829052.6312472.513137.6013780.1242212.8715992.985127.02171100.00 以上是“燃油種類(fuel)”的“6種水平代碼(分別為0、1、2、3、4、5)”及其各水平的頻數(shù)分布。 Fuel頻數(shù)百分比累積頻數(shù)累積百分比82rongas95.2695.2694%Eth2514.623419.88Ethanol9052.6312472.51Gasohol137.6013780.12Indolene2212.8715992.98Methanol127.02171100.00 以上是“燃油種類(fuel)”的“6種水平的真實含義”及其各水平的頻數(shù)分布。 Fuel TransformationtiFuel頻數(shù)百分比累積頻數(shù)累積百分比1.3742357556127.02127.021.55720545762514.623721.641.59754088319052.6312774.274.209335677695.2613679.534.36641535932212.8715892.404.4654059829137.60171100.00 以上是“燃油種類(fuel)” 經(jīng)過“優(yōu)化計分變換”后的取值及其各水平的頻數(shù)分布。 對定量自變量可以進行“樣條變換(spline)”“單調(diào)樣條變換(mspline)”或“變量擴展變換(pspline)(即產(chǎn)生‘派生變量’,例如變量的平方項、立方項等)”。例如:對“eqratio”進行“樣條變換”,可以寫成spline(eqratio)或spline(eqratio/degree=3 nknots=0)。后一種寫法涉及到兩個參數(shù),一是“degree=”,它是指“樣條函數(shù)(通常是多項式)”中多項式的“次數(shù)”,“degree=3”代表“三次多項式”;二是“nknots=”,它是指“結(jié)點(或?qū)懗伞?jié)點’)個數(shù)”,實際上就是“斷點個數(shù)”。其具體含義是:在該定量變量的取值區(qū)間內(nèi),確定“幾個斷點(如nknots=9)”,于是,就將該區(qū)間劃分成“10段”。在每一段上,擬合一個“多項式”。就整體而言,被稱為“分段多項式(a piecewise polynomial)”。 若需要變換的定量變量為“因變量”,除了可以進行“定量自變量”的某些變量變換(注意:“非迭代懲罰B樣條變換”只能用于“定量自變量”)之外,還可以進行其他一些變量變換,如“BOX-COX變換(注意:該變換僅適用于‘定量因變量’)”等。 3.4.1尋找最優(yōu)回歸模型的策略 前面提及的“實際問題”涉及“定量因變量(nox)”“定量自變量(eqratio)”“多值有序自變量(cpratio)”和“多值名義自變量(fuel)”。在擬合多重回歸模型的過程中,涉及到如何對上述4個變量進行“變量變換”。當(dāng)然,對每個變量選取不同的“變量變換”方法,將會得到不同的擬合效果。下面將固定“多值名義自變量(fuel)”的變量變換方法為:優(yōu)化計分變換,即“OPSCORE(fuel)”,對其他變量依次采取各種可能的變量變換,呈現(xiàn)出各種情況下回歸模型的擬合優(yōu)度,最終選擇“擬合優(yōu)度最好”的回歸模型。 3.4.2不同變量變換方法組合下的回歸模型擬合優(yōu)度 在對“多值名義自變量(fuel)”做優(yōu)化計分變換“OPSCORE(fuel)”的前提條件下,對其他變量依次采取各種可能的變量變換,對應(yīng)的回歸模型擬合優(yōu)度見表1。 表1 在OPSCORE(fuel)的前提下對“nox”“eqratio”和“cpratio”進行不同變量變換后對應(yīng)的回歸模型的擬合優(yōu)度 注:“Root MSE”代表“均方誤差的平方根”,此值越小越好;“Nox_B”代表對“Nox”進行變量變換的方法,表中第3列和第4列符號代表的含義與“Nox_B”相同,此處從略 由表1可知:模型1和模型7的結(jié)果相同且擬合優(yōu)度最好,對應(yīng)的SAS過程步程序如下: proc transreg data=sashelp.Gas solve test nomiss plots=all; ods exclude where=(_path_ ?' MV' ); model spline(NOx / nknots=9) =spline(EqRatio / nknots=9) monotone(CpRatio) opscore(Fuel); run; 3.4.3 模型1的主要輸出結(jié)果 由于“模型1”是一個“非參數(shù)模型”,只能給出總模型的有關(guān)信息和與“擬合優(yōu)度”有關(guān)的結(jié)果,不能給出“回歸系數(shù)”等信息。但可以以“圖形”方式給出各變量變換的結(jié)果,以提示對不同變量采取什么變量變換方法最有效。下面先給出與“擬合優(yōu)度”有關(guān)的結(jié)果: Root MSE0.30935R-Square0.9586Dependent Mean2.34593Adj R-Sq0.9527Coeff Var13.18661 表1中的第1行就是摘錄了以上的結(jié)果。 在前面給出的“SAS過程步”的“MODEL語句”中,對定量因變量(nox)、定量自變量(EqRatio)、多值有序自變量(CpRatio)和多值名義自變量(Fuel)分別做了“樣條變換”“樣條變換”“單調(diào)變換”和“優(yōu)化計分變換”,變換的結(jié)果用“圖形”呈現(xiàn)出來,見圖1。 圖1 基于MODEL語句中的要求對4個變量做的變換結(jié)果 在圖1中,左上方一圖表明:定量因變量(nox)呈現(xiàn)出“對數(shù)”變化趨勢,提示可對其進行“對數(shù)變換”;右上角一圖表明:定量自變量(EqRatio)呈現(xiàn)出“二次拋物線”變化趨勢,提示可對其進行“二次多項式變換”;左下角一圖表明:多值有序自變量(CpRatio)呈現(xiàn)出近似“直線”變化趨勢,提示可對其進行“恒等變換(即不做變換)”;右下角一圖表明:多值名義自變量(Fuel)呈現(xiàn)出“兩種級別”,其中,水平代碼為“1、2、5”的燃料種類對結(jié)果的影響數(shù)量大約在“1與2”之間;而水平代碼為“0、3、4”的燃料種類對結(jié)果的影響數(shù)量大約在“4與5”之間,縱軸上的“數(shù)量”就是“優(yōu)化計分”的結(jié)果。 基于以上的“非參數(shù)多重非線性回歸分析結(jié)果(以圖1反映其變量變換的效果)”,可以構(gòu)建參數(shù)多重非參數(shù)回歸分析模型如下: log(nox) = b0+b1*EqRatio+b2*EqRatio**2+b3*CpRatio+Sum b(j)*Fuel(j)+Error 此模型的含義是:對定量因變量(nox)進行自然對數(shù)變換,對定量自變量(EqRatio)進行二次多項式變換,對多值有序自變量(CpRatio)進行恒等變換,而對多值名義自變量(Fuel)進行優(yōu)化計分變換并將其6個計分值[Fuel(j)]與其回歸系數(shù)[b(j)]分別相乘后求和(j=0、1、2、3、4、5)(注意:相當(dāng)于把“6值名義自變量”視為6個新變量,就有6個回歸系數(shù);但它們又屬于同一個多值名義變量,因此,其自由度為6-1=5,在本質(zhì)上相當(dāng)于是5個啞變量)。 3.4.4 構(gòu)建參數(shù)多重非線性回歸模型 采用SAS擬合上述參數(shù)多重非線性回歸模型所需要的SAS過程步程序如下: title2' Now fit log(nox) = b0 + b1*EqRatio + b2*EqRatio**2 +'; title3' b3*CpRatio + Sum b(j)*Fuel(j) + Error'; *-Fit the Parametric Model Suggested by the Nonparametric Analysis-; proc transreg data=sashelp.Gas solve ss2 short nomiss plots=all cldetail; model log(nox) = pspline(EqRatio / deg=2) identity(CpRatio) opscore(Fuel); run; 3.4.5 參數(shù)多重非線性回歸模型輸出結(jié)果 Univariate ANOVA Table Based on the Usual Degrees of FreedomSourceDFSum of SquaresMean SquareFPr > FModel879.338389.917298213.09<0.0001Error1607.446590.046541Corrected Total16886.78498 以上是總模型的假設(shè)檢驗結(jié)果,F(xiàn)=213.09,P<0.0001,表明總模型具有統(tǒng)計學(xué)意義。 Root MSE0.21573R-Square0.9142Dependent Mean0.63130Adj R-Sq0.9099Coeff Var34.17294 以上是擬合優(yōu)度評價結(jié)果:均方誤差平方根為0.21573,R2為0.9142,調(diào)整R2為0.9099,相對于表1中“模型1”的三個對應(yīng)的數(shù)值(0.30935、0.9586、0.9527),均方誤差平方根略微變小了一點,然而,R2和調(diào)整R2的數(shù)值卻下降得比較多(因為模型中自變量的項數(shù)減少得很多,現(xiàn)在模型的自由度df=8;而模型1的df=21)。 Univariate Regression Table Based on the Usual Degrees of FreedomVariableDFCoefficientType II Sum of SquaresMean SquareFPr > FLabelIntercept1-15.27464957.133857.13381227.60<0.0001InterceptPspline.EqRatio_1135.10291462.747862.74781348.22<0.0001Equivalence Ratio 1Pspline.EqRatio_21-19.38646864.643064.64301388.94<0.0001Equivalence Ratio 2Identity(CpRatio)10.0320581.44451.444531.04<0.0001Compression RatioOpscore(Fuel)50.1583885.56191.112423.90<0.0001Fuel 以上為回歸模型的回歸系數(shù)與假設(shè)檢驗結(jié)果(因輸出結(jié)果過寬,將其中參數(shù)的置信區(qū)間等信息省略掉了)。其中,第2、3兩行上分別為“定量自變量(EqRatio)”的“一次項”與“二次項”的計算結(jié)果;第4行上為“多值有序自變量(CpRatio)”的計算結(jié)果;而第5行上為“多值名義自變量(Fuel)”的計算結(jié)果,具有5個自由度,但只有一個回歸系數(shù)。 基于上述回歸模型計算出定量因變量的預(yù)測值作為橫坐標(biāo)的數(shù)值,其原始觀測值作為縱坐標(biāo)的數(shù)值,繪制出散布圖,以直觀的方式呈現(xiàn)模型對資料的擬合效果。見圖2。 圖2 反映因變量(nox)的觀測值與預(yù)測值的散布圖 由圖2可知:此回歸模型對資料的擬合效果比較令人滿意。 本文對“多值名義自變量”進行“優(yōu)化計分變換(Opscore)”,取代了經(jīng)典統(tǒng)計學(xué)中常規(guī)的“啞變量變換”。在此基礎(chǔ)上,對“定量因變量(nox)”分別進行了“樣條變換(Spline)”“單調(diào)樣條變換(Mspline)”和“Box-Cox變換”三種變量變換;對“定量自變量(EqRatio)”分別進行了“樣條變換(Spline)”和“單調(diào)樣條變換(Mspline)”兩種變量變換;對“多值有序自變量(CpRatio)”分別進行了“單調(diào)變換(Monotone)”和“優(yōu)化計分變換(Opscore)”。從而,獲得“最優(yōu)模型”的擬合優(yōu)度評價指標(biāo)為“均方誤差平方根=0.30935、R2=0.9586、調(diào)整R2=0.9527”。由于此模型相對比較復(fù)雜(自變量有21項且寫不出具體的回歸系數(shù));簡化后的參數(shù)多重非線性回歸模型僅含8個自變量(注:包括派生變量),R2為0.9142,調(diào)整R2為0.9099。 在本期科研方法專題的前三篇文章中,分別對“多值名義自變量”進行“啞變量變換”“算術(shù)均值變換”和“校正均值變換”,對“定量因變量(nox)”進行了“不變換”“自然對數(shù)變換”“平方根變換”“倒數(shù)變換”“指數(shù)變換”和“Logistic變換”;對“定量自變量(EqRatio)”和“多值有序自變量(CpRatio)(被視為‘定量的’)”引入了“派生變量”。還在模型中假定“包含截距項”和“不含截距項”兩種情形下,分別采取“前進法”“后退法”和“逐步法”篩選自變量。在前三篇文章中的每一篇中,都從大約72個模型中選出了24個“最優(yōu)模型”且以“模型24”為“最佳模型”。 “四組方法”指本期科研方法專題中四篇文章所介紹的方法,即分別對“多值名義自變量”進行“啞變量變換”“算術(shù)均值變換”“校正均值變換”和“優(yōu)化計分變換”。在此基礎(chǔ)上,再對其他變量采用多種變量變換方法,最后以“擬合優(yōu)度評價指標(biāo)”的取值為判斷標(biāo)準(zhǔn),總結(jié)于下表2。 表2 四篇文章中“最佳”回歸模型擬合優(yōu)度評價指標(biāo)匯總 注:前三篇文章原先給出的是“均方誤差(MSE)”,分別為“0.00076”“0.00060”和“0.00061”,將它們分別開平方根,得到“0.02757”“0.02449”和“0.02470”。由表2可知,前三篇文章回歸模型的擬合效果非常接近,均優(yōu)于第4篇文章回歸模型的擬合效果;然而,第一篇文章回歸模型中保留了“部分啞變量”,嚴(yán)格來說,它不是非常理想的結(jié)果 由此可以得出如下結(jié)論:在對“定量因變量”構(gòu)建多重非線性回歸模型中,對“多值名義自變量”是采取“啞變量變換”“優(yōu)化計分變換”還是“均值或校正均值變換”,對結(jié)果的影響都不是特別大,關(guān)鍵在于以下四點:第一,應(yīng)對“定量因變量”“定量自變量”和“多值有序自變量”采取合適的變量變換;第二,應(yīng)盡可能引入較多的“派生變量”;第三,應(yīng)假定“模型中不含截距項”;第四,應(yīng)盡可能多采取一些篩選自變量的策略,如前進法、后退法和逐步法。 最后還需要注意:引入派生變量后獲得的“最佳”回歸模型通常具有“嚴(yán)重多重共線性”,若最終的回歸模型中的某些回歸系數(shù)的“正負號”與專業(yè)知識不符合,此時,再基于“最佳回歸模型”所確定的“自變量組合”進行“嶺回歸分析”。從而獲得滿意的結(jié)果[4]。1.3 對定量變量進行其他變量變換
2 實際問題與數(shù)據(jù)結(jié)構(gòu)
3 解決問題的思路和做法
3.1 對多值名義自變量和多值有序自變量進行變量變換
3.2 對定量自變量進行變量變換
3.3 對定量因變量進行變量變換
3.4 針對實際問題,尋找合適的變量變換方法
4 總 結(jié)
4.1 本文方法的小結(jié)
4.2 其他方法的回顧
4.3 四組方法的總結(jié)