Matlab實現(xiàn)多元的回歸實例

1、Matlab 實現(xiàn)多元回歸實例(一)一般多元回歸一般在生產(chǎn)實踐和科學(xué)研究中，人們得到了參數(shù)x Xi, ,Xn和因變量y的 數(shù)據(jù)，需要求出關(guān)系式y(tǒng) f X，這時就可以用到回歸分析的方法。 如果只考慮 f是線性函數(shù)的情形，當(dāng)自變量只有一個時，即， X Xi, ,Xn中n 1時，稱 為一元線性回歸，當(dāng)自變量有多個時，即， X Xi, ,Xn中n 2時，稱為多元 線性回歸。進(jìn)行線性回歸時，有4個基本假定： 因變量與自變量之間存在線性關(guān)系； 殘差是獨立的； 殘差滿足方差奇性； 殘差滿足正態(tài)分布。在Matlab軟件包中有一個做一般多元回歸分析的命令regeress，調(diào)用格式如下：此時，默認(rèn)alpha =

2、0.05.m 1的矩陣，其中第一列是全1向1列向量對應(yīng)回歸方程的常數(shù)項)， 回歸方程具有如下形式：mx mb, bint, r, rint, stats = regress(y,X,alpha)或者b, bint, r, rint, stats = regress(y,X)這里，y是一個n 1的列向量，X是一個n 量(這一點對于回歸來說很重要，這一個全 一般情況下，需要人工造一個全1列向量y 01X1其中，是殘差。在返回項b,bint,r,rint,stats中, b o im是回歸方程的系數(shù)； bint是一個m 2矩陣，它的第i行表示i的(1-alpha)置信區(qū)間； r是n 1的殘差列向量；

3、 rint是n 2矩陣，它的第i行表示第i個殘差的(1-alpha)置信區(qū)間；注釋：殘差與殘差區(qū)間杠桿圖，最好在0點線附近比較均勻的分布，而不呈現(xiàn)一定的規(guī)律性，如果是這樣，就說明回歸分析做得比較理想。 一般的，stast返回4個值：R2值、F_檢驗值、閾值f，與顯著性概率相關(guān) 的p值(如果這個p值不存在，貝U，只輸出前3項)。注釋：(1) 一般說來，R2值越大越好。(2) 人們一般用以下統(tǒng)計量對回歸方程做顯著性檢驗： F_檢驗、t_檢驗、以及 相關(guān)系數(shù)檢驗法。Matlab軟件包輸出FJ檢驗值和閾值f。一般說來，F(xiàn)_檢驗值 越大越好，特別的，應(yīng)該有F_檢驗值f。(3) 與顯著性概率相關(guān)的p值應(yīng)該

4、滿足p alpha。如果p alpha，則說明回歸 方程中有多余的自變量，可以將這些多余的自變量從回歸方程中剔除 (見下面逐 步回歸的內(nèi)容)。這幾個技術(shù)指標(biāo)說明擬合程度的好壞。這幾個指標(biāo)都好，就說明回歸方程是有意義的。例 1 (Hamilton ， 1987 )數(shù)據(jù)如下:序號YX1X2112.372.239.66212.662.578.94312.003.874.40411.933.106.64511.063.394.91613.032.838.52713.133.028.04811.442.149.05912.863.047.711010.843.265.111111.203.395.051

5、211.562.358.511310.832.766.591412.633.904.901512.463.166.96第一步分析數(shù)據(jù)在Matlab軟件包中分析是否具有線性關(guān)系，并作圖觀察，M 文件opt_ha nmilton_1987 x1=2.23,2.57,3.87,3.10,3.39,2.83,3.02,2.14,3.04,3.26,3.39,2.35,2.76,3.90,3.16;x2=9.66,8.94,4.40,6.64,4.91,8.52,8.04,9.05,7.71,5.11,5.05,8.51,6.59,4.90,6.96; y=12.37,12.66,12.00,11.93

6、,11.06,13.03,13.13,11.44,12.86,10.84,11.20,11.56,10.83,12.63, 12.46;corrcoef(x1,y);corrcoef(x2,y);plot3(x1,x2,y, '*');得到結(jié)果：ans =1.00000.00250.00251.0000ans =1.00000.43410.43411.0000即，corrcoef(x1,y) = 0.0025 , corrcoef(x2,y) = 0.4341，說明沒有非常明顯的單變量線性關(guān)系。圖形如下：也看不出有線性關(guān)系，但是，旋轉(zhuǎn)圖形，可以看出所有點幾乎在一個平面上這說明，

7、y, x1,x2在一個平面上，滿足線性關(guān)系：ai xi a2 X2 by a或者，換成一個常見的形式y(tǒng) a。ai Xi a? x?其中， 是殘差。于是，在Matlab軟件包中做線性多元回歸，寫一個 M 文件opt_regress_hamilt onx1=2.23,2.57,3.87,3.10,3.39,2.83,3.02,2.14,3.04,3.26,3.39,2.35,2.76,3.90,3.16： x2=9.66,8.94,4.40,6.64,4.91,8.52,8.04,9.05,7.71,5.11,5.05,8.51,6.59,4.90,6.96：y=12.37,12.66,12.00

8、,11.93,11.06,13.03,13.13,11.44,12.86,10.84,11.20,11.56,10.83,12.63, 12.46'e=on es(15,1);x=e,x1,x2;b,b in t,r,ri nt,stats=regress(y,x,0.05)rcoplot(r,ri nt)其中，rcoplot ( Residual case order plot )表示畫出殘差與殘差區(qū)間的杠桿圖。執(zhí)行后得到：b =-4.51543.0970-4.38223.12381.03991.0319 bint =-4.64863.07031.0238 r =0.0113-0.0

9、087-0.0102-0.00690.0101-0.0106-0.0037-0.01050.0049-0.01360.00570.0163-0.00230.01100.0071rint-0.00870.0314-0.03030.0128-0.03010.0098-0.02990.0162-0.01060.0308-0.03130.0102-0.02520.0178-0.02990.0089-0.01740.0272-0.03310.0058-0.01610.0275-0.00270.0354-0.02360.0190-0.00790.0299-0.01560.0298stats =1.0e+0

10、040.00013.922200.0000即，y 4.5153.097x11.0319x2。置信度95%,且R2 1.0, F _檢驗值39222 0 ,與顯著性概率 0.05相關(guān)的p 0.0000 0.05，這說明，回歸方程中的每個自變量的選取，都是有意義的殘差杠桿圖:3n0% !l從杠桿圖看出，所有的殘差都在0點附近均勻分布，區(qū)間幾乎都位于0.03,0.03之間，即，沒有發(fā)現(xiàn)高杠桿點，也就是說，數(shù)據(jù)中沒有強(qiáng)影響點、異常觀測點。綜合起來看，以上回歸結(jié)果(回歸函數(shù)、擬合曲線或曲面)近乎完美。(二)逐步回歸假設(shè)已有數(shù)據(jù)X和丫，在Matlab軟件包中，使用stepwise命令進(jìn)行逐步 回歸，得到回

11、歸方程丫 aiXi 32X2anXn ，其中 是隨機(jī)誤差。stepwise命令的使用格式如下：stepwise(X,Y)注意：應(yīng)用stepwise命令做逐步回歸，數(shù)據(jù)矩陣X的第一列不需要人工加 一個全1向量，程序會自動求出回歸方程的常數(shù)項(in tercept )。在應(yīng)用stepwise命令進(jìn)行運算時，程序不斷提醒將某個變量加入(Move in )回歸方程，或者提醒將某個變量從回歸方程中剔除(Move out )。注釋：使用stepwise命令進(jìn)行逐步回歸，既有剔除變量的運算，也有引 入變量的運算，它 是目前應(yīng)用較為廣泛的一種多元回歸方法。在 運行 stepwise(X,Y)命令時，默認(rèn)顯著性

12、概率0.05。例2 (Hald,1960 ) Hald數(shù)據(jù)是關(guān)于水泥生產(chǎn)的數(shù)據(jù)。某種水泥在凝固時放出的熱量Y (單位：卡/克)與水泥中4種化學(xué)成分所占的百分比有關(guān)：Xi : 3Cao AI2O3X2 : 3Cao SiO2x3 :4Cao Al2o3 Fe2o3X4 :2Cao SiQ在生產(chǎn)中測得13組數(shù)據(jù):序號X1X2X3X4Y172666078.52129155274.331156820104.34113184787.6575263395.961155922109.27371176102.78131224472.59254182293.1102147426115.911140233483.

13、8121166912113.3131068812109.4求出關(guān)系式Y(jié) f X解：(1)本問題涉及的數(shù)據(jù)是5維的，不能畫圖觀察。先做異常值分析。X=7,26,6,60;1,29,15,52;11,56,8,20;11,31,8,47;7,52,6,33;11,55,9,22;3,71,17,6;1,31,22,44;2,54,18,22;21,47,4,26;1,40,23,34;11,66,9,12;10,68,8,12;Y=78.5,74.3,104.3,87.6,95.9,109.2,102.7,72.5,93.1,115.9,83.8,113.3,109.4A=X,Y;mahal(A,

14、A)程序執(zhí)行后得到結(jié)果：ans =5.68033.64846.70023.36763.38394.43004.00806.50673.08497.50165.17682.4701可以認(rèn)為數(shù)據(jù)都是正常的(2) 一般多元回歸。在Matlab 軟件包中寫一個 M 文件opt_cement_1:X=7,26,6,60;1,29,15,52;11,56,8,20;11,31,8,47;7,52,6,33;11,55,9,22;3,71,17,6;1,31,22,44;2,54,18,22;21,47,4,26;1,40,23,34;11,66,9,12;10,68,8,12;Y=78.5,74.3,10

15、4.3,87.6,95.9,109.2,102.7,72.5,93.1,115.9,83.8,113.3,109.4'a1=o nes(13,1);A=a1,X;b,b in t,r,ri nt,stat=regress(Y,A)rcoplot(r,ri nt)程序執(zhí)行后得到：b =62.40541.55110.51020.1019-0.1441bint =-99.1786223.9893-0.16633.2685-1.1589-1.6385-1.7791 r =0.00481.5112-1.6709-1.72710.25083.9254-1.4487-3.17501.37830.28

16、151.99100.9730-2.2943 rint =-4.0390-3.2331-5.31262.17921.84231.49104.04856.25551.97073.1061-6.5603-4.5773-0.5623-6.0767-6.8963-3.5426-3.0098-2.2372-4.1338-6.9115stat =5.07888.41323.17940.54636.29933.57296.21916.07972.3228111.47920.00005.9830以及殘差杠桿圖:于是，我們得到:0.9824Y 62 .40541 .5511 x1 0.5102 x2 0.1019

17、 x3 0.1441 x4并且，殘差杠桿圖顯示，殘差均勻分布在 0點線附近，在stat返回的4個值中,R2 = 0.9824，說明模型擬合的很好。F檢驗值=111.4792>0.000 ，符合要求。 但是，與顯著性概率相關(guān)的p值=5.9830>0.05，這說明，回歸方程中有些變量 可以剔除。(3)逐步回歸在Matlab 軟件包中寫一個 M 文件opt_cement_2:X=7,26,6,60;1,29,15,52;11,56,8,20;11,31,8,47;7,52,6,33;11,55,9,22;3,71,17,6;1,31,22,44;2,54,18,22;21,47,4,26

18、;1,40,23,34;11,66,9,12;10,68,8,12;Y=78.5,74.3,104.3,87.6,95.9,109.2,102.7,72.5,93.1,115.9,83.8,113.3,109.4'stepwise(X,Y)程序執(zhí)行后得到下列逐步回歸的畫面：X1X2X3X41.868753.55000.0046.90.7891254.68620.0007-1.25578-2.09840.0598-0.738162-4.77480.0006» >1.Nexl stepMove X4 InExport.Coeff. t-stat p-valCoefficie

19、nts with Error Bars-3-2-10123Intercepts 95.4231RMSEs 15.0437R-squre - 0F -丹赳NAdj R-sq - -0 .D833333p-NeriN程序提示：將變量x4加進(jìn)回歸方程(Move x4 in )，點擊Next Step按鈕，即，進(jìn)行下一步運算，將第4列數(shù)據(jù)對應(yīng)的變量X4加入回歸方程。點擊Next Step按健后，又得到提示：將變量x1加進(jìn)回歸方程(Move x1 in )，點擊Next Step按鈕，即，進(jìn)行下一步運算，將第1列數(shù)據(jù)對應(yīng)的變量加入回歸方程。點擊Next Step按健后，又得到提示：Move No terms ，即，沒有需要加入(也沒 有需要剔除)的變量了。注意：在Matlab7.0軟件包中，可以直接點擊“ All Steps ”按鈕，直接求 出結(jié)果(省略中間過程)