回歸建模的基礎與要領(Ⅱ)
——偏態分布計量資料的變換

胡良平

(1.軍事科學院研究生院，北京 100850；2.世界中醫藥學會聯合會臨床科研統計學專業委員會，北京 100029*通信作者：胡良平，E-mail：lphu812@sina.com)

1 計量資料的分布類型

1.1 計量資料的概念

測定“身高”“體重”“胸圍”“血脂”等指標的數值所得到的資料，在統計學上都被稱為“計量資料”。它們有一個共同特點：數值可以帶小數且有度量衡單位。嚴格地說，任何一個計量變量的取值可以充滿其取值區間，只是在實際中，滿足一定的精度要求就可以了。例如，當研究者測量人的身高時，若以“厘米”為單位，通常保留到小數點后一位即可，如：165.3厘米，其中的“0.3”是估計出來的，沒有必要寫成“165.324568厘米”。若是計算的“中間結果”，可以保留到小數點后第6位(目的是盡可能減少計算過程中的舍入誤差)，但對于最終結果，一般只保留到測量工具能準確測到的下一位。

1.2 計量資料分布類型的概念[1]

將某實際問題中的一個計量變量的全部取值由小到大排序，再將它們按相等的間隔劃分為若干組，然后，統計出各組中的數據個數，即“頻數”。若用表格形式呈現此時的資料，它就被稱為“頻數分布表資料”。所謂“頻數分布”，就是“頻數”在各組段上是如何“分配”的；若用圖形形式呈現此時的資料，它就被稱為“頻數分布直方圖”。例如：圖1所顯示的直方圖就被稱為“正偏態分布的計量資料”。

圖1 30只老鼠腫瘤發展到特定尺寸所用時間(d)的頻數分布直方圖

由圖1可知：頻數最多的組位于橫坐標軸上偏向“左邊”的位置，右邊出現了較長的“尾巴”。若采用一個叫“偏度系數g1”的公式計算，得到的結果為“g1>0”，故稱具有這樣頻數分布的計量資料為“正偏態分布計量資料”。而圖2所顯示的直方圖被稱為“負偏態分布的計量資料”。

圖2 某地110名健康男性體重(kg)的頻數分布直方圖

由圖2可知：頻數最多的組位于橫坐標軸上偏向“右邊”的位置，左邊出現了較長的“尾巴”。若采用一個叫“偏度系數g1”的公式來計算，得到的結果為“g1<0”，故稱具有這樣頻數分布的計量資料為“負偏態分布計量資料”。

實際計量資料中，還有一些的頻數分布為“基本對稱”分布，即“頻數最多的組位于橫坐標軸上基本居中的位置”。若采用一個叫“偏度系數g1”的公式來計算，得到的結果為“g1≈0”，故稱具有這樣頻數分布的計量資料為“對稱分布計量資料”，見圖3。

圖3 由101名正常成年男子血清總膽固醇數據繪制的頻數分布直方圖

若圖3中的光滑曲線可用公式(1)描述，則該“曲線”被稱為“正態分布曲線”。

(1)

滿足式(1)的曲線被稱為“一般正態分布曲線”，其“均值為μ、標準差為σ”；它是單峰分布的，高峰位于橫坐標軸的正中位置；它的“偏度系數g1=0”，同時，它的“峰度系數g2=0”。若“均值為0、標準差為1”，此時的正態分布就被稱為“標準正態分布”。

1.3 將偏態分布計量資料變換為正態分布計量資料的必要性

在經典統計學中，無論是對計量資料進行假設檢驗(如Z檢驗、t檢驗、方差分析)、區間估計，還是進行簡單線性回歸分析或多重線性回歸分析，首選的方法是“參數法”。而參數法的重要前提條件之一是來自結果變量的計量資料必須服從“正態分布”。例如，文獻[2]中專門用一章篇幅詳細介紹“正態分布的統計方法”。

