SA  >> Vol. 8 No. 4 (August 2019)

    光滑樣條回歸及應用研究
    The Study on Smooth Spline Regressionand Its Application

  • 全文下載: PDF(765KB) HTML   XML   PP.604-612   DOI: 10.12677/SA.2019.84069  
  • 下載量: 85  瀏覽量: 156  

作者:  

王鳳雪:北方工業大學理學院,北京

關鍵詞:
光滑樣條算法光滑度 Smooth Spline Algorithm Smoothness

摘要:

本文通過介紹非參數估計方法光滑樣條回歸的模型、算法、案例分析對光滑樣條回歸進行研究。以R語言中ISLR包中Wage數據集為例,分別運用樣條回歸和多項式進行實證分析,擬合研究工資與年齡的關系,結果顯示樣條回歸結果優于多項式回歸。

This paper studies the smooth spline regression by introducing the model, algorithm and the case for analysis of the nonparametric estimation method. Taking the Wage data set in the ISLR package in R language as an example, spline regression and polynomial are respectively used for empirical analysis, and the relationship between Wage and age is fitted. The results show that spline regression results are superior to polynomial regression.

1. 引言

在實際應用中,我們常常對總體進行某種分布的假設,抽樣得到樣本信息,去估計總體參數,這種方法稱為參數估計方法。但當對總體信息一無所知,或不假定總體分布形式,只通過樣本信息對總體參數進行估計,此時,非參數估計就展現了很強的靈活性。

非參數回歸分為局部回歸、光滑樣條回歸、正交回歸。光滑樣條回歸,因其在抽取樣本對總體進行回歸時,不必依賴總體分布形式,在減小誤差、提高預測精確度、提高擬合曲線的光滑度上都體現了良好的特性。

在諸多方法中,三次光滑樣條因其計算簡便、結果有良好的統計性質而應用廣泛。在光滑樣條的學術研究上,陳長生、徐勇勇(1997)主要重點研究了兒童成長曲線的光滑樣條擬合 [1] ,回歸結果優于參數估計方法。2000年其在此基礎上進行加權光滑樣條的改進 [2] ,改進后模型效果更好,穩定性更強。盧一強、陳中威(2010)主要研究了基于光滑樣條的選擇方法,提出了部分線性模型中的非參數函數部分的假設檢驗是否為多項式的假設檢驗方法 [3] 。本文通過介紹光滑樣條回歸模型、算法、案例分析,基于R語言的實現,對光滑樣條進行研究以及光滑樣條回歸與多項式回歸結果的比較分析。

2. 理論模型及方法論

2.1. 非參數回歸模型的一般形式及模型

設Y為因變量 X 1 , X 2 , ? , X p 為自變量,非參數回歸模型的一般形式為

Y = η ( X 1 , X 2 , ? , X p ) + ε

其中對p元回歸函數只作一些連續性或光滑性的要求。由于非參數回歸模型不假定回歸函的具體形式而增加了模型的靈活性和適應性。

( y i , x i ) Math_12#為來自總體 ( Y , X ) 的一個樣本容量為n的獨立同分布的樣本,需要基于觀測值 ( y i , x i ) Math_15#估計 η ( x ) 并進行有關的統計推斷。

非參數回歸的模型為:

E ( Y | X = x ) = g (X)

數據和模型是統計分析的兩個信息來源,數據帶有“噪聲”,但無偏,而模型實際上是種約束,有助于降低噪聲,是響應的。在“偏差–方差”的平衡表上,代表兩個極值的分別是標準參數模型和無約束非參數模型。在兩個極值之間,存在著大量的非參數或半參數模型,其中大多數被稱為平滑方法。非參數估計族可通過懲罰似然法導出各種隨機環境下的模型 [4] 。

2.2. 三次平滑模型樣條

光滑樣條回歸實際上是一種局部建模方法,是按照一定的光滑性連接起來的分段多項式 [3] 。首先先介紹多項式樣條回歸估計,多項式樣條估計的最小化式為:

s ( f ) = i = 1 T ( Y i ? η ( X i ) )

在實際生活中,三階樣條比較常用,原理如下:

設某區間 [ a , b ] 上有實數 t 1 , t 2 , ? , t n 且滿足 a < t 1 < t 2 < ? < t n < b f ( x ) 是定義在區間 [ a , b ] 的函數,如 f ( x ) 滿足以下條件[兒童參考]:

1) 在 [ a , t 1 ] , ( t 1 , t 2 ] , ? , ( t n , b ] 上,函數 f ( x ) 為三次多項式。

