0802 回歸分析簡(jiǎn)介及其在均勻設(shè)計(jì)中的應(yīng)用

1、第二節(jié)回歸分析簡(jiǎn)介及其在均勻設(shè)計(jì)中的應(yīng)用回歸分析是數(shù)據(jù)分析的有力工具，它能揭示變量之間的相互關(guān)系，因此在均勻設(shè)計(jì)的數(shù)據(jù)分析中成為主要的手段，回歸分析方法和理論十分豐富，有關(guān)書籍?dāng)?shù)以百計(jì)，這里僅作一梗概介紹，細(xì)節(jié)可以參看有關(guān)書籍，如26，29，30數(shù)據(jù)處理可使用統(tǒng)計(jì)軟件包SAS，SPSS，MINITAB，BMDP，S等，國(guó)內(nèi)許多部門如中國(guó)均勻設(shè)計(jì)學(xué)會(huì)為均勻設(shè)計(jì)及其數(shù)據(jù)分析制作了專用統(tǒng)計(jì)軟件包，使用更為方便。2.1 一元線性回歸模型 由于均勻設(shè)計(jì)的數(shù)據(jù)分析要利用回歸分析，因此需要對(duì)回歸分析作一扼要介紹。一元線性回歸是處理兩個(gè)變量之間關(guān)系的最簡(jiǎn)單的模型。本章將詳細(xì)討論這個(gè)模型。一元線性回歸雖簡(jiǎn)單，但

2、從中可以了解回歸分析方法的基本思想/方法和應(yīng)用。 我們首先通過(guò)一個(gè)例子說(shuō)明如何建立一元線性回歸方程。 例3 為了估計(jì)山上積雪融化后對(duì)下游灌溉的影響，在山上建立了一個(gè)觀測(cè)站，測(cè)量了最大積雪深度(X)與當(dāng)年灌溉面積(Y)，得到連續(xù)10年的數(shù)據(jù)。 年序最大積雪深度X(尺)灌溉面積Y(千畝)115.228.6210.419.3321.240.5418.635.6526.448.9623.445.0713.529.2816.734.1924.046.71019.137.4 為了研究這些數(shù)據(jù)中所蘊(yùn)含的規(guī)律性，我們把各年最大積雪深度作橫坐標(biāo)，相應(yīng)的灌溉面積作縱坐標(biāo)，將這些數(shù)據(jù)點(diǎn)標(biāo)在平面直角坐標(biāo)圖上，如圖9，

3、這個(gè)圖稱為散點(diǎn)圖。 從圖9看到，數(shù)據(jù)點(diǎn)大致落在一條直線附近，這告訴我們變量X與Y之間的關(guān)系大致可看作是線性關(guān)系，從圖9還看到，這些點(diǎn)又不都在一條直線上，這表明X與Y的關(guān)系并沒(méi)有確切到給定X就可以唯一地確定Y的程度。事實(shí)上，還有許多其他因素對(duì)Y產(chǎn)生影響，如當(dāng)年的平均氣溫，當(dāng)年的降雨量等等，這些都是影響Y取什么值的隨機(jī)因素。如果我們只研究X與Y的關(guān)系，可以假定 有如下結(jié)構(gòu)式：Y=+X+ (2.1) 式中， 稱為回歸系數(shù)，X為自變量，Y為因變量，表示隨機(jī)誤差，常常假定遵從正態(tài)分布N(0，2)，這表示誤差為正和負(fù)的機(jī)會(huì)一樣多，2 表示誤差的大小。式中，2 通常是未知的，它們要通過(guò)數(shù)據(jù)的信息來(lái)估計(jì)。 設(shè)

4、(Xi,Yi)，i=1，n為一組數(shù)據(jù)，若用回歸方程(2.1)來(lái)擬合，則當(dāng)X=時(shí)的估計(jì)值為 (2.2)自然，我們希望求和使與很接近.也就是說(shuō)，我們要決定一條直線，使其與所有的點(diǎn)都比較接近，最流行求， 估計(jì)值的辦法是用最小二乘法，令 (2.3)最小二乘法是求和使Q達(dá)極小，使Q達(dá)極小的和值記為a和b.利用微積分中求極值的辦法求得 (2.4)式中 (2.5)利用這些公式到例3，得:于是 b=415.606/230.656=1.802×18.88=2.511從而回歸方程為 試將該直線畫在圖9上，可以看到擬合的效果是不錯(cuò)的，衡量擬合效果的好壞，如下的方法是十分有用的。 (a) 相關(guān)系數(shù) 相關(guān)系數(shù)

