非線性最小二乘

目錄

1 非線性最小二乘估計

在經(jīng)典最小二乘法估計中，假定被解釋變量的條件期望是關(guān)于參數(shù)的線性函數(shù)，例如
$$E(y|x)=a+bx$$
其中$a,b$為待估參數(shù)，$E(y|x)$是關(guān)于參數(shù)$a,b$的線性函數(shù)。但$E(y|x)$是關(guān)于參數(shù)的非線性函數(shù)，則利用ols求出的正規(guī)方程組沒有解析解。只能通過相關(guān)數(shù)值計算。考慮一個簡單的非線性模型
$$Y_i=\beta X_{1i}+{\beta }^2{X}_{2i}+{\epsilon }_i$$
其中擾動項${\epsilon }_i$滿足$\mathrm{E}\left({\epsilon }_i\right)=0,\mathrm{var}\left({\epsilon }_i\right)={\sigma }^2$,且為獨(dú)立同分布。其殘差平方和為
$$\begin{array}{rl}S(\beta )& =\sum _{i=1}^n{\epsilon }_i^2=\sum _{i=1}^n{\left[Y_i?f\left(X_i,\beta \right)\right]}^2\\ & =\sum _{i=1}^n{\left[Y_i?\beta X_{1i}?{\beta }^2{X}_{2i}\right]}^2\end{array}$$
為了使回歸線盡可能接近觀測值，要求殘差平方和最小。根據(jù)微積分的知識
$$\begin{array}{rl}\frac{\mathrm{d}S}{\mathrm{d}\beta }& =2\sum _{i=1}^n\left[Y_i?f\left(X_i,\beta \right)\right]\left(?\frac{\mathrm{d}f\left(X_i,\beta \right)}{\mathrm{d}\beta }\right)\\ & =2\sum _{i=1}^n\left[Y_i?\beta X_{1i}?{\beta }^2{X}_{2i}\right]\left[?X_{1i}?2\beta X_{2i}\right]=0\end{array}$$
整理得：
$$2{\beta }^3\sum _{i=1}^n{X}_{2i}^2+3{\beta }^2\sum _{i=1}^n{X}_{1i}{X}_{2i}+\beta \left(\sum _{i=1}^n{X}_{1i}^2?2\sum _{i=1}^n{X}_{2i}{Y}_i\right)?\sum _{i=1}^n{X}_{1i}{Y}_i=0$$
這是關(guān)于參數(shù)$\beta$的三次函數(shù)。盡管三次函數(shù)存在解析解（利用卡丹或盛金公式），其結(jié)果極為復(fù)雜。若上述三次方程存在實根${\beta }_i(i=1,2,3)$（最多三個），則將${\beta }_i$代入殘差平方和，取$S(\beta )$最小所對應(yīng)的${\beta }_i$。上述例子中，被解釋變量條件期望是關(guān)于參數(shù)的二次函數(shù)，如果將這種函數(shù)形式改為指數(shù)、對數(shù)或三角函數(shù)形式，則一般不存在解析解。