通常，需要先對計量資料進行正態性檢驗。當正態性檢驗得出該組計量資料服從正態分布時，可以采用相應的“參數法”對計量資料進行處理；反之，則要求采用“非參數法”處理計量資料。然而，在多因素或多自變量的情況下，常沒有合適的“非參數法”可運用。有時，人們習慣借助某種變量變換方法，希望經變換后的計量資料滿足“正態分布”的要求，再對變換后的計量資料采取“參數法”處理。

事實上，并非所有計量資料通過某種變量變換方法變換后都能符合“正態分布”要求。若能將服從偏態分布的計量資料變換為“對稱分布”的計量資料，也就很接近“參數法”的要求了。文獻[2]中給出了“為對稱性而變換”的方法。

1.4 常用的計量資料變換方法[2]

1.4.1 對數變換

當計量資料x呈正偏態分布時，對其進行對數變換可使其偏態狀況有所減弱；有時，取對數變換后，計量資料就接近正態分布了。用y表示變換后的計量資料，見式(2)和式(3)：

y= ln(x),x>0

(2)

y= ln(x+C), 部分x<0或x=0

(3)

在式(3)中，取C>max|x|，這里的x為負值，“max”為取“最大值”之意，應確保所有的“C+x”一定大于0。

1.4.2 平方根變換

當一組計量資料x的算術平均值近似等于其方差(若將此時的計量資料近似視為“計數資料”，則這樣的計數資料被認為近似服從Poisson分布的計數資料[3])時，對其進行平方根變換可使其偏態狀況有所減弱；有時，取平方根變換后，計量資料就接近正態分布了。用y表示變換后的計量資料，見式(4)和式(5)：

(4)

(5)

在式(5)中，取C≥max|x|，這里的x為負值，“max”為取“最大值”之意，應確保所有的“C+x”一定大于或等于0。

1.4.3 倒數變換

當計量資料x呈負偏態分布時，對其進行倒數變換可使其偏態狀況減弱；有時，取倒數變換后，計量資料就接近正態分布了。用y表示變換后的計量資料，見式(6)和式(7)：

y= 1/x,x>0

(6)

y= 1/(C+x), 部分x<0或x=0

(7)

在式(7)中，取C>max|x|，這里的x為負值，“max”為取“最大值”之意，應確保所有的“C+x”一定大于0。

1.4.4 Box-Cox變換(包含“冪變換”與“對數變換”)[4]

當計量資料x呈偏態分布(包括正偏與負偏兩種情形)時，對其進行Box-Cox變換可使其偏態狀況減弱；有時，經此變換后，計量資料就接近正態分布了。用y表示變換后的計量資料，見式(8)和式(9)：

(8)

y=ln(x), λ=0

(9)

Box-Cox變換的一般形式見下面的式(10)與式(11)：

(10)

(11)

在式(10)與式(11)中，取C>max|x|，這里的x為負值或0，“max”為取“最大值”之意，應確保所有的“C+x”一定大于0；g通常取值為1。

1.4.5 變量變換的效果

值得注意的是：對一組計量資料或計數資料做任何變換，都不可能絕對保證一定能使其呈“對稱分布”或“正態分布”。通常，經過合適的變量變換后，會使變換后的資料較原始資料具有更好的“對稱性”。有時，可能需要相繼采取多種變量變換方法。若目的是為了使變換后的資料接近“正態分布”，則必須對變換后的資料進行嚴格的正態性檢驗。只有通過了正態性檢驗(最好，正態性檢驗的結果為P>0.2；通常，P>0.1即可；但至少也應滿足P>0.05)的資料，才適合選用相應的參數統計分析方法(如t檢驗、方差分析、簡單線性回歸分析或多重線性回歸分析，在回歸分析中，應特別強調：因變量應近似服從正態分布；然而，在統計理論上，假定“模型的誤差項服從正態分布”)。

以下基于SAS中的“TRANSREG過程”[4]并采用“Box-Cox變換”將偏態分布計量資料變換為近似呈正態分布的計量資料。

2 將偏態分布計量資料變換為正態分布計量資料

2.1 將正偏態分布計量資料變換為正態分布計量資料

