非線性回歸分析

（曲線曲面擬合）

顧世梁

2011年11月1非線性回歸分析的任務(wù)

非線性關(guān)系是最普遍的變數(shù)間量化關(guān)系，合適的非線性回歸分析對(duì)研明變數(shù)間的數(shù)量關(guān)系有重要作用。非線性回歸分析的廣泛應(yīng)用，將促使試驗(yàn)研究從定性向定量發(fā)展，由粗放向精細(xì)發(fā)展。線性關(guān)系形式單一，而非線性關(guān)系多種多樣，選擇合適的非線性模型并非易事。多項(xiàng)式也是一種（簡(jiǎn)單的一種）非線性關(guān)系，先前已有論述，本章僅討論多項(xiàng)式以外的純非線性關(guān)系。對(duì)于純非線性回歸分析，非線性回歸統(tǒng)計(jì)數(shù)的估計(jì)、假設(shè)測(cè)驗(yàn)等均有很大難度。非線性回歸分析的主要任務(wù)有下列4項(xiàng)：1)建立合適的非線性模型；2)估計(jì)非線性方程的統(tǒng)計(jì)數(shù)——曲線曲面擬合；3)合理的顯著性測(cè)驗(yàn)；4)方程的進(jìn)一步利用（插值與外推）。2非線性回歸方程的選擇主要有3種方法：1）解微分和偏微分方程組dsolve(‘Dy+y+c’,’…’)y=dsolve(‘Dy-b*y+c*y^2’,’y(0)=k/(1+a)’)symscbk;y=subs(y,c

b/k);pretty(y)2）根據(jù)機(jī)理或基本數(shù)量關(guān)系推導(dǎo)每一種函數(shù)關(guān)系都有一些基本特點(diǎn)，可以根據(jù)這些基本要素確定不同的方程。這些基本要素如零點(diǎn)(初值點(diǎn))、峰值點(diǎn)（極大、極?。?、拐點(diǎn)、漸近點(diǎn)等，應(yīng)符合數(shù)據(jù)事實(shí)。

3)試算、比較與選擇當(dāng)變數(shù)間的可能關(guān)系所知甚少，可對(duì)不同方程進(jìn)行試擬合，比較分析后選出最佳關(guān)系模型。除了前述的關(guān)鍵點(diǎn)數(shù)據(jù)應(yīng)與曲線、曲面有好的吻合外，也應(yīng)保證數(shù)據(jù)在前、中、后段都能較好地?cái)M合；另外也應(yīng)保證較高的擬合度（決定系數(shù)）、較小的離回歸平方和以及較好的插值和外推。通常，較少參數(shù)的曲線剛性有余、柔性不足，而參數(shù)較多的方程有較大的柔性。但參數(shù)太多往往會(huì)過參數(shù)化（over-parameterization），擬合的難度大大增加。3參數(shù)估計(jì)目標(biāo)函數(shù)：當(dāng)給定Xi

與Yi

(i=1,2,…,n)時(shí)，Q

也是b的函數(shù)：

Q=F(b)。擬合即為尋找βopt=min(F(b))的過程。發(fā)展穩(wěn)定高效實(shí)現(xiàn)全局最優(yōu)擬合的算法是非線性回歸的關(guān)鍵，難度較大。1）線性化法

對(duì)一些簡(jiǎn)單的方程，我們可以采用數(shù)據(jù)轉(zhuǎn)換的方式將其化成線性方程，然后用一元或多元線性回歸的方式進(jìn)行分析。如：其缺陷是該類方法僅適用于簡(jiǎn)單的方程，而絕大多數(shù)純非線性方程較復(fù)雜，不能用線性化方法進(jìn)行參數(shù)估計(jì)。2）一些通用方法梯度法(快速登山法,Gradient)；給定某一起始參數(shù)點(diǎn)：

若=0,

bj

在該點(diǎn)前后的變化不會(huì)使Q變化<0,

bj

在該點(diǎn)的增加將使Q變小

>0,

bj

在該點(diǎn)的增加將使Q變大令<0>0朝著使Q減小的方向