因此，數(shù)值分析自然成為解決上述問題的有力武器?？紤]一般化的非線性回歸問題，設(shè)總體回歸模型滿足
$$Y=f(X,\beta )+\epsilon$$
對應(yīng)的殘差平方和為
$$S(\beta )=\sum _{i=1}^n{\left[Y_i?f\left(X_i,\beta \right)\right]}^2$$
要使其最小化，需要滿足一階條件
$$\frac{\mathrm{d}S}{\mathrm{d}\beta }=?2\left[\sum _{i=1}^n\left[Y_i?f\left(X_i,\beta \right)\right]\left(?\frac{\mathrm{d}f\left(X_i,\mathbf{\beta }\right)}{\mathrm{d}\beta }\right)\right]=0$$
顯然，上述問題不存在解析解，因此考慮對$f(X_i,\beta )$進(jìn)行一階泰勒展開。設(shè)參數(shù)向量$\beta$的初始值為${\beta }_1$,則可以在${\beta }_1$附近找到函數(shù)$f(X_i,\beta )$使得
$$f\left(X_i,\beta \right)\approx f\left(X_i,{\beta }_1\right)+\frac{\mathrm{d}f\left(X_i,\beta \right)}{\mathrm{d}\beta }{\mid }_{\beta ={\beta }_1}\left(\beta ?{\beta }_1\right)$$
記 d f ( X i , β ) d β ∣ β 1 ≈ f ( X i , β ) ? f ( X , β ) β ? β 1 \left.\frac{\mathrmn5n3t3z f\left(X_{i}, \beta\right)}{\mathrmn5n3t3z \beta}\right|_{\beta_{1}} \approx \frac{f\left(X_{i}, \beta\right)-f(X, \beta)}{\beta-\beta_{1}} dβdf(Xi?,β)? ?β1??≈β?β1?f(Xi?,β)?f(X,β)?，簡記 X ~ i ( β 1 ) = d f ( X i , β ) d β ∣ β 1 \widetilde{X}_{i}\left(\beta_{1}\right)=\left.\frac{\mathrmn5n3t3z f\left(X_{i}, \beta\right)}{\mathrmn5n3t3z \beta}\right|_{\beta_{1}} X i?(β1?)=dβdf(Xi?,β)? ?β1??，則
$$\begin{array}{rl}S(\beta )& =\sum _{i=1}^n{\left[Y_i?f\left(X_i,{\beta }_1\right)?X_i\left({\beta }_1\right)\left(\beta ?{\beta }_1\right)\right]}^2\\ & =\sum _{i=1}^n{\left[Y_i\left({\beta }_1\right)?X_i\left({\beta }_1\right)\beta \right]}^2\end{array}$$
其中
$$Y_i\left({\beta }_1\right)=Y_i?f\left(X_i,{\beta }_1\right)+X_i\left({\beta }_1\right){\beta }_1$$
給定初始值向量${\beta }_i$，則$Y_i\left({\beta }_1\right)$與$X_i\left({\beta }_1\right)$可計算，從而求出最小殘差平方和。$S(\beta )$對應(yīng)的回歸方程為
$$Y_i\left({\beta }_1\right)=X_i(\beta )\beta +{\epsilon }_i$$
最小二乘估計量為
$${\beta }_2={\left[X{\left({\beta }_1\right)}^{\prime }X\left({\beta }_1\right)\right]}^{?1}X{\left({\beta }_1\right)}^{\prime }Y\left({\beta }_1\right)$$
其中
X ~ ( β 1 ) = [ X ~ 1 ( β 1 ) ? X ~ n ( β 1 ) ] , Y ^ ( β 1 ) = [ Y ~ 1 ( β 1 ) ? Y ~ n ( β 1 ) ] \widetilde{X}\left(\beta_{1}\right)=\left[\begin{array}{c} \widetilde{X}_{1}\left(\beta_{1}\right) \\ \vdots \\ \widetilde{X}_{n}\left(\beta_{1}\right) \end{array}\right], \quad \hat{Y}\left(\beta_{1}\right)=\left[\begin{array}{c} \widetilde{Y}_{1}\left(\beta_{1}\right) \\ \vdots \\ \widetilde{Y}_{n}\left(\beta_{1}\right) \end{array}\right] X (β1?)= ?X 1?(β1?)?X n?(β1?)? ?,Y^(β1?)= ?Y 1?(β1?)?Y n?(β1?)? ?
此時我們求出${\beta }_2$,再將${\beta }_2$作為初始值依次迭代計算，得到關(guān)于向量參數(shù)${\beta }_i$的一個序列，當(dāng)且僅當(dāng)
$$||{\beta }^{(k+1)}?{\beta }^{(k)}||<\delta$$
其中$\delta >0$為事先預(yù)定的絕對誤差。不難得到，參數(shù)$\beta$滿足遞推關(guān)系
$$\begin{array}{rl}{\mathbf{\beta }}_{n+1}& ={\left[X{\left({\mathbf{\beta }}_n\right)}^{\prime }X\left({\mathbf{\beta }}_n\right)\right]}^{?1}X{\left({\mathbf{\beta }}_n\right)}^{\prime }Y\left({\mathbf{\beta }}_n\right)\\ & ={\left[X{\left({\mathbf{\beta }}_n\right)}^{\prime }X\left({\mathbf{\beta }}_n\right)\right]}^{?1}X{\left({\mathbf{\beta }}_n\right)}^{\prime }\left[\mathbf{Y}?f\left(X,{\mathbf{\beta }}_n\right)+X\left({\mathbf{\beta }}_n\right){\mathbf{\beta }}_n\right]\\ & ={\mathbf{\beta }}_n+{\left[X{\left({\mathbf{\beta }}_n\right)}^{\prime }X\left({\mathbf{\beta }}_n\right)\right]}^{?1}X{\left({\mathbf{\beta }}_n\right)}^{\prime }\left[Y?f\left(X,{\mathbf{\beta }}_n\right)\right]\end{array}$$
通過證明，隨著樣本容量$n\to \infty$,參數(shù)$\beta$估計量服從漸進(jìn)正態(tài)分布，即
$$\beta ～N\left(\beta ,{\hat{\sigma }}^2{\left[X(\beta {)}^{\prime }X(\beta )\right]}^{?1}\right),{\hat{\sigma }}^2=\frac{S(\beta )}{n?1}$$