2.1.1問題與數據結構

【例1】給30只老鼠注射給定的腫瘤接種物，腫瘤發展到特定尺寸所用的時間(d)如下：

1.7、3.7、5.0、5.1、5.3、5.9、6.0、6.0、7.4、8.0、8.3、8.3、8.3、9.1、9.6、11.3、12.1、12.3、13.1、13.4、14.0、15.9、16.1、16.7、17.0、21.0、22.7、30.0

2.1.2 所需的SAS程序

利用下面的SAS程序創建SAS數據集并進行Box-Cox變換：

/*以下的SAS數據步程序用于創建SAS數據集a1*/

data a1;

inputy@@;

z=0;

cards;

1.7 3.7 5.0 5.1 5.3 5.9 6.0 6.0 7.4 8.0

8.3 8.3 8.3 9.1 9.6 11.3 12.1 12.3 13.1

13.4 14.0 15.9 16.1 16.7 17.0 21.0 22.7 30.0

;

run;

/*以下程序繪制a1的頻數分布直方圖并進行正態性檢驗*/

proc univariate data=a1 normal;

vary;

histogramy;

run;

/*以下程序對a1進行Box-Cox變換，求出合適的lambda值*/

/*經過Box-Cox變換后的數據存儲在數據集aaa中*/

ods graphics on;

proc transreg details data=a1 maxiter=0 nozerocon

stant plots=(transformation(dependent) obp);

model BoxCox(y/ convenient lambda=-10 to 10 by 0.01)=identity(z);

output out=aaa approximations;

run;

/*以下程序繪制aaa的頻數分布直方圖并進行正態性檢驗*/

proc univariate data=aaa normal;

varty;

histogramty/normal;

run;

2.1.3 輸出結果及解釋

2.1.3.1 反映原始數據分布狀況的結果

本例中原始數據的頻數直方圖見前面的圖1(呈正偏態分布)，此處從略。對原始數據進行正態性檢驗的結果為：W=0.925501、P=0.0475，說明原始數據不服從正態分布。

偏度系數與峰度系數分別為g1=1.097與g2=1.415，表明原始數據具有正偏態(偏度系數明顯大于0)和尖翹峰(峰度系數明顯大于0)分布。

2.1.3.2 對原始數據作Box-Cox變換的結果

求得公式(8)中的參數λ為“0.29”，對經Box-Cox變換后的數據作正態性檢驗，得到：W=0.971648、P=0.6254，說明經Box-Cox變換后的數據服從正態分布。繪制經Box-Cox變換后數據的頻數直方圖，見圖4。

圖4 經Box-Cox變換后的30只老鼠腫瘤發展到特定尺寸所用時間(d)的頻數分布直方圖

2.2 將負偏態分布計量資料變換為正態分布計量資料

2.2.1 問題與數據結構

【例2】某研究者收集到某地110名健康成年男性的體重(kg)數據如下：

43.5、70.0、45.0、45.0、46.5、69.5、58.0、68.0、66.5、70.1、67.0、66.5、68.0、59.0、66.0、68.0、69.8、68.8、67.0、55.5、51.5、52.5、61.0、58.0、47.5、53.0、53.0、54.0、59.0、54.0、46.0、54.0、55.0、57.0、52.0、52.0、50.0、54.0、62.5、54.5、65.0、61.5、60.5、60.0、68.5、67.0、70.0、67.0、75.0、70.5、64.5、68.0、72.0、63.0、63.5、64.0、65.0、74.5、72.5、67.5、72.0、69.0、61.5、69.0、60.0、40.0、71.2、74.0、71.0、69.5、69.0、61.0、70.3、68.5、64.5、70.5、73.0、65.0、67.5、71.0、79.0、80.0、72.5、79.0、75.2、81.0、82.0、75.0、73.0、77.0、80.0、81.5、42.0、77.0、75.3、81.0、83.0、74.5、80.0、75.5、80.0、77.0、75.0、61.0、79.0、81.5、78.0、73.5、81.9、85.0

