第三章：概率論與數理統計模型與實驗3.2多元線性回歸分析3.2多元線性回歸分析

某種水泥在凝固時放出的熱量y（卡/克）與水泥中下列4種化學成份有關。

x1：3CaO·Al2O3的成份（%）

x2：3CaO·SiO2的成份（%）

x3：4CaO·Al2O3·Fe3O3的成份（%）

x4：2CaO·SiO2的成份（%）現記錄了13組數據，列在表3—4中，根據表中的數據，試研究y與x1，x2，x3，x4四種成份的關系。

3.2.1水泥凝固時放出熱量問題

表3—4表3—4編號x1(%)x2(%)x3(%)x4(%)y(卡/克)172666078.52129155274.331156820104.34113184787.6575263395.961155922109.27371176102.78131224472.59254182293.1102147426115.911140233483.8121166912113.3131068812109.43.2.2生豬生產決策問題

影響生豬生產的原因很多，為尋找影響生產的原因，現搜集整理了某市1980～1990年共11年的糧食產量、生豬外調量、鮮魚產量、家禽產量、豬肉銷售價格和生豬出欄量的統計數據（見表3—5），試研究生豬出欄量與糧食產量、生豬外調量、鮮魚產量、家禽產量、豬肉銷售價格之間的關系，并據此對未來的生豬生產進行預測。

表3—5

指標年份糧食產量（萬噸）x1豬外調量（千頭）x2鮮魚產量（百噸）x3家禽產量（萬只）x4豬肉售價(元/公斤)x5生豬出欄（萬頭）Y1980215.52190.94132.777701.60210.291981219.91270.42153.308341.607222.091982229.53230.72120.738171.627218.811983258.44250.58188.059062.014237.021984234.25305.76241.7211611.978263.191985211.57358.45255.3811592.262277.741986243.43426.07255.8512272.437272.421987239.79414.77304.6411842.827313.251988245.94345.53305.1711484.185321.911989253.17380.18353.9613584.547330.091990232.08333.60335.8812894.344335.92在現實生活中，變量與變量之間經常存在一定的關系，一般來說，變量之間的關系可以分為兩大類，一類是確定性的關系，這種關系通常用函數來表示。例如，已知圓的半徑r，那么圓的面積S與半徑r的關系就可用函數關系：來表示，這時如果取定了r的值，S的值就會完全確定了。另一類是非確定性關系，例如，人的體重與身高之間的關系就是非確定性關系，一般來說，身高越高，體重越大，但是身高相同的人體重往往是不相同的。再如，鋼材的強度與鋼材中含某種元素的含量，纖維的拉伸倍數與強度，降雨量、氣溫、施肥量與農作物的產量等均屬于這種關系。變量之間的這種非確定性關系通常稱為相關關系。

3.2.3多元線性回歸分析模型回歸分析就是數理統計中研究相關關系的一種數學方法，它就是通過大量的試驗或觀測，發現變量之間關系的統計規律。為了研究方便，我們考慮一個變量受其它變量影響時，仍把這變量稱為因變量，記為Y，其它變量稱為自變量，記為X，這時相關關系可記作Y=f（x）+ε

（1）其中f（x）為當X=x時，因變量Y的均值，即f（x）=E（Y|X=x）

稱f（x）為Y對X的回歸函數，ε為Y與f（x）的偏差，它是一個隨機變量，并假定E（ε）=0。“回歸”一詞始于英國統計學家蓋爾登（Galton）和皮爾遜（Pearson），他們研究兒子的身高Y與其父母的平均身高X之間的關系，通常雙親的平均身高較矮時，其子的身高也較矮；雙親的平均身高較高時，其子的身高也較高，但是他們發現，如果雙親屬于高個類（高于1000多對夫婦的平均身高）時，其子比他父親更高的概率會較小；相反地，如果雙親屬于低個類時，其子比他父親更矮的概率也較小，所以，平均身高偏高或偏矮的雙親，其子的身高都有一種“向中心回歸”的現象，基于這個事實，蓋爾登把求描述兒子與雙親身高關系的直線稱為回歸直線。

