R語言時間序列應用(decompose、Holt-Winters初步)


對於明顯的周期性時間序列,可以使用decompose函數對數據進行分解成季節部分、趨勢部分、隨機部分三種。decompose函數有兩種type,即“additive”以及“multiplicative”兩種,還有一個fliter選項,表示是否加入線性濾波,一般fliter選擇NULL即可。下面的例子展現了使用decompose分析含有季節因素時間序列數據的例子

將某地區1962-1970年平均每頭奶牛的月度產奶量數據導入outcome內。對於時間序列數據,常常還要使用ts函數將其形式進行轉換成時間序列專用的數據形式:

outcome<-ts(B,frequency = 12,start = c(1962,1))

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Result1<-decompose(outcome,type = "additive")

plot(Result1)

 

Result2<--decompose(outcome,type = "multiplicative")

plot(Result2)

 

 

選擇乘性模型繼續進行分析。通過以下命令提取模型中的趨勢項,並使用線性模型進行擬合。

results<-Result1

plot.ts(results$trend)

 

abline(lm(results$trend~time(outcome)),col="red")

 

reg<-lm(results$trend~time(outcome))

summary(reg)

par(mfrow=c(2,2))

plot(reg,which=c(1:4))

 模型回歸的摘要,以及統計圖如下

Ver

Coefficients:

                Estimate Std. Error t value Pr(>|t|)    

(Intercept)   -4.269e+04  4.678e+02  -91.25   <2e-16 ***

time(outcome)  2.207e+01  2.379e-01   92.75   <2e-16 ***

---

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

 

Residual standard error: 5.383 on 94 degrees of freedom

  (12 observations deleted due to missingness)

Multiple R-squared:  0.9892, Adjusted R-squared:  0.9891

F-statistic:  8603 on 1 and 94 DF,  p-value: < 2.2e-16

 

顯然回歸模型擬合效果一般,殘差出現了明顯的自相關性。使用該方法對模型進行預測需要提取出序列的季節因子,並與線性擬合的趨勢項預測數據相乘。首先我們提取出季節指數,得到季節指數sea,注意要先把results$seasonal項化為向量型,以便后面的計算。

注意如果前面分解采用加性模型,應該采用原始數據減去季節因素。

sea<-as.vector(results$seasonal)

sea<-sea[1:12]

par(mfrow=c(1,1))

plot(sea,type = "l")

pre<-1971+(0:11)/12

pre<--42690.59+22.07*pre

pred<-pre*sea

 

 最后得到pred項即為預測的結果

> sea [1] 0.9752900 0.9239188 1.0481092 1.0718562 1.1567964 1.1178806

 [7] 1.0416054 0.9783800 0.9233615 0.9286178 0.8919580 0.9422262

> pred [1] 789.3802 749.5006 852.1739 873.4530 944.7980 915.0700

 [7] 854.5487 804.4770 760.9361 766.9756 738.3376 781.6811

 

 含季節因素時間序列預測的另一個方式是采用x11方法,與前面decompose不同的地方在於提取季節因子,x11方法提取季節因子是通過某季節該變量的平均數/所有季節該變量的平均數得到。使用x11方法得到季節指數sea2

xbar<-rep(0,12)

for (i in 1:108){

  xbar[i%%12+1]<-xbar[i%%12+1]+B[i]

}

xbar<-xbar/9

ex<-sum(B)/108

sea2<-xbar/ex

 

> sea2 [1] 0.9551790 0.9607222 0.9125752 1.0381691 1.0643016 1.1536269

 [7] 1.1165664 1.0429205 0.9841622 0.9309471 0.9385493 0.9022806

 

從模型中去除季節因素,並對趨勢項使用線性擬合

trend<-c(0)

for (i in 1:108){

  trend[i]<-B[i]/sea2[i%%12+1]

}

trend<-ts(trend,frequency = 12,start = c(1962,1))

ts.plot(trend)

reg2<-lm(trend~time(trend))
abline(lm(trend~time(trend)),col="red")

 

顯然x11方法提取的趨勢項沒有decompose函數提取的平滑,但是其趨勢擬合效果更好,殘者基本符合正態分布

Call:

lm(formula = trend ~ time(trend))

 

Coefficients:

(Intercept)  time(trend)  

  -41832.05        21.63  

 

 

時間序列預測的有效方法為Holt-Winters方法平滑預測,該方法有三個參數來控制:alpha,beta和gamma,分別對應當前時間點上的水平,趨勢部分的斜率和季節性部分。參數alpha,beta和gamma的取值都在0和1之間,並且當其取值越接近0意味着對未來的預測值而言最近的觀測值占據相對較小的權重。如果我們事先采用decompose對季節因素進行提取

將某城市1980年1月至1995年8月每月屠宰生豬數量導入data集中,並用ts函數轉化

datats<-ts(data,frequency = 12,start =c(1980,1) )

ts.plot(datats)

datatshw<-HoltWinters(datats)

 

#需要使用forecast包中forecast.HoltWinters以及plot.forecast函數

datatshwforecast<-forecast.HoltWinters(datatshw,h=48)
plot.forecast(datatshwforecast)

 

  datatshw中儲存了平滑參數的值如下:

Smoothing parameters:

 alpha: 0.3074266

 beta : 0.00366508

 gamma: 0.4374631

 

參考資料:http://blog.csdn.net/jiabiao1602/article/details/43153139

              http://www.dataguru.cn/article-3235-1.html

              http://bbs.pinggu.org/thread-4116352-1-1.html

 


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM