李艷嬌，李瑞敏，陳經(jīng)偉

（1.沈陽建筑大學(xué) 建筑設(shè)計(jì)研究院，遼寧 沈陽 110015；2.沈陽機(jī)床（集團(tuán)）有限責(zé)任公司 設(shè)計(jì)研究院，遼寧 沈陽 110142）

在社會(huì)生活及生產(chǎn)實(shí)踐中經(jīng)常會(huì)遇到這樣一種問題，即我們非常關(guān)注一個(gè)量的變化，而這個(gè)量受到另一個(gè)或是多個(gè)因素的影響.我們想要了解這些因素是如何影響我們最為關(guān)注的這個(gè)量的以及這些因素對我們最為關(guān)注的這個(gè)量的影響權(quán)重分別有多大.知道了這些，我們就可以對該量變化所反應(yīng)的相關(guān)問題做出分析和評價(jià)，并對其未來發(fā)展趨勢進(jìn)行預(yù)測和控制.這里就要用到數(shù)理統(tǒng)計(jì)中一個(gè)非常重要而普遍的分析方法，即回歸分析法.

MATLAB是一個(gè)非常強(qiáng)大的通用數(shù)學(xué)工具，利用其可以輕而易舉地解決各學(xué)科領(lǐng)域中的一般性數(shù)學(xué)問題.下面我們首先簡單介紹一下回歸分析法，在此基礎(chǔ)上結(jié)合一具體實(shí)例，詳細(xì)介紹利用MATLAB進(jìn)行多元線性回歸分析的幾種方法.

1 回歸分析及多元線性回歸模型

回歸分析是研究一個(gè)因變量與另一個(gè)或一組自變量間定量關(guān)系的一種統(tǒng)計(jì)分析方法.只有一個(gè)自變量的回歸分析叫做一元回歸分析，有多個(gè)自變量的回歸分析叫做多元回歸分析.按照自變量及因變量間的關(guān)系類型，回歸分析又有線性回歸及非線性回歸.所以就有一元線性回歸、多元線性回歸、一元非線性回歸和多元非線性回歸四種基本回歸分析類型[1-2].

如果一個(gè)因變量y與k個(gè)自變量x1,x2,...,xk存在線性相關(guān)關(guān)系，那么就可以用多元線性回歸模型

對其進(jìn)行描述，其中未知常量a0,a1,...,ak稱為回歸模型系數(shù).若n次抽樣，第i次抽樣數(shù)據(jù)為(yi,x1i,x2i,…,xki)，那么就有

其中 ε0,ε1,...,εn為隨機(jī)誤差項(xiàng)[1-2].回歸分析的主要任務(wù)就是以誤差ε0,ε1,...,εn的平方和最小為原則，求取多元線性回歸模型的回歸系數(shù)a0,a1,...,ak.

某單位生產(chǎn)部門有一具體問題，生產(chǎn)人員已明確他們關(guān)心的變量y與兩個(gè)因素x1,x2具有線性相關(guān)關(guān)系，根據(jù)式（1）知可以用二元線性回歸模型對其進(jìn)行描述.車間給出的n=16個(gè)測試數(shù)據(jù)見表1.

現(xiàn)在需要技術(shù)部門幫助定量地確定兩個(gè)因素x1,x2對他們關(guān)心的變量y的影響權(quán)重并對變量y的變化趨勢作出預(yù)判以便控制.

2 基于回歸算法的MATLAB編程實(shí)現(xiàn)

根據(jù)式（2）及表1中的16個(gè)測試數(shù)據(jù)有

上式可以寫成形式

其中：

式（5）是適定方程組，在MATLAB中用 A=inv(X)*Y即可求解，其中inv(X)是X的逆陣.也可以用左除函數(shù)“”求解，即A=XY或 A=mldiv ide(X,Y).對于適定方程組來說，XY與inv(X)*Y具有相同的意義[3-4].

基于回歸算法求解上述問題的MATLAB程序如下：

x1=[0，0，0，0，0.1，0.1，0.1，0.1，0.2，0.2，0.2，0.2，0.3，0.3，0.3，0.3]；

x2=[0，0.1，0.2，0.3，0，0.1，0.2，0.3，0，0.1，0.2，0.3，0，0.1，0.2，0.3]；

y=[0，2*10^-6，10^-6，-2*10^-6，10^-6，2*10^-6，4*10^-6，-10^-6，-2*10^-6，-10^-6，0，3*10^-6，4*10^-6，-4*10^-6，2*10^-6，5*10^-6]；

X（3，3）=0；Y（3，1）=0；

for i=1：16

X（1，1）=X（1，1）+1；