2) 函數 f ( x ) 及其二階導數在 t i ( i = 1 , 2 , ? , n ) 處處連續。

則稱這樣的分段多項式函數為三次樣條函數, t i 為樣條函數的節點。令 t 0 = a , t n + 1 = b

η ( x ) = d i ( x ? t i ) 3 + c i ( x ? t i ) 2 + b i ( x ? t i ) + a i , t i < x < t i + 1

為三次光滑樣條的表達式。光滑樣條估計的基本思想是尋找一個光滑函數使得殘差平方和最小,所以引入懲罰函數,使得懲罰平方和最小,估計方法為懲罰最小二乘估計,表達式為:

min s ( f ) = i = 1 T { y i ? η ( x i ) } 2 + λ a b { η ( x ) } 2 d x

損失函數 i = 1 T { y i ? η ( x i ) } 2 是使 η 能很好的擬合數據,而 λ a b { η ( x ) } 2 d x 則對函數 η 的波動性進行懲罰,二

階導數 η ( x ) 對應了斜率的變化程度,衡量的是函數的粗糙度。 λ 為需要選擇的光滑參數也稱拉格朗日乘子,此方法可以避免多項式樣條估計的節點選擇對非參數擬合曲線的光滑程度影響。 λ 值越大,函數 η 越光滑。

3. 光滑樣條基本算法

1) 目標函數均方誤差最小

min s ( f ) = i = 1 T { y i ? η ( x i ) } 2 + λ a b { η ( x ) } 2 d x

2) 寫成矩陣形式,其中 η ( x ) = j = 1 N N j ( x ) θ j

s ( θ , λ ) = ( y ? N θ ) T ( y ? N θ ) + λ θ T Ω N θ

其中 { N } i j = N j ( x i ) { Ω N } i j = N j N k d t

3) 運用最小二乘法估計

θ ^ = ( N T N + λ Ω N ) ? 1 N T y

4) 帶入擬合函數

η ^ = j = 1 N N j ( x ) θ ^ j

設B為N*M的矩陣,

η ^ = B ( B T B ) ? 1 B T y = H y

5) 求光滑參數 λ

d f λ = t r a c e (Sλ)

S λ = N ( N T N + λ Ω N ) N T = N ( N T [ I + λ N ? T Ω N N ? 1 ] N ) ? 1 N T = ( I + λ N ? T Ω N N ? 1 ) ? 1

矩陣 S λ 可以寫成 S λ = ( I ? λ K ) ? 1

此時 min s ( η ) = ( y ? η ) T ( y ? η ) + λ η T K η η ^ = S λ y

矩陣S具有對稱半定性質,對其進行特征分解:

S λ = k = 1 N ρ k ( λ ) u k u K T

其中, ρ k ( λ ) = 1 1 + λ d k ,這里的 d k 是矩陣K的特征值

6) 估計樣條函數

η ^ = S λ y = k = 1 N ρ k ( λ ) u k u K T y

4. 基于R語言的隨機模擬比較研究

光滑樣條的擬合用R語言的smoothing.spline()實現。隨機模擬實驗運用R語言,隨機生成變量X,Y,運用光滑樣條方法回歸擬合X與Y之間的關系。

#模擬實驗

>set.seed(1)

>ei<-rnorm(4000,0,3) #生成均值為0,標準差為1的4000條殘差序列

>x<-rnorm(4000,1,0.5) #生成均值為1,標準差為0.5的4000條x序列

>y<-5*x^2+*x^3+ei #設置x、y關系式,生成y

>plot(x,y,cex=.8,col=darkgrey) #畫出x、y的散點圖

>title(Smoothing Spline) #題頭命名“光滑樣條”

>fit=smooth.spline(x,y,df=10) #設置初始自由度為10,進行光滑樣條擬合

>fit2=smooth.spline(x,y,cv=TRUE)#運用交叉驗證法,調整自由度,再次利用光滑樣條擬合

>fit2$df

[1] 9.793812 #讀出調整的自由度值為9.79

>fit2$lambda

[1] 0.004389605##讀出調整的光滑參數為0.0044

>lines(fit2,col=blue,lwd=2) #畫出光滑樣條曲線

>attr(bs(x,df=9.8)knots)#讀出自由度為9.8時,樣條結點的位置

12.5% 25% 7.5% 50% 62.5% 75%100%

0.4115038 0.6566871 0.8289012 0.9908485 1.1491730 1.3207943 2.8121807

運用光滑樣條方法擬合的圖形如圖1所示。