5、 用于描敘變量X和Y的線性相關(guān)的程度，并常用r來(lái)表示，r的值介于-1，1之間，它的意義由圖10可以知道。r的絕對(duì)值越接近于1表示X和Y之間的線性關(guān)系越密切；r0，兩者呈正比關(guān)系，叫正相關(guān)；r 0兩者呈負(fù)相關(guān)。r的值接近于0，兩者沒(méi)有線性相關(guān)關(guān)系。圖10中(c)表示X和Y沒(méi)有任何關(guān)系，(d)表示X和Y有非線性相關(guān)關(guān)系，r的計(jì)算公式為 (2.6)式中 (2.7)對(duì)例3 =764.861r=415.605/ =0.9894 后者很接近于1，故最大積雪深度與灌溉面積有很密切的線性相關(guān)關(guān)系，且是正相關(guān).但是，相關(guān)系數(shù)有一個(gè)缺點(diǎn)，就是它接近1的程度與樣本的組數(shù)n是有關(guān)的，當(dāng)n較小時(shí)，相關(guān)系數(shù)的絕對(duì)值容易接

6、近于1，當(dāng)n較大時(shí)，相關(guān)系數(shù)的絕對(duì)值容易偏小。特別當(dāng)n=2時(shí)，因?yàn)閮牲c(diǎn)決定一條直線，所以相關(guān)系數(shù)的絕對(duì)值總為1，在許多統(tǒng)計(jì)書中29給出相關(guān)系數(shù)的起碼值，當(dāng)相關(guān)系數(shù)的絕對(duì)值大于表中之值時(shí)才可以認(rèn)為X和Y有線性關(guān)系。此例當(dāng)顯著性水平=1%時(shí)，表中的起碼值為0.765，今計(jì)算r=0.9894 0.765，故最大積雪深度與灌溉面積有高度的線性關(guān)系。在有些統(tǒng)計(jì)軟件中，常給出，這時(shí)便于區(qū)別記為。 (b)方差分析和F檢驗(yàn) 因變量的波動(dòng)可用來(lái)表達(dá)，這個(gè)波動(dòng)是由兩個(gè)因素造成的；一個(gè)是X的變化引起Y相應(yīng)的變化，另一個(gè)是隨機(jī)誤差。前者造成Y的波動(dòng)可用回歸平方和來(lái)表達(dá)，后者用殘差平方和來(lái)度量。它們分別用 和來(lái)表示，從

7、數(shù)學(xué)上可以導(dǎo)出 - (2.8)當(dāng)X和Y為線性回歸模型(2.1)時(shí)，它們有如下更方便的計(jì)算公式- (2.9) 利用統(tǒng)計(jì)量F (2.10)可以來(lái)檢驗(yàn)回歸方程(2.1)是否可信.當(dāng)方程可信時(shí)F ，這里為F表中的臨界值，1和n-2為自由度，為顯著水平.對(duì)例3可以算得=1.802×415.606=748.922=764.961-748.922=16.039 F=8×748.922/16.039=373.55當(dāng)=1%時(shí)。用F值和F表上的臨界值相比，若F ，表明Y的變化主要是由X的變化造成的，回歸方程(2.1)可信；若F值小于，回歸方程不可信?？尚诺某潭纫部煞殖刹煌燃?jí)，在本書中，=5%