X（1，2）=X（1，2）+x1（i）；

X（1，3）=X（1，3）+x2（i）；

X（2，1）=X（2，1）+x1（i）；

X（2，2）=X（2，2）+x1（i）^2；

X（2，3）=X（2，3）+x1（i）*x2（i）；

X（3，1）=X（3，1）+x2（i）；

X（3，2）=X（3，2）+x1（i）*x2（i）；

X（3，3）=X（3，3）+x2（i）^2；

Y（1，1）=Y（1，1）+y（i）；

Y（2，1）=Y（2，1）+x1（i）*y（i）；

Y（3，1）=Y（3，1）+x2（i）*y（i）；

end

A=inv（X）*Y%同A=XY

解得結(jié)果為a0=-1.0×10-7,a1=3×10-6,a2=3.5×10-6.

3 利用MATLAB預(yù)定義函數(shù)實(shí)現(xiàn)

式（4）可以寫成形式

其中：

式（6）是一個(gè)超定方程組，可以用 A=regress(Y,X)來求解.在MATLAB中，regress函數(shù)返回的是Y=XA的最小二乘擬合解[5-6].式（6）方程組也可以用左除函數(shù)“”來求解，即 A=XY或 A=mldiv ide(X,Y).對于超定方程組來說，左除函數(shù)“”與regress函數(shù)相同，返回的也是方程Y=XA的最小二乘擬合解.

利用預(yù)定義函數(shù)regress求解上述問題的MATLAB程序如下：

x1=[0，0，0，0，0.1，0.1，0.1，0.1，0.2，0.2，0.2，0.2，0.3，0.3，0.3，0.3]'；

x2=[0，0.1，0.2，0.3，0，0.1，0.2，0.3，0，0.1，0.2，0.3，0，0.1，0.2，0.3]'；

Y=[0，2*10^-6，10^-6，-2*10^-6，10^-6，2*10^-6，4*10^-6，-10^-6，-2*10^-6，-10^-6，0，3*10^-6，4*10^-6，-4*10^-6，2*10^-6，5*10^-6]'；

X=[ones（16，1），x1，x2]；

A=regress（Y，X）%同A=XY

解得結(jié)果為a0=-1.0×10-7,a1=3×10-6,a2=3.5×10-6.

與上節(jié)基于回歸算法編程所解得的結(jié)果相同.將其代入式（3）知生產(chǎn)人員關(guān)心的變量y與兩個(gè)因素x1,x2之間的關(guān)系為y=-1.0×10-7+3×10-6x1+3.5×10-6x2.

4 結(jié)束語

（1）多元線性回歸問題中一般會(huì)有大量數(shù)據(jù)，它的解算涉及到數(shù)理統(tǒng)計(jì)、線性代數(shù)、矩陣分析及微積分等諸多高等數(shù)學(xué)知識(shí).MATLAB作為一個(gè)強(qiáng)大的通用數(shù)學(xué)工具，可以在沒有專用軟件的情況下，快速方便地解決諸如多元線性回歸等一般性數(shù)學(xué)問題.

（2）在MATLAB中，有利用預(yù)定義函數(shù)及自行編程實(shí)現(xiàn)算法兩種解決數(shù)學(xué)問題的基本方法.其中利用預(yù)定義函數(shù)更為簡單快捷，不需要了解具體算法，但不一定所有問題都可以找到最適合的預(yù)定義函數(shù)；自行編寫M文件需要理解原理、熟悉算法，較復(fù)雜，但可能更靈活，更適合具體問題的解決.所以實(shí)際應(yīng)用中應(yīng)采用二者相結(jié)合的方法以快速有效地解決問題.

[1]Chatterjee S，Hadi A S，Price B.例解回歸分析[M].3版.鄭明，徐勤豐，胡瑾瑾，譯.北京：中國統(tǒng)計(jì)出版社，2004：41-43.

[2]王松桂，陳敏，陳立萍.線性統(tǒng)計(jì)模型：線性回歸與方差分析[M].北京：高等教育出版社，2004：28-33.

[3]陳杰.MATLAB寶典[M].3版.北京：電子工業(yè)出版社，2007：82-83.

[4]Recktenwald G.數(shù)值方法和MATLAB實(shí)現(xiàn)與應(yīng)用[M].伍衛(wèi)國，萬群，張輝，等，譯.北京：機(jī)械工業(yè)出版社，2004：318-320.

[5]Chatterjee S,Hadi A S.Influential Observations,High Leverage Points,and Outliers in Linear Regression[J].Statistical Science,1986,1(3):379-416.

[6]許以超.線性代數(shù)與矩陣論[M].2版.北京：高等教育出版社，2008：116-120.