0
$\begingroup$

I am running ARIMA model in R and I used auto.arima(X) function to decide appropriate model.After using this function I found that the order of my model is ARIMA(2,1,0). The problem is I run the same ARIMA(2,1,0) model using arima(X,order=c(2,1,0)) and I got AIC as AIC=832.16. but for same model by using auto.arima(X) as AIC=805.29. I dont know why for the same model AIC is different. Please hep me to over come this problem. Thank you in advance.

$\endgroup$
1
$\begingroup$

Please, give further details of what you are doing. Are you using the default options for these functions?. Posting the data (the output from dput(series.name)) will help as well.

In the example below I get the same AIC ($330.97$) for model ARIMA(2,1,0) with forecast::auto.arima and stats::arima:

require(forecast)
set.seed(127)
x <- arima.sim(n=120, model=list(order=c(2,1,0), ar=c(0.7,-0.3)))
forecast::auto.arima(x)
# Series: x 
# ARIMA(2,1,0)                    
# Coefficients:
#          ar1      ar2
#       0.8619  -0.4422
# s.e.  0.0820   0.0816
# sigma^2 estimated as 0.8719:  log likelihood=-162.48
# AIC=330.97   AICc=331.17   BIC=339.33
stats::arima(x, order=c(2,1,0))
# Series: x 
#ARIMA(2,1,0)                    
#Coefficients:
#         ar1      ar2
#      0.8619  -0.4422
#s.e.  0.0820   0.0816
#sigma^2 estimated as 0.8719:  log likelihood=-162.48
#AIC=330.97   AICc=331.17   BIC=339.33
$\endgroup$
1
$\begingroup$

I have a model i use to forecast with ARIMA which i find to be much better than auto arima:

N<-c(#Input your data#)
deltaT<-~set observations per year, 1/4 for quarterly, 1/12 for monthly
horiz<-#set the amount of forecasts you want to do
startY<-c(#,#) #set the start date c(8,1) for 2008 quarter 1 etc
force.log <- "log" # takes values "s" (lets script choose logging) "level" (forces levels) or "log" (forces logs)
Nu <- ts(N, deltat=deltaT, start = startY)
NbMAT <- matrix(0, ncol= (1/deltaT), nrow = (length(Nu)+(1/deltaT)), byrow=TRUE)
for (i in 1:(1/deltaT)) 
{NbMAT[,i] <- c(rep(0, length = i-1), lag(Nu, k=-i+1), rep(0, length = (1/deltaT)-i+1))}
Nbind <- NbMAT[c(-1:(-(1/deltaT)+1), -(length(NbMAT[,1])-1/deltaT+1):-(length(NbMAT[,1]))),]
Nbind2 <- data.frame(Nbind)
lm1 <- lm(X1 ~., data=Nbind2)
lm2 <- lm(X1 ~ X2 + Nbind2[, length(Nbind2[1,])], data=Nbind2)
library(lmtest)
library(car)
bptest1 <- bptest(lm1)
bptest2 <- bptest(lm2) 
gqtest1 <- gqtest(lm1)
if (force.log == "level") {aslog <- "n"} else {{if (force.log == "log") {aslog <- "y"} else {if (bptest1$p.value < 0.1 | bptest2$p.value < 0.1 | gqtest1$p.value < 0.1) {aslog <- "y"} else {aslog <- "n" }}}}
    if(aslog == "y") {N <- log(Nu)} else {N <- Nu}
    plot(N)
    pacf(N)
    cat("log series, y/n? : ", aslog)
    library(tseries)
    library(forecast)
    max.sdiff <- 3 # set maximum seasonal differences allowed. For typical series, 0 is fine.
    arima.force.seasonality <- "n"
    kpssW <- kpss.test(N, null="Level")
    # alt ndiffs(N, alpha = 0.05, test = "kpss") #or test ="adf" / "pp"
    ppW <- tryCatch({ppW <- pp.test(N, alternative = "stationary")},  error = function(ppW) {ppW <- list(error = "TRUE", p.value = 0.99)})
    # if > 0.05, assume unit root. Automated error catching throws in unit root if pp test matrix is singular
    adfW <- adf.test(N, alternative = "stationary", k = trunc((length(N)-1)^(1/3)))
    if(kpssW$p.value < 0.05 | ppW$p.value > 0.05 | adfW$p.value > 0.05) {ndiffsW = 1} else {ndiffsW = 0}
aaw <- auto.arima(N, max.D= max.sdiff, d=ndiffsW, seasonal=TRUE, 
              allowdrift=FALSE, stepwise=FALSE, trace=TRUE, seasonal.test="ch")
orderWA <- c(aaw$arma[1], aaw$arma[6] , aaw$arma[2])
    orderWS <- c(aaw$arma[3], aaw$arma[7] , aaw$arma[4])
if(sum(aaw$arma[1:2])==0) {orderWA[1] <- 1} else {NULL}
    if(arima.force.seasonality == "y") {if(sum(aaw$arma[3:4])==0) {orderWS[1] <- 1} else {NULL}} else {NULL}
    Arimab <- Arima(N, order= orderWA, seasonal=list(order=orderWS), method="ML")
    fArimab <- forecast(Arimab, h=8, simulate= TRUE, fan=TRUE)
    if (aslog == "y") {fArimabF <- exp(fArimab$mean[1:horiz])} else {fArimabF <-fArimab$mean[1:horiz]}
    # View(fArimabF, title = "Arima forecast") # gives the forecast
    plot(fArimab, main = "ARIMA Forecast", sub="blue=fitted, red=actual") # ylim=c(17, 20)
    lines(N, col="red", lwd=2)
    lines(ts(append(fitted(Arimab), fArimab$mean[1]), deltat=deltaT, start = startY), 
  col= "blue", lwd = 2) # makes the graph look nicer
if (aslog == "y") {Arimab2f <- exp(fArimab$mean[1:horiz])} else {Arimab2f <- fArimab$mean[1:horiz]}
    start(fArimab$mean) -> startARIMA
ArimaALTf <- ts(prettyNum(Arimab2f, big.interval = 3L, big.mark = ","), deltat = deltaT , start= startARIMA)
View(ArimaALTf, title = "ARIMA forecast")
summary(Arimab)
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.