8、時(shí)可信用“*” 表示，=1%時(shí)可信用“*” 表示。上述計(jì)算結(jié)果常列成方差分析表，如表10所示。表10 方差分析表 方差來(lái)源平方和自由度均方F顯著性回歸748.9221748.922373.550*誤差16.03982.005總和764.9619 (c) 殘差分析稱為殘差，它能提供許多有用的信息，表11給出了例3的10個(gè)殘差，利用殘差可以提供如下信息： (i)之估計(jì) (2.11)給出了回歸方程的精度，它稱為殘差標(biāo)準(zhǔn)差，若隨機(jī)誤差遵從正態(tài)分布N(0，)，則Y的預(yù)報(bào)落在之內(nèi)的概率大約為95%，對(duì)例3可以算得=1.416，且10個(gè)均落于2×1.416之內(nèi)。 (ii)數(shù)據(jù)和模型之診斷 由殘差之

9、大小，可以發(fā)現(xiàn)異常(或叫離群)數(shù)據(jù)，可以發(fā)現(xiàn)模型(2.1)是否合適，是否要用非線性回歸模型等，這些已形成一整套理論，稱為回歸診斷，有興趣的讀者可參見文獻(xiàn)31。2.2多元線性回歸模型 當(dāng)影響因變量Y的自變量不止一個(gè)時(shí)，比如有m個(gè)，這時(shí)Y和X之間的線性回歸方程為 (2.12)其中為回歸系數(shù)，為隨機(jī)誤差，常假定 。 設(shè)為觀測(cè)值，回歸分析的首要任務(wù)是利用它們來(lái)估計(jì)和，它們的最小二乘估計(jì)記作求估計(jì)值 需要解下面的線性方程組 (2.13)其中 當(dāng) 求得后，計(jì)算 (2.14)回歸方程(2.12)建立后，檢驗(yàn)其是否可信可用方差分析，這時(shí)公式(2.8)依然有效，但方差分析表(參看表10)將成為表12之形式，其中

10、 (2.15) 表12 方差分析表方差來(lái)源平方和自由度均方F顯著性回歸m殘差n-m-1總和n-1 它將與F的臨界值來(lái)比較，其比較的結(jié)果和結(jié)論請(qǐng)參見上節(jié)的討論，反映回歸精度的的估計(jì)公式為 (2.16) 類似于一元回歸相產(chǎn)系數(shù)r，可以定義適用于多元回歸的全關(guān)系數(shù)R，R定義)為和的相關(guān)系數(shù)，或定義為 (2.17) 例4 試用線性回歸模型(2.10)來(lái)擬合表9的試驗(yàn)數(shù)據(jù)。 解：這時(shí)n=7，7組觀察值為(0.330，1.0，13，1.5)，(0.336，1.4，19，3.0)， (0.482，3.4，28，3.5)，它們的均值和為由于，故它們不必全部列出，將它們代入到方程級(jí)(2.13)中可以解得從而&#

11、215;2.2+0.00343 ×19-0.077×2.0 =0.201的估計(jì)為.于是回歸方程為 (2.18)進(jìn)一步對(duì)它作方差分析，其方差分析表列于表13.表13 方差分析表方差來(lái)源自由度平方和均方F回歸30.0487700.0162573.29誤差30.0148380.004946總和60.063608當(dāng) =0.05 時(shí)F表的臨界值，回歸方程(2.18)不可信.這時(shí)，是否Y和三個(gè)因素之間不可能建立回歸關(guān)系呢？不是的，我們還應(yīng)作進(jìn)一步探討，在下節(jié)我們將繼續(xù)討論該例。2.3 二次型回歸模型與變量篩選 由于因變量常常有交互作用，回歸模型(2.12)不足以反映實(shí)際，于是二次型回歸

