#加载程序包
library(tseries)
library(forecast)


#加的可能是错的
library(Metrics)
#加的可能是错的



#读取数据
data <-read.table("yichang.txt",header=FALSE)
data <- as.matrix(data)
#画出data变量图
plot(data)

n<-length(data)
m<-floor(n*0.8)
L<-n-m

#将数据拆分为训练集
#取data的80%进行训练

train<-data[1:m,1]

#tsdisplay用来看ACF和PACF

tsdisplay(train)
auto.arima(train)

#diff():求差和求导

###  differential   ###
s1<-diff(train,12)
adf.test(s1)
tsdisplay(s1)

ndiffs(train)
#######################

###  forecast ###
#nrow:行 ncol：列 matrix:矩阵
fore<- matrix(NA, nrow=L, ncol=1)
#循环拟合，for(var in seq),var:变量名，seq：向量表达式
for (i in 1:L) {

fit1<-arima(data[1:(m+i-1),1],order=c(1,0,1))
tsdiag(fit1)
f.p1<-forecast(fit1,h=3,level=c(99.5))
fore<-Vectorize(fore)
fore[i,1]<-f.p1$mean
}


f.p1$mean
# 7/25 #

plot(fore)

lines(fore,col="green")
lines(data[m+1:n],col="red")

write.table(fore,"fore_mean_h=3.csv",sep=",")
obs<-data[m+1:n]
write.table(obs,"obs_h=3.csv",sep=",")