5. 實證分析

本文以R語言中ISLR包中Wage數據集為例,此數據集是以美國中部大西洋地區男員工收入水平為背景的調查數據,通過此數據集運用多項式回歸和光滑樣條回歸兩種方法分析比較研究工資水平和年齡之間的關系。

Figure 1. Smoothing spline

圖1. 光滑線條

#R語言加載ISLR包,調出數據集,進行數據可視化(如圖2)

>library(ISLR)

>attach(Wage)

>plot(age,wage,xlim=agelims,cex=.5,col=darkgrey)

Figure 2. Scatter diagram

圖2. 散點圖

5.1. 多項式回歸模型

由wage-age的散點圖可以,兩變量之間不存在線性關系,最好是擬合一條曲線,使散點均勻地散落在曲線兩側,首先嘗試構造多項式回歸模型。構造多項式回歸模型時,應該在足以解釋自變量和因變量關系的前提下,回歸次數越低越好,避免變量間產生多重共線性。運用交叉驗證法選擇合適的多項式次數。

#運用交叉驗證法選擇多項次最佳回歸次數

>nrow(Wage)

[1] 3000

>train = sample(1:3000,2000)

>cv.err = vector(numeric,5)

>for(i in 1:5){

+ fit = lm(wage~poly(age,i),data = Wage,subset = train)

+ pred = predict(fit,newdata=Wage[-train,])

+ cv.err[i] = mean((pred-wage[-train])^2)

+}

>plot(1:5,cv.err,type=l)

圖3為交叉驗證的結果圖,橫軸表示多項式的次數,縱軸表示均方誤差,從圖中可以看出,多項式次數為4、5時,均方誤差較小,考慮到模型的簡單化,多項式回歸次數得到最佳選擇是4次。

Figure 3. Cross validation diagram

圖3. 交叉驗證圖

#構造4次多項式回歸模型

>fit = lm(wage~poly(age,4),data = Wage)

>agelims = range(age)

>age.grid = seq(from=agelims[1],to=agelims[2])

>preds = predict(fit,newdata = list(age=age.grid),se=T)

>se.bands=cbind(preds$fit+2*preds$se.fit,preds$fit-2*preds$se.fit)# 構建預測值的置信區間

>plot(age,wage,xlim=agelims,cex=0.5,col=darkgrey)

>title(Degree-4 Polynomial,outer = T)

>lines(age.grid,preds$fit,lwd=2,col=blue) # 多項式回歸預測曲線

>matlines(age.grid,se.bands,lwd=2,col = red,lty=3) # 置信區間曲線

圖4中,紅色虛線代表置信區間曲線,由擬合的圖形可以看出,擬合的效果還不錯。

Figure 4. Degree-4 polynomial regression

圖4. 4次多項式回歸

#方差分析

> anova(fit)Analysis of Variance Table

Response: wage

Df Sum Sq Mean Sq F value Pr(>F)

poly(age, 4) 4 450482 112620 70.689 <2.2e?16 ***

Residuals 2995 4771604 1593

---

Signif. codes:0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

方差分析結果中,F值較大,回歸模型顯著,同樣說明4次多項式回歸擬合Wage-age的關系效果良好。

5.2. 光滑樣條回歸

光滑樣條的擬合用R語言的smoothing.spline()實現。和多項式回歸要確定次數一樣,光滑樣條回歸的關鍵是要確定樣條結點的個數、位置和光滑參數。

#光滑樣條回歸

>plot(age,wage,xlim=agelims,cex=.8,col=darkgrey)

>title(Smoothing Spline)

>fit=smooth.spline(age,wage,df=18)#第一次設定自由度為18

>fit2=smooth.spline(age,wage,cv=TRUE)#交叉驗證法調整自由度

>fit2$df

[1] 6.794596#調整后的自由度為6.7946

>lines(fit,col=red,lwd=2)

>lines(fit2,col=blue,lwd=2)

>legend(topright,legend=c(18 DF6.8 DF),col=c(redblue),lty=1,lwd=2,cex=.8)

此段代碼中進行了兩次擬合,使用了兩次smoothing.spline()函數,第一次設定初始自由度df = 18,此時函數值確定自由度df為18對應的光滑參數λ的值。第二次調整,迭代計算,通過交叉驗證選擇合適的值,最終結果就是由λ值選出的自由度為6.7946。

> library(splines, lib.loc=F:/R-3.5.1/library)

> attr(bs(age,df=6.7946)knots)#樣條結點位置

20%40%60% 100%

32 39 46 80