因而一個(gè)實(shí)例：b0=[3，20，0.5]XYf0Y-f0df/dKdf/dadf/db20.300.35896-0.058960.11965-0.01580.6320140.860.809340.050660.26978-0.029552.3639961.731.503200.226800.50107-0.03754.4999882.202.195690.004310.73190-0.029434.70937102.472.64373-0.173730.88124-0.01573.13958122.672.85830-0.188300.95277-0.006751.62009142.802.94627-0.146270.98209-0.002640.73879

-.35275-4.813e-3.16480(2)高斯法(Gauss)；(3)高斯-牛頓法(Gauss-Newton)；以新的b值再運(yùn)行前述過程，反復(fù)迭代，直至delta趨于0，或Q已不再變小。

f按多元Taylor級(jí)數(shù)展開（略去二次及二次以上各項(xiàng)）：

則目標(biāo)函數(shù)可以轉(zhuǎn)化為:令得新的優(yōu)化點(diǎn):當(dāng)b與b(0)有差異時(shí)，應(yīng)令b替代b(0)重新計(jì)算

由

△=0，或Q的前后差異小于某一定值。一個(gè)實(shí)例：b0=[3，20，0.5]XYf0Y-f0df/dKdf/dadf/db20.300.35896-0.058960.11965-0.01580.6320140.860.809340.050660.26978-0.029552.3639961.731.503200.226800.50107-0.03754.4999882.202.195690.004310.73190-0.029434.70937102.472.64373-0.173730.88124-0.01573.13958122.672.85830-0.188300.95277-0.006751.62009142.802.94627-0.146270.98209-0.002640.73879df/dKdf/dadf/dbY-f00.1197-0.01580.63201-0.058960.2698-0.02962.363990.050660.5011-0.03754.499980.22680.7319-0.02944.709370.004310.8812-0.01573.13958-0.173730.9528-0.00681.62009-0.18830.9821-0.00260.73879-0.14627XY得新的優(yōu)化點(diǎn):反復(fù)迭代…(4)改良高斯牛頓法(Levenberg-Marquardt)

這是梯度法和高斯-牛頓法相結(jié)合的一種方法。A很可能是奇異的，需對(duì)此陣進(jìn)行調(diào)整:作用:一可解決A陣奇異，無法求解Δ之困；二是A陣對(duì)角線元素包含了較大的與求解Δ相關(guān)的信息量，加快趨于全局優(yōu)的進(jìn)程。

(5)極大似然法(maximumlikelihood)。 大多數(shù)著名的統(tǒng)計(jì)軟件如SAS,Matlab,Sigmaplot等包含了基于這些算法的非線性方程擬合模塊。3）上述通用算法存在的問題：(1)需提供方程的導(dǎo)數(shù)或偏導(dǎo)數(shù)；(2)需提供合適的初值；(3)一般難于實(shí)現(xiàn)全局最優(yōu)擬合。最后一點(diǎn)往往是最主要、最致命的缺陷。4）曲線、曲面擬合新算法

(Contraction-ExpansionAlgorithm)CE算法包含三個(gè)基本步驟：(1)收縮步，縮小步長(zhǎng)的搜索過程；(2)擴(kuò)張步，擴(kuò)大步長(zhǎng)的搜索過程；(3)調(diào)整步，中心點(diǎn)、臨界值的重新調(diào)整。(1)收縮步(2)擴(kuò)張步(3)中心點(diǎn)和步長(zhǎng)的確定全局最優(yōu)擬合的能力和效率很大程度上取決于初始點(diǎn)和步長(zhǎng)，初始步長(zhǎng)一般總不是很合適的，必須由尋優(yōu)過程的信息反饋調(diào)整。記錄在尋優(yōu)搜索過程中的度點(diǎn)（即滿足一定要求的參數(shù)點(diǎn)）的數(shù)量和位置，算出它們平均數(shù)和標(biāo)準(zhǔn)差（Sj為bj的二階原點(diǎn)矩）：(4)調(diào)整臨界值C

若