回歸函數可以是一元函數，也可以是多元函數，即Y=f（x1，x2，…，xm）+ε

（2）其中f（x1，x2，…，xm）=E（Y|X1=x1，X2=x2，…，Xm=xm）為m元回歸函數，統稱為多元回歸函數。若回歸函數f（x1，x2，…，xm）中的m=1,且是一元線性函數，則稱為是一元線性回歸；m>1且是多元線性函數，則稱為是多元線性回歸；若回歸函數f（x1，x2，…，xm）是非線性函數，則稱為是非線性回歸。對非線性回歸，經常采用線性化的方法來處理。所以，目前研究最多的是線性回歸，且假定X1，X2，…，Xm和Y均服從正態分布。回歸分析的任務就是要求出滿足式（2）的回歸函數f（x1，x2，…，xm），從而對所研究的相關關系作出所需的預測和控制。多元回歸模型的應用是相當廣泛的，例如，某種商品的銷售量可能受收入水平、風俗習慣、產品質量、價格、宣傳廣告等多種因素的影響；某種產品的質量可能受生產該產品時的溫度、濕度、壓力、原材料的質量和有害成份的含量等影響；工人的勞動生產率可能受學歷、智力水平、情緒的穩定性和才能等因素的影響；某城市的用水量可能與該城市的人口數及工業總產值有關。諸如此類的關系，可以通過多元回歸分析模型進行研究。例如，在水泥凝固時放出熱量問題中，可建立線性回歸模型其中E（ε）=0，D（ε）=σ2，b0，b1，，b2，b3，b4和σ2是未知參數，為了估計這些參數，將表3—4的值代入模型(3），其中，x1，x2，…，xm是自變量，b0為常數，b1，b2，…，bm為回歸系數，b0，b1，b2，…，bm皆為未知，統稱b0，b1，b2，…，bm為回歸參數，一旦回歸參數確定，則多元線性回歸模型就完全確定，一般假定隨機誤差ε～N（0，σ2）。得線性模型一般地，多元線性回歸模型可表示為為了得到回歸參數的估計值，就要對變量進行觀測，假設對變量的n（n>m）次獨立觀測數據為（yi，xi1，xi2，…，xim）,i=1～n,則這些觀測數據應滿足式（5），即有其中則多元線性回歸的數學模型式（6）可以寫成矩陣形式

若記其中。3.2.4參數的最小二乘估計

為了獲得參數β的估計，我們采用最小二乘法，即選擇β，使達到最小。將Q（β）對β求導數并令其為零，得

方程（9）稱為正規方程，其中X為n×（m+1）階矩陣，一般假定rank（X）=m+1，由線性代數理論可知，L=XTX為滿秩矩陣，它的秩rank（L）=m+1，則正規方程（9）有唯一解，記作

即記，則我們來證明（10）式中的為參數向量β的最小二乘法估計量，現用矩陣形式來敘述其證明步驟。從式（8）知，對任意的β則有上述證明過程中應用了如下結果：至此，在時，證明了式（10）中的是β的最小二乘法估計量。在實際工作中，常稱為經驗線性回歸方程。3.2.5最小二乘法估計量的性質

首先我們在假定的條件下，探討一下由式（10）確定β的最小二乘法估計量的性質。

1、是β的線性無偏估計量。證：由于，每一個都是y1，…，yn的線性組合，因而是的線性估計量，此時稱的線性估計量。

證：記

=E{[BY—BE(Y)][BY—BE(Y)]T}=B·E{[Y—E(Y)][Y—E(Y)]T}BT=B·σ2In·BT

=L-1XT·σ2In·（L-1XT）T

=σ2L-1即2、的協方差矩陣為，即其中性質2告訴我們，用最小二乘法求出的諸回歸系數之間存在相關性，進一步可以證明.3、是β的最小方差線性無偏估計，即在所有線性無偏估計類中，有且只有使其方差達到最小。

3.2.6多元線性回歸方程的顯著性檢驗

從上面的參數估計過程可以看出，對于一批觀察數據不論它們是否具有線性關系，總可以利用最小二乘法建立起多元線性回歸方程但是Y與x1，x2，…，xm

是否確實存在相關關系呢？回歸方程的效果如何呢？這就要進行“整個回歸效果是否顯著”的檢驗。當時，沒有關系，回歸模型沒有意義，于是我們要檢驗是否成立。若H0成立，則x1，x2，…，xm對y沒有影響；反之，若H0不成立，則x1，x2，…，xm對y有影響，此時y與x1，x2，…，xm的線性關系顯著，也稱為整個回歸效果顯著。但要注意，即使整個回歸效果是顯著的，y也可能只與某幾個xi關系密切（相應的bi顯著不為零），而與另幾個xi關系不密切（相應的bi為零）。這就是說，多元線性回歸除了首先要檢驗“整個回歸是否顯著”外，還要逐個檢驗每一個bi是否為零，以便分辨出哪些xi對y并無顯著影響，最后，還要對各個bi作出區間估計。為了進行檢驗和區間估計，可以證明以下結論成立：則稱lyy為總變差或稱為y的離差平方和。lyy可進行如下分解：這時稱為殘差平方和。稱為回歸平方和。記稱為剩余標準差或估計的標準差。由于lyy不變，當然希望Q越小越好，即U越大越好，因此，定義復相關系數（complexcorrelationcoefficient）

當觀察值全都與回歸值

吻合時，Q=0，R=1；當時，，，在一般情況下，R的數值在0與1之間。復相關系數R的定義，類似于兩個變量時的相關系數的定義，但要注意，復相關系數R只取正值。在兩個變量時，有正相關與負相關之分，在多個變量時，就沒有這一說了，所以復相關系數R只取正值。（2）在的條件下且U與Q獨立，因此（3）

i=1，2，…，m

這里cii為C=L-1中第i個對角線元素。利用上述幾條結論，可進行下列檢驗、估計和預測。（1）回歸顯著性檢驗（F檢驗）即整個回歸效果是否顯著性檢驗。若整個回歸效果不顯著，即全部回歸系數為零，因此，設原假設H0：b1=b2=…=bm=0。若H0為真，則而且在H0不成立時，F值有變大的趨勢，因此應取右側否定域，故檢驗法是當時拒絕原假設，認為回歸效果顯著；否則認為回歸效果不顯著。（2）單個回歸系數為零的檢驗（t檢驗），即某個自變量是否對因變量有顯著性影響的檢驗。在多元回歸分析中可能出現y與所有自變量的總體是有相關關系的，但y與某個特定的xi則可能無關，即xi對y并不起作用或者已被其它的xi的作用所代替，為此設m個原假設若H0i為真，統計量而當H0i不成立時，有變大的趨勢，因而應取雙側拒絕域，故當時，否定H0i，即認為xi對y是有作用的，若某幾個xi是有作用的，而另幾個xi是不起作用的，則應從回歸方程中刪除那些不起作用的自變量。單個回歸系數是否為零，也可以用F檢驗，即若H0i為真，統計量故當時，拒絕原假設，即認為xi對y的影響是顯著的；否則認為xi對y的影響是不顯著的。（4）y0的95%預測區間近似為，其中（3）對bi的區間估計

由于因而bi的置信區間為

其中3.2.7多元線性回歸分析模型的推廣

（1）多項式回歸分析模型類似于模型（5），由自變量多項式和隨機項組成的回歸模型稱為多項式回歸模型，它的一般形式為初看模型（13）不是線性回歸，因自變量中含有冪函數，但由于未知參數b0，b1，…，bm都是線性出現的，因此，令則模型（13）就變成為多元線性回歸模型：從而多項式回歸模型可以用多元線性回歸模型的計算公式和檢驗方法。多項式回歸還有許多推廣的形式，

例如：①②③④⑤

上述模型的共同特點是未知參數都是以線性形式出現，所以都可以采用恒等變換像模型（13）化為模型（14）一樣化為多元線性回歸模型。（2）廣義線性回歸模型廣義線性回歸模型的一般形式為

其中：是一個不含未知數參數的一元函數，且有反函數：

的不含未知參數的多元函數。

廣義線性回歸模型的回歸系數的確定主要是從自變量和因變量y的n組觀察值出發，用最小二乘法求出，使得達到最小。此時也就是令則這樣就把廣義線性回歸模型化為多元線性回歸模型。例如：這里為已知的常數，就屬于廣義線性回歸模型。3.2.8最優回歸模型的選擇

很多實際問題中，例如氣象預報、地震預報等，自變量個數多得驚人，過多的自變量不僅給數據的收集、計算帶來很大麻煩，而且多個自變量之間往往會產生線性相關或近似線性相關關系。這種相關關系，稱為“復共線性”，復共線性的存在，會使得回歸分析的誤差增大，得到不合理的結果，甚至使計算溢出，無法求得回歸分析問題的解。所以，我們總是希望從眾多自變量中挑選較少的一些變量，建立“最優”的回歸模型。所謂“最優”的回歸模型有兩個方面的含義：一方面指方程中要包含所有的有顯著作用的自變量；另一方面，希望自變量的個數盡可能地少，因為自變量過多，使用起來不方便，而且引入了無意義的量，反而會使回歸方程的穩定性降低。為了解決這些問題，就產生了一種重要的新的方法，即逐步回歸法，它是一種能避免由復共線性產生的不良影響的回歸分析方法，其基本思想是：不取全部自變量作回歸，而是從一個含自變量的方程出發，通過逐步引入和刪除，從全部自變量中選取一部分自變量作回歸，既要使殘差平方和盡量小，又要避免引入過多的自變量，以致產生復共線性現象。

來表示它，得回歸方程及回歸平方和U，再用回歸模型

來表示它，得回歸方程

及回歸平方和Ui，則稱為偏回歸平方和。下面簡單介紹一下數學原理。若xi可用其余自變量線性表示，則U與Ui應當差不多，因而Si應當很小。所以當從F檢驗中發現Si很小時，我們就可以從P個自變量中刪除它。類似地，我們也可以檢驗xp+1，…，xm中的自變量是否應當添到回歸模型中去。逐步回歸的具體做法是：事先給定一個引入變量的F水平界限Fin和一個刪除變量的F水平界限Fout，首先在所有尚未引進方程的自變量中，找出Fi值最大的一個，如果Fi>Fin就引入自變量xi，然后在所有已引進方程的自變量中，找出Fi值最小的一個，如果Fi≤Fout，就刪除自變量xi。這樣逐步引入、刪除，直到方程外沒有一個自變量的Fi大于Fin，方程內沒有一個自變量的Fi小于Fout為止。在逐步回歸中，檢驗作用顯著的臨界值F像一個門坎，它的值取大些，進入回歸方程的變量就會少些，反之就會多些。3.2.9MATLAB和Mathematica求解

（1）MATLAB命令命令格式[b,bint,r,rint,stats]=regress(Y,X,alpha)其中輸入向量X，Y的排列方式分別為

alpha為顯著性水平（缺省時設定為0.05）。輸出向量b為回歸系數的估計值，即輸出向量bint為回歸系數估計值的置信區間；

輸出向量r為殘差向量；輸出向量rint為殘差向量的置信區間；輸出向量，它是一個3維向量，用于檢驗回歸模型的統計量，其中第一個分量R2中的R是相關系數，第二個分量是F統計量，第三個分量是與統計量F對應的概率P，當P<alpha時拒絕原假設H0，說明回歸模型成立。（2）Mathematica命令在Mathematica中健入命令<<Statisti\Linearre.m按Shift+Enter健，即可調入線性回歸軟件包。輸入（3）實際問題的求解①水泥凝固時放出熱量問題在MATLAB編輯器中輸入以下程序：%水泥放出熱量問題ch321%文件名：ch321.mx1=[7111117113122111110]';x2=[26295631525571315447406668]';x3=[615886917221842398]';x4=[6052204733226442226341212]';y=[78.574.3104.387.695.9109.2102.772.5...93.1115.983.8113.3109.4]';x=[ones(13,1)x1x2x3x4];[b,bint,r,rint,stats]=regress(y,x,0.05);disp('回歸系數估計值')bdisp('回歸系數估計值的置信區間')bintdisp('殘差平方和')r'*rdisp('相關系數的平方')stats(1)disp('F

統計量')stats(2)disp('與統計量F對應的概率p')stats(3)

執行后輸出回歸系數估計值b=62.40541.55110.51020.1019-0.1441回歸系數估計值的置信區間bint=-99.1786223.9893-0.16633.2685-1.15892.1792-1.63851.8423-1.77911.4910殘差平方和ans=47.8636相關系數的平方ans=0.9824F統計量ans=111.4792與統計量F對應的概率pans=4.7562e-007

從計算結果可知，回歸方程查表得：，易見統計量，進一步可得，所以回歸效果是極高度顯著的。

②生豬生產決策問題1°數學模型首先選取生豬出欄量作因變量Y，用直接影響生豬生產的糧食產量x1，豬外調量x2，鮮魚產量x3，家禽產量x4，豬肉價格x5作為自變量。

假設上一年度各自變量的變化影響下一年度的

t——年份，i=1，2，3，4，5

對于某些變量，這一點是顯然的。例如x1（糧食產量），上一年的豐歉程度，對下一年的糧食消耗至關重要，若其對生豬飼養有影響，自然應是上一年而不是本年度的糧食產量在起主要作用，其它因素也是如此。為尋找出影響生豬出欄量的關鍵變量，設計多元線性逐步回歸模型來擬合生豬出欄量的關鍵變量，設計多元線性逐步回歸模型來擬合生豬出欄量與各自變量xi的數量關系，

其模型形式如下：

式中：Yt——第t年的生豬出欄量；

——第（t-1）年的糧食產量；

——第（t-1）年的生豬外調量；

——第（t-1）年的鮮魚產量；

——第（t-1）年的家禽產量；

——第（t-1）年的豬肉售價；

bi——待估計參數，i=1，2，3，4，5。2°應用逐步回歸求解模型（a）先考慮一個自變量人們普遍認為糧食產量x1是生豬生產的氣候表，故首先引入x1，此時（18）變為按最小二乘法，求得回歸方程這表明，對于生豬出欄量來說，糧食產量只解釋了總變差40%。以上結果并不象人們一般地認識那樣，生豬生產最主要地取決于糧食生產狀況，嚴格地說，可以認為當糧食減少到某一臨界范圍時，即只能滿足人們生活的基本必需時，才有可能對生豬生產產生決定性的影響。（b）為簡便起見，同時引入其余5個變量，（18）變為

下面給出模型（19）的MATLAB計算程序。%生豬生產決策問題ch322%文件名：ch322.mx1=[215.52219.91229.53258.44243.25211.57243.43239.79245.94253.17232.08]';x2=[190.94270.42230.72250.58305.76358.45426.07414.77345.53380.18333.60]';x3=[132.77153.30120.73188.05241.72255.38255.85304.64305.17353.96335.88]';x4=[7708348179061161115912271184114813581289]';x5=[1.601.6071.6272.0141.9782.2622.4372.8274.1854.5474.344]';y=[210.29222.09218.81237.02263.19277.74272.42313.25321.91330.09335.92]';x=[ones(11,1)x1x2x3x4x5];[b,bint,r,rint,stats]=regress(y,x,0.05);disp('回歸系數估計值')bdisp('回歸系數估計值的置信區間')bintdisp('殘差平方和')r'*rdisp('相關系數的平方')stats(1)disp('F

統計量')stats(2)disp('與統計量

F對應的概率

p')stats(3)執行后輸出回歸系數估計值b=192.1008-0.16240.05680.5233-0.048810.3357回歸系數估計值的置信區間

bint=67.7756316.4260-0.65810.3333-0.14620.25980.03781.0087-0.18630.0887-7.413928.0854

殘差平方和ans=350.5414相關系數的平方ans=0.9841F統計量ans=62.0595與統計量

F對應的概率

pans=1.6912e-004從計算結果易見，回歸方程