經光滑參數選出df = 6.7946的樣條回歸模型中,R語言給出的樣條結點為32、39、46、80。下圖是分別為自由度為18和6.8的回歸樣條。

圖5中可以看出,自由度為6.8的樣條(藍線)比自由度為18的樣條(紅線)要光滑的多。這就是光滑參數起到了重要作用。

Figure 5. Smoothing spline

圖5. 光滑樣條

5.3. 多項式回歸與光滑樣條比較

圖6中,綠線是4次多項式的擬合曲線,藍線是利用光滑樣條擬合的曲線,在題中,兩者的擬合效果均較好,擬合曲線相差不大,圖中頭尾兩部分可以看出,采用光滑樣條擬合的曲線在首尾兩端較多項式擬合的曲線平穩,結合圖上散點的位置,得出光滑樣條回歸擬合的結果要優于多項式回歸的結果。

Figure 6. Degree-4 polynomial VS smoothing spline

圖6. 4次多項式 VS 光滑樣條

6. 結論

基于美國中部大西洋地區男員工收入水平數據,通過實證分析,可以得出結論如下:

工資會隨著年齡的增長呈先上升后下降的趨勢;人從20歲開始工作到40歲收入是上升的,此階段事業處于上升期,工資水平的增速是緩慢下降的;40歲到60歲工資的水平趨于平穩;60歲之后工資水平緩慢下降。從圖中可以看出60歲到80歲收入水平雖然在下降,依然處于較高的水平,這和美國的退休制度有關系,美國政府沒有硬性規定的統一退休年齡,工薪族往往依照自己的身體和財務狀況作出選擇,可以選擇提前退休或正常退休。基本上在美國65歲被普遍認為是一個正常的退休年齡,圖中在65歲左右的高、中、低的收入分布較40~60歲階段沒有顯著變化。美國實行的是一個“彈性退休制度”,所以在美國經常還有年紀大的老人堅持工作,這一點不同于中國,我國有明確的法定退休年齡,男性60歲,女性50歲。

從回歸擬合的結果角度看,光滑樣條的擬合要優于多項式的擬合。懲罰平方和綜合考慮了曲線擬合的兩個方面:擬合優度和光滑度。同時光滑樣條回歸不必事先對結點進行選擇,克服了以往利用樣條函數進行曲線擬合時存在的缺點,避免了結點選擇的盲目性,既提高擬合程度又在一定程度保證曲線光滑,使得擬合曲線精確美觀。非參數回歸應用廣泛,使用性強,它所需要假定比參數回歸弱的多,適用任何分布的數據,尤其當反應變量與解釋變量間函數關系不清楚時,參數模型難以進行擬合處理,此時非參數估計可作為擬合曲線的一個非常有效的方法。

综合分布图 新浪爱彩网 彩票开奖号码是多少 彩票开奖代码 广东十一选5开奖信息 福彩3d惊人的规律 电脑怎样赚钱 体彩山西11选5开奖查询结果 江西快三在线计划 高频彩吧 中原风釆22选5走势图 3d087开组三还是组六 排列三走势图 二肖中特碼 澳洲幸运10是什么颜色 甘肃十一选五怎么选 吉林快3走势图淘宝1

文章引用:
王鳳雪. 光滑樣條回歸及應用研究[J]. 統計學與應用, 2019, 8(4): 604-612. https://doi.org/10.12677/SA.2019.84069

參考文獻

[1] 陳生長, 徐勇勇, 夏結來. 光滑樣條非參數回歸方法及醫學應用[J]. 中國衛生統計, 1999(6): 23-26.
[2] 陳長生, 伍稚萍, 吳冰. 兒童生長曲線的光滑樣條非參數回歸方法構建[J]. 北華大學學報(自然科學版), 2001, 2(3): 194-196.
[3] 盧一強, 陳中威. 部分線性模型的光滑樣條推斷[J]. 山西大同大學學報(自然科學版), 2010, 26(1): 1-4.
[4] 宋向東. 非參數回歸的罰樣條算法[J]. 燕山大學學報, 2007, 31(3): 263-265.

新浪爱彩网 彩票开奖号码是多少 彩票开奖代码 广东十一选5开奖信息 福彩3d惊人的规律 电脑怎样赚钱 体彩山西11选5开奖查询结果 江西快三在线计划 高频彩吧 中原风釆22选5走势图 3d087开组三还是组六 排列三走势图 二肖中特碼 澳洲幸运10是什么颜色 甘肃十一选五怎么选 吉林快3走势图淘宝1