I have written the following code to generate 500 data points from a $SARIMA$ model, use $400$ as training data and then predict the following $100$, while estimating the model with AIC. It appeared to me that I did everything correctly as I could see the AIC correctly choosing the model with certain phi values, etc. etc., however my plot outputs for the estimations is very very dense and incomprehensible and I'm not sure why. I checked the number of data points, the window size, etc. and am not sure what I have done wrong in my implementation.
library(CombMSC) library(forecast) ##################### #generate data sdat1<- sarima.Sim(n=500, period=12, model = list(order = c(1,0,0), ar=0.5),list(order= c(1,0,0), ar = 0.5)) #procure training data x.tr <- window(sdat1, start=1, end=400) #candidate models op1 <- arima(x.tr, order=c(0,0,0), list(order= c(0,1,0))) op2 <- arima(x.tr, order=c(1,0,0), list(order= c(1,0,1))) op3 <- arima(x.tr, order=c(1,0,0), list(order= c(0,2,0))) op4 <- arima(x.tr, order=c(0,1,0), list(order= c(0,0,0))) op5 <- arima(x.tr, order=c(1,0,1), list(order= c(0,0,0))) op6 <- arima(x.tr, order=c(1,0,0), list(order= c(1,0,0))) models <- c(op1,op2,op3,op4,op5,op6) models.AIC <- c(op1$aic,op2$aic,op3$aic,op4$aic,op5$aic,op6$aic) mod.best = NULL if (min(models.AIC) == op1$aic){ mod.best=op1 } else if (min(models.AIC) == op2$aic){ mod.best=op2 } else if (min(models.AIC) == op3$aic){ mod.best=op3 } else if (min(models.AIC) == op4$aic){ mod.best=op4 } else if (min(models.AIC) == op5$aic){ mod.best=op5 } else if (min(models.AIC) == op6$aic){ mod.best=op6 } models.AIC mod.best modpred <- predict(mod.best, n.ahead=100) vld.data1<- sdat1[401:500] plot.ts(sdat1, ylim=c(floor(min(sdat1)),ceiling(max(sdat1)))) plot.ts(sdat1, xlim=c(0,400),ylim=c(floor(min(sdat1)),ceiling(max(sdat1)))) lines(modpred$pred, col='blue') lines(modpred$pred-(1.96*modpred$se), col='red') lines(modpred$pred+(1.96*modpred$se), col='red') 
Contents
hide
Best Answer
Seems you are simulating 500 years * 12 observations by years, not 500 observations. You must do:
sarima.Sim(n=500/12, period=12, model=list(order=c(1,0,0), ar=0.5), list(order=c(1,0,0), ar=0.5))