#清理环境，加载包
rm(list=ls())
library(forecast)
library(tseries)
# as.ts()与is.ts()
Data <- read.table("F:\\桌面\\temp/arima_data.txt", header = TRUE)[, 2]
is.ts(Data)
video1 <- ts(Data)
is.ts(video1)
video2 <- as.ts(Data)
is.ts(video2)



# 示例：时序图

plot.ts(video1, xlab = "时间", ylab = "学习播课视频时长（分钟）")



#平稳性检验
#自相关图
#  acf（）使用示例
acf(video1) # lag.max = 50
tsdisplay(video1) # 查看自先关与偏自相关图


# 单位根检验
library(fUnitRoots)
# unitrootTest(video1)
adf.test(video1)



# 非平稳时间序列差分
difvideo <- diff(video1)
plot.ts(difvideo, xlab = "时间", ylab = "学习播课视频时长（分钟）")
acf(difvideo, lag.max = 30)
adf.test(difvideo)




# 一阶差分后序列的偏自相关图
pacf(difvideo, lag.max = 30)



#  BIC图
library(TSA)
res <- armasubsets(y = difvideo, nar = 5, nma = 5, y.name = 'test', ar.method = 'ols')
plot(res)



#  ARIMA模型
library(forecast)
arima <- Arima(video1, order = c(1, 1, 0))
arima

# auto.arima通过选取AIC和BIC最小来选取模型，与根据acf和pacf图建立的模型进行比较
(mod_auto=auto.arima(video1))

# 残差正态性检验
qqnorm(arima$residuals)
qqline(arima$residuals)

qqnorm(mod_auto$residuals)
qqline(mod_auto$residuals)



# 残差白噪声检验
Box.test(arima$residuals, lag = 5, type = "Ljung-Box")
Box.test(mod_auto$residuals, lag = 5, type = "Ljung-Box")
# ！！！检验结果越大越不好

#模型预测
pre <- forecast(arima, h = 5, level = c(80, 95))
forecast

plot(video1,col='pink')
par(new=T)
plot(pre,col='green')


# 绘制原始值与预测值图形
plot(pre)



#  decompose()函数
video <- ts(Data, start = c(2013, 1), frequency = 12)
video.de <- decompose(video, type = "additive")
video.de
plot(video.de)



#  stl()函数 
video.stl <- stl(video, s.window = "periodic")
video.stl
plot(video.stl)



#  HoltWinters()函数建模
hw.video <- HoltWinters(video, alpha = TRUE, beta = TRUE, gamma = TRUE)
hw.video
plot(hw.video)



# 示例：建模并对模型的残差进行自相关检验与白噪声检验。
library(forecast)
hw.model <- forecast(hw.video, h = 6, level = c(80, 95))
hw.model
hw.model$residuals[which(is.na(hw.model$residuals))] = 0
acf(hw.model$residuals)
Box.test(hw.model$residuals, lag = 10, type = 'Ljung-Box')



# 画出原始值与预测值的图形

plot(hw.model)