張子賢
(江蘇建筑職業技術學院,江蘇 徐州 221116)
指數函數擬合公路隧道工程沉降規律的方法研究
張子賢
(江蘇建筑職業技術學院,江蘇 徐州 221116)
對于指數函數回歸,只當采用乘積隨機誤差時才能夠線性化。導出了采用乘積隨機誤差及采用線性化回歸方法時,指數函數因變量的數學期望的表達式,該式表明,該因變量的估計值并非是其數學期望的估值。分析表明,采用線性化回歸方法所求指數函數的回歸系數不滿足該因變量的殘差平方和為最小。基于上述不合理現象,對指數函數的回歸計算應采用非線性回歸方法求解。文中給出了采用高斯-牛頓法或借助MATLAB軟件中nlinfit函數求解指數函數非線性回歸的方法。實例進一步表明,采用非線性回歸方法擬合效果顯著優于線性化的回歸方法,且借助MATLAB軟件易于實現。
指數函數;線性化回歸方法;非線性回歸方法;沉降規律擬合;擬合精度
指數回歸函數(以下簡稱指數函數)常用于擬合路基、隧道等建筑物的沉降量隨時間的變化規律[1-3]。對該函數的回歸計算,數理統計教科書、以往科技文獻以及生產實際中常用的方法是:首先,通過變量代換轉化為線性模型,并利用線性回歸方法推求線性回歸系數;然后,根據線性回歸系數反求指數函數的回歸系數。這種方法稱為非線性模型線性化回歸方法。該方法看起來是合理的,其實不然。實際上,只當采用乘積隨機誤差時指數函數才能夠線性化,本文將分析指數函數采用線性化回歸方法的統計特性、分析指數函數非線性回歸與其線性化回歸二者的殘差平方和之間的關系式,進而指出指數函數采用線性化回歸方法的不合理之處。給出了指數函數采用高斯-牛頓法進行非線性回歸的方法或借助MATLAB軟件中nlinfit函數的回歸方法。實例表明,對于指數函數,采用非線性回歸方法的擬合效果顯著優于線性化的回歸方法。
采用回歸分析方法,對因變量的值進行估計,實質上是對其數學期望進行估計。因此,因變量的估計值應等于其數學期望的估計值,這是回歸方法應具備的統計特性。為分析指數函數采用線性化回歸方法時的統計特性,必須引入隨機誤差項ε。
當采用乘積隨機誤差時[4],指數函數回歸模型:

式中:A、B為待估參數。
對式(1)進行自然對數變換后,并令v=lny,a=lnA,得線性模型

傳統的做法是假設 ε~N(0,σ2),利用線性最小二乘法求解式(2)的線性回歸系數,進而得式(1)的非線性回歸系數,為簡便起見且不致混淆,將也簡記為ˇ(下同)。因此,將作為y的估計值。然而并非是y的數學期望的估計值。因為:

文獻[5]已證明:E(eε)=1+σ2/2+σ4/8
由式(3)可知ˇ,對數變換后y的數學期望不為AeBx,而估計值=AeBx不等于y的數學期望的估值,顯然不是好的估計。
對于指數函數,若采用加性隨機誤差:

對于式(4),設ε~N(0,σ2),則E(y)=AeBx。因此,指數函數因變量的估計值:

一般地,設可線性化的曲線回歸變量代換后的線性回歸的因變量v的殘差平方和為曲線回歸因變量的殘差平方和為,文獻[5]導出二者的近似關系式為

對于指數函數,v=lny,即y=ev,dy/dv=ev,則有:

綜上所述,為取得好的擬合效果以及提高擬合精度,對于指數函數應采用非線性回歸方法求解。
需要指出,反映線性化回歸密切程度的指標常用相關系數r或r2:

反映非線性回歸模型的相關指數R2(也稱為決定系數)

顯然,當對非線性待估變量y作了變換時,r2≠R2。
對于指數函數回歸模型式(4)、式(5),引入參數向量 θ=(A,B)。設 y和 x具有 n組觀察值(xi,yi),i=1~n。采用高斯-牛頓法[6]求解指數函數回歸參數“最小二乘”估計的參數遞推公式,寫成矩陣形式為:


采用高斯-牛頓法對指數函數非線性回歸計算的步驟如下。
(1)分別對式(5)中參數 A、B 求偏導數,得:

(2)確定待估參數初值 θ(0)=(A0,B0)′,可利用實測值中任兩組關系值求得。
(3)利用 θ(0)、式(12)、式(13)及 n 組實測值(ti,y)i,i=1~n,計算偏導數矩陣J(θ(0))及θ(0),進而根據式(10)進行迭代計算,直到 θ(k)收斂穩定,即小于或等于預先指定的小正數(例如δ=0.000 5),從而得到非線性回歸系數A,B的估計值。
實例1:柳山隧道左洞K39+630拱頂及該隧道斷面內空收斂1#、2#、3#測線的累計沉降量y的觀測數據見表1[1]。
根據實例1各測點的散點圖,采用指數函數估計累計沉降量:

