趙奮軍,胡遠(yuǎn)新
(浙江省第七地質(zhì)大隊(duì), 浙江麗水市 323000)
基于MATLAB的測(cè)量數(shù)據(jù)回歸分析研究
趙奮軍,胡遠(yuǎn)新
(浙江省第七地質(zhì)大隊(duì), 浙江麗水市 323000)
闡述了回歸模型的形式,從測(cè)量數(shù)據(jù)處理角度運(yùn)用實(shí)例在MATLAB下從回歸方程的假設(shè)檢驗(yàn)和殘差分布規(guī)律 2方面進(jìn)行詳細(xì)的分析,同時(shí)預(yù)測(cè)未觀測(cè)的數(shù)據(jù)及其置信區(qū)間,減少了大量煩瑣的計(jì)算,編程代碼簡(jiǎn)單,從而使測(cè)量數(shù)據(jù)回歸分析問題變的簡(jiǎn)單易行。
回歸分析;假設(shè)檢驗(yàn);MATLAB;測(cè)量數(shù)據(jù)
回歸分析是研究 1個(gè)變量 Y與其它若干變量 X之間相關(guān)關(guān)系的 1種統(tǒng)計(jì)推斷法。它是在一組試驗(yàn)或觀測(cè)數(shù)據(jù)的基礎(chǔ)上,利用數(shù)理統(tǒng)計(jì)方法建立因變量與自變量之間的回歸關(guān)系函數(shù)表達(dá)式(稱回歸方程式),尋找被隨機(jī)性掩蓋了的變量之間的依存關(guān)系。粗略地講,可以理解為用 1種確定的函數(shù)關(guān)系去近似代替比較復(fù)雜的相關(guān)關(guān)系,這個(gè)函數(shù)稱為回歸函數(shù),在實(shí)際問題中稱為經(jīng)驗(yàn)公式?;貧w分析所研究的主要問題就是如何利用變量 X,Y的觀察值(樣本),對(duì)回歸函數(shù)進(jìn)行統(tǒng)計(jì)推斷,包括進(jìn)行估計(jì)及檢驗(yàn)與有關(guān)的假設(shè)等。
在一元線性回歸中,有 2個(gè)變量,其中 x是可觀測(cè)、可控制的普通變量,稱為自變量或控制變量,y為隨機(jī)變量,稱為因變量或響應(yīng)變量。通過散點(diǎn)圖或計(jì)算相關(guān)系數(shù)判定 y與 x之間存在著顯著的線性相關(guān)關(guān)系,即 y與 x之間存在如下關(guān)系:
其中未知參數(shù) a,b及σ2都不依賴于 x,稱為一元線性回歸模型。ε為隨機(jī)誤差或隨機(jī)干擾,是 1個(gè)分布與 x無關(guān)的隨機(jī)變量,常假定其為均值為 0的正態(tài)變量。
建立一元線性回歸模型的過程,就是利用一組觀測(cè)數(shù)據(jù) xi,yi(i=1,2,…,n)確定參數(shù) a,b的最小二乘估計(jì)值 過程,進(jìn)而得到 y關(guān)于 x的經(jīng)驗(yàn)回歸方程一元線性回歸分析的任務(wù)就是要利用這組數(shù)據(jù)求出回歸系數(shù),,并對(duì)參數(shù)和方差進(jìn)行估計(jì),并對(duì)回歸的效果進(jìn)行顯著性檢驗(yàn),從而接受回歸模型,最后在把模型用于預(yù)測(cè)和控制。
在實(shí)際的問題中,影響變量 y的因素往往不只1個(gè),而包含多種影響的多個(gè)自變量 x。通常要研究 1個(gè)因變量 y與多個(gè)自變量之間的相互關(guān)系稱為多元回歸分析,其回歸模型為:
其中 b0,b1,…,bm,σ2都是與 x1,x2,…,xm無關(guān)的未知參數(shù),ε為互相獨(dú)立的服從均值為 0,方差為σ2的正態(tài)隨機(jī)變量。
建立多元線性回歸模型的過程,就是利用一組觀測(cè)數(shù)據(jù) xi,yi(i=1,2,…,n),在最小二乘法原則下確定m+1個(gè)回歸參數(shù) b0,b1,…,bm的估值的過程,即得到m元經(jīng)驗(yàn)線性回歸方程。多元線性回歸分析的過程與一元線性回歸分析類似,即把 b0,b1,…,bm作為未知數(shù),令 X=[I xi],i=1,2,…,m作為已知系數(shù),把多元回歸模型表示成線性方程組的形式 y=X*[b0b1… bm]T,然后采用一元線性回歸分析的方法進(jìn)行參數(shù)的估計(jì)以及回歸效果的假設(shè)檢驗(yàn)。
自變量與因變量之間的關(guān)系并非都是線性的,常常會(huì)出現(xiàn)非線性關(guān)系。解決這種非線性回歸問題,一般都是通過變量的變換化為線性回歸問題:即把曲線方程化為直線方程。當(dāng)把非線性模型化為線性形式以后,就可以采用線性回歸分析方法。建立非線性回歸模型的過程:通過適當(dāng)?shù)淖兞刻鎿Q將非線性關(guān)系線性化;用線性回歸分析方法分析新變量下的線性回歸模型,求出未知參數(shù)的估計(jì)值,得到非線性回歸方程,并對(duì)其做相應(yīng)的顯著性檢驗(yàn),從而驗(yàn)證模型的嚴(yán)密性;通過新變量之間的線性相關(guān)關(guān)系反映原來變量之間的非線性相關(guān)關(guān)系。
結(jié)合文獻(xiàn)[1]中的實(shí)驗(yàn)數(shù)據(jù),就測(cè)量數(shù)據(jù)回歸分析進(jìn)行討論。按觀測(cè)數(shù)據(jù) (xi,yi)的對(duì)應(yīng)關(guān)系在MATLAB中繪出其散點(diǎn)圖,確定回歸模型。由MATLAB中的散點(diǎn)圖可知這些點(diǎn)成直線關(guān)系,故可以用線性模型進(jìn)行回歸分析。
首先在MATLAB的M文件中將觀測(cè)數(shù)據(jù) xi和yi表示成向量的形式,把回歸模型 y=a+bx表示成線性方程組的形式,令 X=[I xi],系數(shù) a和 b為未知數(shù)?,F(xiàn)將主要步驟的部分代碼及結(jié)果分析如下:
用回歸分析函數(shù) regress:[B,B int,R,R int,Stats]=regress(y′,X,0.05),其中函數(shù) regress采用的是最小二乘法進(jìn)行的回歸分析,B為返回回歸模型的系數(shù) a和 b的最小二乘估值,滿足無偏性;B int得到回歸系數(shù) a和 b的置信區(qū)間;R為觀測(cè)數(shù)據(jù)的殘差值,R int為各殘差值的置信區(qū)間;Stats得到回歸分析擬合優(yōu)度系數(shù)值和 F檢驗(yàn)值以及其對(duì)應(yīng)的概率 P值,顯著水平α=0.05。
F檢驗(yàn)值是按照 F=U/(Q/(n-2))進(jìn)行計(jì)算得到的。F≈Fα(1,n-2),在顯著水平α下,若 F>Fα(1,n-2),則認(rèn)為回歸方程效果在此水平下顯著;反之則認(rèn)為方程效果不明顯。
計(jì)算后得到的結(jié)果見表1。實(shí)驗(yàn)數(shù)據(jù)的殘差值:R=[-1.0909 1.4727-0.9636 0.6000 2.1636-1.2727 -1.7091 0.8545 0.4182 -2.0182 1.5455]。實(shí)驗(yàn)殘差值的置信區(qū)間:R int=[(-4.0065,1.8247), (-1.5249,4.4704),(-4.2261, 2.2989), (-2.8096, 4.0096),(-0.8747, 5.2019), (-4.6362, 2.0908),(-4.9274, 1.5092), (-2.5201, 4.2292),(-2.9174, 3.7538), (-4.8053, 0.7689),(-1.2356,4.3265)]
表1 參數(shù)估計(jì)及假設(shè)檢驗(yàn)計(jì)算結(jié)果(α =0.05)
殘差及其置信區(qū)間如圖1所示,由此可以確定殘差落在其置信區(qū)間內(nèi)的大致位置,也可以觀察殘差的分布變化的趨勢(shì),殘差圖越散亂代表模型的適配越好。
圖1 殘差及其置信區(qū)間
方差σ2的無偏估計(jì):,反映了回歸直線擬合的程度。
上面的計(jì)算過程都是在假定 y與 x呈現(xiàn)線性相關(guān)關(guān)系的前提下進(jìn)行的,若這個(gè)假定不成立,則建立的回歸直線方程也失去意義,為此必須對(duì) y與 x之間的線性相關(guān)關(guān)系作假設(shè)檢驗(yàn)。
(1)F檢驗(yàn) (方程顯著性檢驗(yàn))。其以方差分析為基礎(chǔ),是對(duì)回歸總體線性關(guān)系是否顯著的一種假設(shè)檢驗(yàn),是解釋模型中因變量與所有自變量之間的線性關(guān)系在總體上是否顯著的方法。返回 F=96.1798,作原假設(shè) H0∶b=0的檢驗(yàn)統(tǒng)計(jì)量,當(dāng) H0為真時(shí) F的值不應(yīng)太大,故對(duì)選定的顯著性水平α=0.05下查表 ,Fα(1,9)=5.12,則 F>Fα(1,9),故拒絕原假設(shè) H0,認(rèn)為建立的回歸方程有顯著意義。
通過 F檢驗(yàn)得到回歸方程有顯著意義,只能說明 y與 x1,x2,…,xi之間存在顯著的線性相關(guān)關(guān)系,但還不能確定影響 y的因素除了 x外是否還有 1個(gè)或幾個(gè)不可忽視的其他因素,也不能表明這個(gè)回歸方程擬合的很好,而衡量回歸方程與觀測(cè)值之間擬合好壞常用回歸分析擬合優(yōu)度系數(shù)檢驗(yàn)。
通過 F檢驗(yàn)和回歸分析擬合優(yōu)度檢驗(yàn)都驗(yàn)證了 x與 y間存在顯著的線性相關(guān)關(guān)系,而線性相關(guān)程度有多大則需要用γ相關(guān)系數(shù)檢驗(yàn) (γ檢驗(yàn))。
(3)γ相關(guān)系數(shù)檢驗(yàn) (γ檢驗(yàn))。γ相關(guān)系數(shù)反映因變量與自變量的本質(zhì)聯(lián)系,即 x與 y間線性相關(guān)程度。用γ對(duì)原假設(shè) H0∶b=0進(jìn)行檢驗(yàn),用函數(shù)corrcoef實(shí)現(xiàn) [r,P,rlow,rup]=corrcoef(x,y),按照公式 Cov(x,y)/sqrt(D(x)×D(y))計(jì)算得到,返回相關(guān)系數(shù)值γ =0.9563,P為 x和 y不相關(guān)的概率值,P=0,rlow,rup分別代表相關(guān)系數(shù)γ在 95%置信區(qū)間上限值和下限值,rlow=0.8359,rup=0.9889。γ值與γα(n-2)比較 ,γ >γα(9)=0.6022,拒絕原假設(shè),故認(rèn)為回歸效果顯著。
從殘差的分布規(guī)律及對(duì)其均值的檢驗(yàn)的角度出發(fā),對(duì)回歸效果進(jìn)行分析。由殘差分布的規(guī)律知道殘差必須服從零均值,與樣本同方差的正態(tài)分布,故需對(duì)計(jì)算得到的殘差是否服從正態(tài)分布進(jìn)行測(cè)試。
(1)殘差正態(tài)分布測(cè)試。[h,p,j,cv]=jbtest(R),在顯著性水平為 5%下,h=0表示接受殘差的分布為正態(tài)分布的假設(shè),h=1表示拒絕原假設(shè),P為接受假設(shè)的概率值,P越接近于 0,則可以拒絕是正態(tài)分布的原假設(shè);測(cè)試統(tǒng)計(jì)量的值 j大于接受假設(shè)的臨界值 cv表示拒絕假設(shè)。返回結(jié)果為 h=0,p=0.5039,正態(tài)分布檢測(cè)值 j=1.3708,檢測(cè)的臨界值 cv=5.9915。由測(cè)試的結(jié)果可以看出接受原假設(shè),即殘差的分布符合正態(tài)分布。
(2)參數(shù)估計(jì)。因?yàn)闅埐罘恼龖B(tài)分布,故可以進(jìn)一步利用殘差值求出其均值μ和標(biāo)準(zhǔn)差σ的點(diǎn)估計(jì)值和區(qū)間估計(jì):[muhat,sigmahat,muci,sigmaci]=nor m fit(R,0.05),表示在顯著性水平為5%下求殘差 R的均值 muhat和標(biāo)準(zhǔn)差 sigmahat的估計(jì)值,muci,sigmaci分別為其對(duì)應(yīng)的置信區(qū)間。返回均值 muhat=1.1304e-015≈0,置信區(qū)間 [-0.9790,0.9790],標(biāo)準(zhǔn)差 sigmahat=1.4573,置信區(qū)間為 [1.0182,2.5574]。
(3)假設(shè)檢驗(yàn)。因?yàn)棣拧玁(0,σ2),則必須對(duì)殘差ε的均值進(jìn)行假設(shè)檢驗(yàn),看是否能夠接受其原假設(shè) H0∶μ=0,備擇假設(shè) H1∶μ≠0的檢驗(yàn)。殘差的方差未知,故利用 T檢驗(yàn)法:[h,sig,ci,tstat]=ttest(R,0,0.05),在顯著性水平為 5%下對(duì)殘差 R的均值μ進(jìn)行 T檢驗(yàn),h=0表示接受原假設(shè),h=1表示拒絕原假設(shè);sig在假設(shè) H0下殘差均值出現(xiàn)的概率,sig越小 H0越值得懷疑;ci為真正均值μ的 1-α置信區(qū)間;tstat返回 3個(gè)值:T統(tǒng)計(jì)量的值、自由度和殘差標(biāo)準(zhǔn)差。返回計(jì)算的結(jié)果如下:h=0,sig=1,ci=[-0.9790,0.9790],tstat:T統(tǒng)計(jì)量的值 T=2.5727e-015≈ 0,自由度 df=10,標(biāo)準(zhǔn)差 sd=1.4573。說明接受殘差的均值為 0的假設(shè)。
最后在求得的殘差均值和標(biāo)準(zhǔn)差估計(jì)值下繪制殘差的正態(tài)分布,如圖2所示,可以更加直觀的確定殘差的分布規(guī)律,它是符合線性回歸模型所要求的殘差序列須服從與樣本等方差的正態(tài)分布。
圖2 殘差ε~N(0,σ2)正態(tài)分布
經(jīng)過假設(shè)檢驗(yàn)驗(yàn)證了 x與 y的回歸效果顯著后,就可以把回歸方程運(yùn)用于實(shí)際生產(chǎn)的預(yù)測(cè)與控制。在實(shí)際應(yīng)用中,若因變量 y比較難觀測(cè),而控制變量 x卻比較容易觀察或測(cè)量,那么根據(jù)觀測(cè)資料得到經(jīng)驗(yàn)公式后,只要觀測(cè) x就能求得 y的估計(jì)和預(yù)測(cè)值,這是回歸分析最重要的應(yīng)用之一。
該預(yù)測(cè)問題即對(duì) x的可取值范圍內(nèi)的任一個(gè)x0,作出 y的相應(yīng)估計(jì)值 y0。所謂控制,是指通過控制 y的值以便確定 x的范圍,是預(yù)測(cè)的反問題,即觀測(cè)值 y在某區(qū)間 [y1,y2]內(nèi)取值時(shí),x的控制范圍。
現(xiàn)選取一些預(yù)測(cè)點(diǎn) x0=[6,7,8,9,10],在MATLAB中很快的可以得到其相應(yīng)的預(yù)測(cè)區(qū)間,并對(duì)回歸值 y0精度給出 1個(gè)預(yù)測(cè)值的置信區(qū)間。計(jì)算數(shù)據(jù)見表2。
表2 數(shù)據(jù)預(yù)測(cè)
上面的線性回歸分析分別從回歸參數(shù)的假設(shè)檢驗(yàn)和殘差分布規(guī)律的驗(yàn)證角度 2個(gè)方面考慮分析,更全面的檢驗(yàn)了回歸模型的可行性,同時(shí)對(duì)未觀測(cè)的數(shù)據(jù)也給出了預(yù)測(cè)和區(qū)間估計(jì)。通常線性回歸分析法是最基本的分析方法,遇到非線性回歸問題都可以借助數(shù)學(xué)手段化為線性回歸問題來處理。
在回歸分析中,回歸模型采取何種形式,在沒有對(duì)所討論問題進(jìn)行全面考察的情況下是很難肯定的。通?;貧w模型受到各種因素的限制,但是模型選取的原則一定是最優(yōu)的。具體選取時(shí)應(yīng)首先要結(jié)合具體的專業(yè)理論和經(jīng)驗(yàn)給出因變量可能受影響的自變量,也可以在相關(guān)的軟件中把數(shù)據(jù)點(diǎn)描繪在坐標(biāo)系內(nèi),根據(jù)觀測(cè)數(shù)據(jù)的散點(diǎn)圖分析其大致變化趨勢(shì),然后確定回歸模型。而在MATLAB中對(duì)回歸模型進(jìn)行回歸分析擬合驗(yàn)證是非常方便的,其合理性也在本文得到論證。
[1] 劉大杰,陶本藻.實(shí)用測(cè)量數(shù)據(jù)處理方法[M].北京:測(cè)繪出版社,2000.
[2] 盛 驟,謝式千,潘承毅.概率論與數(shù)理統(tǒng)計(jì) (第三版)[M].北京:高等教育出版社,2000.
[3] 姚 東,王愛民,馮 峰,等.MATLAB命令大全 [M].北京:人民郵電出版社,2000.
[4] 張志涌.精通MATLAB 6.5[M].北京:北京航空航天大學(xué)出版社,2003.
2011-06-15)
趙奮軍 (1977-),男,陜西鳳翔人,工程師,主要從事基礎(chǔ)測(cè)繪和工程測(cè)量的研究和生產(chǎn)。