12、模型常常為人們所采用.若有m個(gè)因素則二次型回歸模型為 (2.19)其中為回歸系數(shù)，為隨機(jī)誤差.我們看到，這時(shí)除了常數(shù)項(xiàng) 以外，方程有m(m+3)/2 項(xiàng)，當(dāng)m=1，2， 時(shí)項(xiàng)數(shù)為m12345678910項(xiàng)數(shù)25914202735445465若使回歸系數(shù)的估計(jì)有可能，必要條件為n>1+m(m+3)/2.當(dāng)m 較大時(shí)，通常不能滿足這個(gè)必要條件.于是有必要從方程(2.19)中選擇貢獻(xiàn)顯著的項(xiàng)，刪除不重要的項(xiàng)。有時(shí)，實(shí)際問(wèn)題需要考慮高階的交互作用，如 等，這時(shí)篩選變量的任務(wù)就更為重要。在回歸分析中，有許多有效的篩選變量的技術(shù)，如a) 前進(jìn)法，b) 后退法，c) 逐步回歸法，d) 最優(yōu)子集法(參看

13、25)。本章僅僅采用逐步回歸技術(shù)來(lái)篩選變量，這并不意味著逐步因歸在上述四項(xiàng)技術(shù)中最好的。 逐步回歸是回歸分析中的一種篩選變量的技術(shù)。開始它將貢獻(xiàn)最大的一個(gè)變量選入回歸方程，并且預(yù)先確定兩個(gè)閾值 和 ，用于決定變量能否入選或剔除。逐步回歸在每一步有三種可能的功能： a) 將一個(gè)新變量引進(jìn)回歸模型，這時(shí)相應(yīng)的F統(tǒng)計(jì)量必須大于 b) 將一個(gè)變量從回歸模型中剔除，這時(shí)相應(yīng)的F統(tǒng)計(jì)量必須小于 c) 將回歸模型內(nèi)的一個(gè)變量和回歸模型外的一個(gè)變量交換位置。 執(zhí)行功能a)和b)時(shí)要注意如下原則： 設(shè)在當(dāng)前步驟中有s個(gè)變量不在回歸模型中，有t個(gè)變量在回歸模型中。今欲從s個(gè)變量中挑選一個(gè)加入回歸模型之中，顯然應(yīng)挑

14、選使回歸效果最好的變量。這里回歸的效果可用方差分析表(見表10，表13)中F值來(lái)衡量，顯然我們要從s個(gè)變量中挑選一個(gè)變量使F值達(dá)到極大。類似地，若欲從t 個(gè)變量中刪除一個(gè)變量使其離開回歸模型，我們就是要選擇刪除后使回歸效果最好的變量，或選擇對(duì)當(dāng)前回歸模型貢獻(xiàn)最小的變量。如果在某一步中，既能實(shí)現(xiàn)a)又能實(shí)現(xiàn)b)，兩者之和就是功能c)。大部分統(tǒng)計(jì)軟件包均有逐步回歸之功能，例如中國(guó)均勻設(shè)計(jì)學(xué)會(huì)推薦的軟件包。 現(xiàn)在我們對(duì)例4繼續(xù)進(jìn)行討論。設(shè)先用后退法來(lái)選變量。所謂后退法，就是開始將所有的變量全部采用，然后逐步剔除對(duì)方程沒(méi)有顯著貢獻(xiàn)的變量，直到方程中所有的變量都有顯著貢獻(xiàn)為止。 仍考慮線性模型，開始三個(gè)

15、因素全部進(jìn)入方程，得(2.18)。統(tǒng)計(jì)軟件包通常還會(huì)提供每個(gè)變量的t值，t值越大(按絕對(duì)值計(jì))表示該因素越重要。對(duì)例2有這表明三個(gè)因素中以(反應(yīng)時(shí)間)對(duì)得率(Y)影響最大，配比次之，吡啶量最小。這些t 值都是隨機(jī)變量，它們遵從分布。若取=0.05 ，這時(shí)n=7，m=3， = 的臨界值(0.05)=3.18。t 值大于該值的因素表示對(duì)方程有顯著貢獻(xiàn)，否則表示不顯著。今 均小于(0.05)=3.18 ，說(shuō)明回歸方程(2.18)的三個(gè)變量至少有一個(gè)不起顯著作用。于是我們將貢獻(xiàn)最小的刪去，重新建立Y和及的線性回歸方程，得 (2.20)，三個(gè)t 值分別為這時(shí)這三個(gè)t值遵從含四個(gè)自由度的t 分布，臨界值為