C很小，產(chǎn)生的度點(diǎn)數(shù)量太少，若C很大，產(chǎn)生的度點(diǎn)數(shù)量太多，這些情形都將使算法的能力和效率受損。臨界值C必須有反饋調(diào)節(jié)機(jī)制。若N是每一輪次的試算節(jié)點(diǎn)總數(shù)，nE是擴(kuò)張步一個(gè)循環(huán)（由3~7個(gè)輪次組成）的度點(diǎn)數(shù)量，在一次循環(huán)后重新計(jì)算臨界值（包括步長(zhǎng)）。前后兩次循環(huán)(v,v+1)使用不同的公式是為了減少循環(huán)過程波浪形C值的發(fā)生。在mod(v,2)=0時(shí)，需將nE清零。（5）CE算法的主要優(yōu)缺點(diǎn)：

不必提供導(dǎo)數(shù)與偏導(dǎo)數(shù)，利于通用程序的編制；

無需提供適合的初值；

實(shí)現(xiàn)最優(yōu)擬合的能力較強(qiáng)；

搜索效率不高，對(duì)多參數(shù)非線性問題難于實(shí)施。

(6)縮張算法的一些改進(jìn)：

1每一輪次的試算節(jié)點(diǎn)數(shù)(z)隨p的增加而指數(shù)（爆炸）式增長(zhǎng)。5步點(diǎn)時(shí)，z=5^p=exp(1.60944p)；在3步點(diǎn)時(shí)，z=3^p=exp(1.09861p)。因此在p>7（5步點(diǎn)）或p>13（3步點(diǎn)）時(shí)，算法負(fù)荷量已超出普通pc機(jī)的上限（每輪次試算節(jié)點(diǎn)數(shù)以1m計(jì)），該法不適宜用于參數(shù)數(shù)p>15的非線性方程的擬合。在參數(shù)較多(p>8)時(shí)，只在p

維參數(shù)空間中均勻隨機(jī)布點(diǎn)，試算節(jié)點(diǎn)數(shù)在基礎(chǔ)條件下隨p的增加而呈多項(xiàng)式（二次式）增長(zhǎng)（z=300+25*p^2），這比指數(shù)式增長(zhǎng)大為減少，使多參數(shù)復(fù)雜非線性問題的擬合成為可能。

2與解析法中的改良高斯牛頓法相結(jié)合，在給定的參數(shù)初值（或中間值）點(diǎn)處，利用參數(shù)微小差量Δ的差分方程獲得方程對(duì)某一參數(shù)的近似偏導(dǎo)函數(shù)值，再將各（觀察值）點(diǎn)的偏導(dǎo)函數(shù)值的乘積累加，得到近似的Jaccobi矩陣(A，或A*)和常數(shù)陣K，再由AΔ=K，解出Δ=A-1K（Δ=A*-1K），當(dāng)Δ接近0或RSS(Q)小于收斂標(biāo)準(zhǔn)時(shí)結(jié)束。

f依第j個(gè)參數(shù)bj的近似偏導(dǎo)數(shù)為：是Xi及參數(shù)點(diǎn)bj(0)處僅第bj

參數(shù)具微小差值時(shí)的回歸值。是bj

微小差值參數(shù)增量；bi(0)基于數(shù)值微分基礎(chǔ)的改良高斯-牛頓法：當(dāng)b與b(0)有差異時(shí)，應(yīng)令b替代b(0)重新計(jì)算Δ，當(dāng)Δ接近0或小于收斂標(biāo)準(zhǔn)時(shí)結(jié)束。構(gòu)建矩陣4非線性回歸統(tǒng)計(jì)數(shù)的假設(shè)測(cè)驗(yàn)

Jaccobi陣A的逆陣C(C=A-1)對(duì)角線元素為相應(yīng)回歸統(tǒng)計(jì)數(shù)標(biāo)準(zhǔn)化的方差，所謂標(biāo)準(zhǔn)化的方差是指離回歸誤差方差為1時(shí)的方差。因此，第j個(gè)回歸統(tǒng)計(jì)數(shù)bj（與0的差異顯著性）