3 非線性最小二乘的實現(xiàn)

在R語言中，可以適用nls函數(shù)實現(xiàn)非線性最小二乘法。以C-D函數(shù)為例，

設(shè)一國產(chǎn)出取決于資本、勞動與全要素的投入，即
$$Y=A{K}^{\alpha }{L}^{\beta }\mu$$
下面通過R代碼運(yùn)行實現(xiàn)對參數(shù)$\alpha ,\beta$的估計

t = 1:12 #時間設(shè)定
Y=c(26.74, 34.81, 44.72, 57.46, 73.84, 88.45, 105.82,126.16, 150.9, 181.6, 204.3, 222.8) #產(chǎn)出序列
K=c(23.66,30.55,38.12,46.77,56.45,67.15,78.92,91.67,105.5, 121.3, 128.6, 132.5) #資本序列
L=c(26, 28, 32, 36, 41, 45, 48, 52, 56, 60, 66, 70) #勞動投入序列
Cdnls <- nls(Y~A*K^a*L^b,start = list(A = 0.1,a = 0.5,b = 0.5)) #非線性最小二乘,start為參數(shù)初始值向量
summary(Cdnls)
#-------------------運(yùn)行結(jié)果---------------------------
#Formula: Y ~ A * K^a * L^b
Parameters:
  Estimate Std. Error t value Pr(>|t|)    
A   0.1129     0.0159    7.12  5.6e-05 ***
a   0.6568     0.0652   10.07  3.4e-06 ***
b   1.0298     0.1044    9.86  4.0e-06 ***
---
Signif. codes:  
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.7 on 9 degrees of freedom

Number of iterations to convergence: 9 
Achieved convergence tolerance: 7.55e-06

結(jié)果顯示，參數(shù)$\alpha =0.6568$,$\beta =1.0298$。對比直接取對數(shù)的OLS,即估計
$$lnY=lnA+\alpha lnK+\beta lnL+e$$

CDlm <- lm(log(Y)~log(K)+log(L))  #對數(shù)形式
summary(CDlm)
#--------------------運(yùn)行結(jié)果--------------
Call:
lm(formula = log(Y) ~ log(K) + log(L))

Residuals:
     Min       1Q   Median       3Q      Max 
-0.02714 -0.00595 -0.00118  0.00764  0.02557 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -2.0737     0.2355   -8.80  1.0e-05 ***
log(K)        0.6258     0.0916    6.83  7.6e-05 ***
log(L)        1.0379     0.1621    6.40  0.00012 ***
---
Signif. codes:  
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.0173 on 9 degrees of freedom
Multiple R-squared:     1,	Adjusted R-squared:  0.999 
F-statistic: 9.16e+03 on 2 and 9 DF,  p-value: 1.29e-15

結(jié)果顯示，參數(shù)$\alpha =0.6268$,$\beta =1.0379$。因此，CD函數(shù)對數(shù)化的結(jié)果回歸與非線性最小二乘回歸的參數(shù)基本一致。但一些不能對數(shù)化的方程，非線性最小二乘的作用更為明顯?？紤]真實模型
$$y=2sin(x)+4cos(x)$$
接下來我們進(jìn)行仿真模擬

set.seed(123) #隨機(jī)種子
x <- seq(1,100,by = 0.1) #1-100，步長為0.1
e <- rnorm(length(x),0,1) #長度為序列x的長度，服從標(biāo)準(zhǔn)正態(tài)分布的誤差
y <- 2*sin(x)+4*cos(x)+e #實際觀測的被解釋變量
plot(x,y,type = "o") #打印散點(diǎn)圖

nls1 <- nls(y~a*sin(x)+b*cos(x),
            start = list(a = 0,b =0)) #非線性最小二乘，初始值設(shè)定為0，0
nls1
#-------------運(yùn)行結(jié)果------------------
Nonlinear regression model
  model: y ~ a * sin(x) + b * cos(x)
   data: parent.frame()
   a    b 
1.92 4.03 
 residual sum-of-squares: 974

Number of iterations to convergence: 1 
Achieved convergence tolerance: 6.73e-10

結(jié)果顯示估計量$a=1.92$,$b=4.03$，與總體參數(shù)$a=2,b=4$即為接近

參考文獻(xiàn)

王斌會(2015).計量經(jīng)濟(jì)學(xué)建模及R語言應(yīng)用[M].北京大學(xué)出版社文章來源地址http://www.zghlxwxcb.cn/news/detail-440166.html

到了這里，關(guān)于非線性最小二乘的文章就介紹完了。如果您還想了解更多內(nèi)容，請在右上角搜索TOY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關(guān)文章，希望大家以后多多支持TOY模板網(wǎng)！