16、(0.05)=2.78，從而 應(yīng)從方程中剔除。然后對(duì)Y和建立回歸方程 (2.21)相應(yīng)的。因此，回歸方程(2.21)為“最終”的回歸模型。這里最終加上引號(hào)，表示并非真正的最終模型，而是在線性模型框架下的最終產(chǎn)物。 上述的分析只發(fā)現(xiàn)對(duì)Y有顯著作用，其它兩個(gè)因素均沒(méi)有顯著作用，該結(jié)論與實(shí)際經(jīng)驗(yàn)不吻合，因此，猜想用線性模型不一定符合實(shí)際。于是進(jìn)一步考試二次回歸模型(2.19)。這時(shí)方程中有9項(xiàng)(不算)。利用逐步回歸技術(shù)求得回歸方程如下： (2.22)其相應(yīng)的 。顯然，回歸方程(2.22)的效果優(yōu)于回歸方程(2.21)。方程(2.22)表明，因素和交互作用對(duì)Y有顯著的影響。值得注意的是，有些人對(duì)回歸分

17、析沒(méi)有足夠的理解，片面追求大的(或小的)，致使選進(jìn)方程中的項(xiàng)過(guò)多，使誤差自由度為1或甚至為0，這時(shí)有關(guān)的結(jié)可靠性是很差的。因此，不應(yīng)片面追求大的，應(yīng)選擇n 稍大的均勻設(shè)計(jì)表，使得誤差有足夠的自由度5。2.4 應(yīng)用實(shí)例 均勻設(shè)計(jì)和正交設(shè)計(jì)以及其他試驗(yàn)設(shè)計(jì)方法一樣，在工農(nóng)業(yè)生產(chǎn)和科學(xué)實(shí)驗(yàn)中有廣闊的應(yīng)用前景，本文的文獻(xiàn)中列舉了部分應(yīng)用成果，其中有些成果成績(jī)顯著。更多的成果有待搜集。均勻設(shè)計(jì)不但在國(guó)內(nèi)得到廣泛應(yīng)用，在彼得格勒和香港也已開花結(jié)果。本節(jié)選擇香港浸會(huì)學(xué)院生物系的一項(xiàng)試驗(yàn)，供讀者參考。例5 為了研究環(huán)境污染對(duì)人體的危害，今考核六種金屬的含量：鎘(Cd)，銅(Cu)，鋅(Zn)，鎳(Ni)，鉻(

18、Cr)，鉛(Pb)，每種金屬含量分別取了17個(gè)水平(百萬(wàn)分之一，ppm)：0.01，0.05，0.1，0.2，0.4，0.8，1，2，4，5，8，10，12，14，16，18，20。今欲考慮這些金屬含量(包括它們的交互作用)對(duì)老鼠壽命的影響，該試驗(yàn)考核老鼠身上某種細(xì)胞的死亡率。它們選用表，根據(jù)使用表的指示，它們選用了表中1，4，6，10，14，15列來(lái)安排六個(gè)因素，其試驗(yàn)方案如表14所示。試驗(yàn)的結(jié)果為死亡率。為了了解試驗(yàn)誤差，提高結(jié)論的精度，他們?cè)谕辉囼?yàn)條件下將試驗(yàn)重復(fù)三次，三次結(jié)果()列于表15，三次死亡率的均值為，列于表15的最后一列。我們看到第17號(hào)試驗(yàn)的死亡率為最高，因?yàn)檫@時(shí)六種金屬

19、都是最高含量，表明這些金屬對(duì)老鼠細(xì)胞確有致命作用。 現(xiàn)進(jìn)一步用回歸分析來(lái)分析數(shù)據(jù)。由于每種金屬的含量由0.01(ppm)變到20(ppm)，最大得出小相差2000倍，于是直接用各因素的水平值作回歸不易獲得好的結(jié)果，通常要對(duì)水平值先作變換，用變換后的數(shù)據(jù)進(jìn)行回歸。最常見的變換是取對(duì)數(shù)。于是回歸分析中的自變量成為logCd，logCu，logZn，logNi，logCr和logPb。根據(jù)以往經(jīng)驗(yàn)，知道六種金屬間有交作用，故應(yīng)選用二次型回歸模型，并用逐步回歸來(lái)篩選變量。用同樣的 和 ，對(duì) 和分別進(jìn)行逐步回歸，發(fā)現(xiàn)四組數(shù)據(jù)的結(jié)果非常吻合，表明試驗(yàn)誤差不大，該試驗(yàn)可以獲得可靠結(jié)論。為節(jié)省篇幅，我們僅列出