試呈現原始數據的分布情況，并對其進行變量變換，使其接近正態分布。

2.2.2 所需要的SAS程序

利用下面的SAS程序創建SAS數據集并進行Box-Cox變換：

/*以下的SAS數據步程序用于創建SAS數據集a1*/

/*在原始數據中增加一個新變量z，它是原始數據*/

/*取倒數變換后再乘以1000得到的結果*/

data a1;

inputy@@;

z=(1/y)*1000;

cards;

43.5 70.0 45.0 45.0 46.5 69.5 58.0 68.0 66.5 70.1 67.0 66.5 68.0 59.0 66.0 68.0 69.8 68.8 67.0 55.5 51.5 52.5 61.0 58.0 47.5 53.0 53.0 54.0 59.0 54.0 46.0 54.0 55.0 57.0 52.0 52.0 50.0 54.0 62.5 54.5 65.0 61.5 60.5 60.0 68.5 67.0 70.0 67.0 75.0 70.5 64.5 68.0 72.0 63.0 63.5 64.0 65.0 74.5 72.5 67.5 72.0 69.0 61.5 69.0 60.0 40.0 71.2 74.0 71.0 69.5 69.0 61.0 70.3 68.5 64.5 70.5 73.0 65.0 67.5 71.0 79.0 80.0 72.5 79.0 75.2 81.0 82.0 75.0 73.0 77.0 80.0 81.5 42.0 77.0 75.3 81.0 83.0 74.5 80.0 75.5 80.0 77.0 75.0 61.0 79.0 81.5 78.0 73.5 81.9 85.0

;

run;

/*以下程序繪制a1中原始數據y和倒數變換后z的頻數分布直方圖并進行正態性檢驗*/

proc univariate data=a1 normal;

varyz;

histogramyz;

run;

/*以下程序在數據集a1中增添w=0的一列，形成數據集a2*/

data a2;

set a1;

w=0;

run;

/*以下程序對a2進行Box-Cox變換，求出合適的lambda值*/

/*經過Box-Cox變換后的數據存儲在數據集aaa中*/

ods graphics on;

proc transreg details data=a2 maxiter=0 nozerocon

stant plots=(transformation(dependent) obp);

model BoxCox(z/ convenient lambda=-10 to 10 by 0.05)=identity(w);

output out=aaa approximations;

run;

/*以下程序繪制aaa的頻數分布直方圖并進行正態性檢驗*/

proc univariate data=aaa normal;

var tz;

histogram tz/normal;

run;

2.2.3 輸出結果及解釋

2.2.3.1 反映原始數據分布狀況的結果

本例中的原始數據y的頻數直方圖見前面的圖2(呈負偏態分布)，此處從略。對原始數據y進行正態性檢驗的結果為：W=0.967065、P=0.0080，說明原始數據y不服從正態分布。

偏度系數與峰度系數分別為g1=-0.495，g2=-0.385，表明原始數據具有負偏態(偏度系數明顯小于0)和平闊峰(峰度系數明顯小于0)分布。

本例中的原始數據經倒數變換后的數據z的頻數分布直方圖見圖5(呈正偏態分布)。

圖5 本例數據經倒數變換后的數據z的頻數分布直方圖

對經倒數變換后的數據z進行正態性檢驗的結果為：W=0.892042、P<0.0001，說明經倒數變換后的數據z不服從正態分布。

偏度系數與峰度系數分別為g1=1.237，g2=1.242，表明經倒數變換后的數據z具有正偏態(偏度系數明顯大于0)和尖翹峰(峰度系數明顯大于0)分布。

2.2.3.2對經倒數變換后的數據z作Box-Cox變換的結果

求得公式(8)中的參數λ為“-2.2”，對經Box-Cox變換后的數據作正態性檢驗，得到：W=0.980521、P=0.1079，說明經Box-Cox變換后的數據服從正態分布。

繪制經Box-Cox變換后的數據的頻數分布直方圖，見圖6。

圖6 經Box-Cox變換后的某地110名健康男性體重(kg)的頻數分布直方圖