[Rcode]时间序列（整理自R in action）【待更新】

一 生成时序对象

sales<-c(18,33,41,7,34)
tsales<-ts(sales,start=c(2003,1),frequency=12)
###常规操作得到基础信息
plot(tsales)                        #绘制时序图
start(tsales);end(tsales)           #起止点
frequency(tsales)                   
###取子集
tsales.subset<-window(tsales,start=c(2003,5),end=(2003,7));tsales.subset

二 时间序列的平滑化与季节性分解

#转化为对数数据可降低波动性同时可以使用相加模型方便做季节性分解
 plot(AirPassengers)       
 lAirPassengers<-log(AirPassengers)
 plot(lAirPassengers)
#分解时间序列
 fit<-stl(lAirPassengers,s.window = "period")
 plot(fit)
 fit$time.series                                   #查看每个观测值各分解项的值
                                                   #注意这里s.window="period",从而每年的相同月份季节效应相同
#把对数数据转化为原始尺度
 exp(fit$time.series)                              #数据解读:7月的季节性效应为1.24，即乘客数增加了24%

#季节的可视化
 install.packages("forecast")
 library(forecast)
 opar<-par(no.readonly = T)
 par(mfrow=c(2,1))  
 monthplot(AirPassengers,xlab="",ylab="")
 seasonplot(AirPassengers,years.labels="TRUE",main="seasonplot")

三 指数预测模型

 #指数模型预测(单，双，三指数模型)
 #单指数平滑
 library(forecast)
 fit<-ets(nhtemp,model="ANN")
 fit
  #一步向前预测并预测
  forecast(fit,1)              
  plot(forecast(fit,1),xlab="Year",ylab=expression(paste("Temperature(",degree*F,")",)),main="New Haven Annual Mean Temperature")
  accuracy(fit)                  #得到准确性度量
  
 #Holt指数平滑和Holt-Winters指数平滑
 library(forecast)
 fit<-ets(log(AirPassengers),model="AAA")
 fit
 pred<-forecast(fit,5)
 pred   
 #ets()函数与自动预测
 library(forecast)
 fit<-ets(JohnsonJohnson)        #我们没有指定模型时，软件自动搜索了一系列模型并且选择了最小化拟合标准的模型。
 fit    

四 ARIMA预测模型 （保证数据的平稳性）

  #回归前的操作：
   #平稳性：
   #diff(ts,differences=d),forecast的ndiffs()函数可以帮我找到最优的d值
   #保证方差为常数：
   #作对数变换或者Box-Cox变换
 
  #example:Nile数据：
   library(forecast)
   library(tseries)  
   plot(Nile)  
   ndiffs(Nile)   
   dNile<-diff(Nile)   
   plot(dNile)   
   adf.test(dNile)                 #tseries包的adf.test做ADF检验判断序列平稳性
   #选择模型：
   Acf(dNile)
   Pacf(dNile)    
   #拟合模型：
   library(forecast)
   fit<-arima(Nile,order=c(0,1,1))
   fit
   accuracy(fit)
   #模型评价：(主要是残差的独立正态性)
   qqnorm(fit$residuals)
   qqline(fit$residuals)   
   Box.test(fit$residuals,type="Ljung-Box")   
   #预测：
   forecast(fit,3)
   plot(forecast(fit,3),xlab="Year",ylab="Annual Flow")   
  
  #ARIMA的自动预测
   library(forecast)
   fit<-auto.arima(sunspots)
   fit
   forecast(fit,3)   
   accuracy(fit)