令x=1/t,式(14)則轉化為式(5)。
對實例1分別采用線性化回歸方法、非線性回歸的高斯-牛頓法,擬合拱頂及該隧道斷面內空收斂1#、2#、3#測線的累計沉降量,其結果見表2。擬合曲線分別見圖1~圖4。可見,高斯-牛頓法擬合效果顯著優于線性化回歸方法。
對指數函數非線性回歸計算也可直接調用MATLAB軟件中nlinfit函數,以實例1拱頂為例,具體方法如下[7,8]。
function yhat=volumsq(beta,t)
yhat=beta(1)*exp(beta(2)./t)

表1 實例1累計沉降量觀測數據

表2 實例1不同方法求得指數函數回歸的計算結果
(2)在命令窗口輸入
在命令窗口輸入
t=[1.007,2.313,3.021,4.007,4.986,6.149,7.18,7.975,9.027,10.037,10.999,12.02,13.117,13.978,15.159,15.996,17.006,18.201,19.218,20.017,21.221,22.256]y=[0.432,0.829,1.182,1.577,2.028,2.393,2.584,2.833,3.033,3.215,3.314,3.462,3.578,3.661,3.75 3.799,3.799,3.908,3.906,3.912,3.907,3.92]beta0=[3,-1]′[beta]=nlinfit(t′,y′,′volumsq′,beta0);beta

圖1 實例1的拱頂由不同方法所得擬合曲線比較

圖2 實例1的1#測線由不同方法所得擬合曲線比較

圖3 實例1的2#測線由不同方法所得擬合曲線比較

圖4 實例1的3#測線由不同方法所得擬合曲線比較
得結果:beta=[4.887 8,-4.303 7]。
實例2:根據文獻[2]坡頭公路隧道拱頂沉降觀測數據,分別采用線性化回歸方法、非線性回歸高斯-牛頓法進行指數函數回歸計算,結果見表3。擬合曲線見圖5。可見高斯-牛頓法擬合效果顯著優于線性化回歸方法。利用高斯-牛頓法求得的回歸方程,計算回歸線的均方誤為0.169 1 mm、殘差絕對值和為4.927 6,分別小于文獻[2]相應的計算結果 0.191 2 mm、6.021 6。

表3 實例2不同回歸方法求得指數函數回歸的計算結果
對于指數函數回歸,只當采用乘積隨機誤差時才能夠線性化。導出了采用乘積隨機誤差及采用線性化回歸方法時,指數函數因變量的數學期望的表達式,該式表明,該因變量的估計值并非是其數學期望的估值。得出了指數函數非線性回歸與其線性化回歸二者的殘差平方和之間的關系式,該式表明,采用線性化回歸方法所求得的非線性回歸系數不滿足指數函數因變量的殘差平方和為最小。因此,對指數函數的回歸計算應采用非線性回歸方法求解。

圖5 實例2由不同方法所得擬合曲線比較
文中給出了采用高斯-牛頓法或借助MATLAB軟件中nlinfit函數求解指數函數非線性回歸的方法。結合兩則實例5個測點的沉降觀測數據,采用非線性回歸方法進行指數函數的回歸計算,結果進一步表明,非線性回歸方法的擬合效果顯著優于線性化回歸方法。
[1]朱建宇,鄢志輝.回歸分析在隧道量測數據處理中的應用[J].湖南城市學院學報(自然科學版),2009,18(1):28-31.
[2]郭云開,李亮,崔曉如.公路隧道監測數據非線性回歸的穩健估計[J].公路工程,2010,35(3):140-144.
[3]胡達,劉杰,楊慶光,等.公路隧道圍巖變形規律及回歸分析[J].公路工程,2010,35(1):5-8.
[4]盛驟,謝式千,潘承毅.概率論與數理統計(第二版)[M].北京:高等教育出版社,1994.
[5]張子賢.可線性化的非線性回歸的有關問題與幾種回歸方法的比較[J].數學的實踐與認識,2015,45(18):167-173.
[6]袁志發,周靜芋.多元統計分析[M].北京:科學出版社,2002.
[7]周品.MATLAB概率與數理統計 [M].北京:清華大學出版社,2012.
[8]謝中華,李國棟,劉煥進,等.MATLAB從零到進階[M].北京:北京航空航天大學出版社,2012.
U456.3
A
1009-7716(2017)12-0157-04
10.16799/j.cnki.csdqyfh.2017.12.045
2017-08-07
住房與城鄉建設部項目(2015-K7-009)。
張子賢(1958-),女,河北豐南人,教授,從事應用概率統計等方面的研究工作。