20、對(duì)的回歸方程=32.68+5.03LogCd+3.84LogCu+2.03LogNi+0.55(LogCu) -0.63(LogZn)+0.94(LogNi)+0.53(LogCd)(LogCu) -0.70(LogCd)(LogCr)+0.92(LogCu)(LogPb)方程中每一項(xiàng)的t 值分別為(常數(shù)項(xiàng)t 值未列)11.5，7.8，4.9，2.6，-3.4，4.1，2.4，-2.8，5.3，它們均遵從 分布，因(0.05)=2.365 小于上述所有t值之絕對(duì)值，故方程可信。表15 死亡率17.9517.6518.3317.922.0922.8522.6222.531.7432.7932.8

21、732.439.3740.6537.8739.331.9031.1833.7532.231.1430.6631.1831.039.8139.6140.8040.042.4841.8643.7942.724.9724.6525.0524.850.2951.2250.5450.660.7160.4359.6960.267.0171.9967.1268.732.7730.8633.7032.429.9428.6830.6629.767.8769.2567.0468.055.5655.2856.5255.779.5779.4378.4879.1 由方程我們可以給出如下結(jié)論：a)Cd，Cu 和Ni含量過(guò)

22、高，對(duì)老鼠細(xì)胞的死亡率有顯著作用，b)金屬Cd和Cu，Cd和Cr，Cu和Pb有交互作用，其中Cd和Cu，Cu和Pb對(duì)死亡率起正交互作用，而Cd和Cu對(duì)死亡率起負(fù)交互作用，c)Zn可能會(huì)中和其它金屬的破壞作用，降低老鼠細(xì)胞的死亡率。2.5 尋求最優(yōu)工藝條件 試驗(yàn)設(shè)計(jì)的目的通常主要有二個(gè)，一是揭示變量(Y)與各因素之間的定性關(guān)系，二是尋求最優(yōu)工藝條件，回歸方程的建立可以達(dá)到一箭雙雕的目的。 現(xiàn)以例2來(lái)說(shuō)明如何尋求最好的工藝條件，表9告訴我們，第7號(hào)試驗(yàn)是7次試驗(yàn)中最好工藝條件，即配比3.4，吡啶量28，反應(yīng)時(shí)間3.5 ，這個(gè)工藝條件和最優(yōu)工藝條件常常是很接近的。 在上述討論中，我們最終建立了回歸模

23、型(2.22)。該方程一般僅在試驗(yàn)范圍內(nèi)成立，即配比，吡啶量10-28，反應(yīng)時(shí)間。尋求最優(yōu)模型等價(jià)于在這個(gè)范圍內(nèi)求方程(2.22)中的極大值。如果回歸方程比較復(fù)雜，可以用任何一個(gè)優(yōu)化算法(參見文獻(xiàn)33，34)來(lái)求最佳工藝條件，許多軟件包都含有優(yōu)化算法。數(shù)論方法也可以用來(lái)求的極大值，方開泰和王元提出了一個(gè)序貫算法SNTO，可以方便地求得的極大值，鑒于篇幅，這里就不詳細(xì)介紹了，有興趣的讀者可以參看文獻(xiàn)16。 對(duì)例2來(lái)講，可以用簡(jiǎn)單的微積分求得極值，由于X在試驗(yàn)范圍內(nèi)恒正，故由(2.22)知X 越大，越高，故X應(yīng)取試驗(yàn)范圍內(nèi)極大值3.4。將X=3.4 代入(2.22)得令，解得=0，=2.7575，這時(shí)的極大值為51.85%。工藝條件=3.4，=2.7575 并未出現(xiàn)在原有試驗(yàn)方案中，故應(yīng)在這個(gè)條件追加試驗(yàn)，由于的最佳條件在試驗(yàn)范圍邊界，故應(yīng)擴(kuò)大試驗(yàn)范圍。 對(duì)于許多實(shí)際工作者，不一定熟悉優(yōu)化方法，手邊沒(méi)有優(yōu)化的軟件。他們也不一定知道SNTO，也可能不會(huì)用微積分去求解極值。下面介紹一種“笨”辦法，其計(jì)算量較大，但程序好編。這種方法將每個(gè)因素的試驗(yàn)