phaser phaser - 1 year ago 136
R Question

ARIMA modelling, prediction and plotting with CO2 dataset in R

I am working with

. I would like to plot
model over my data. I have tried
with no success.

Here is my code:

###### Time Series

# format: time series

# format: matrix
dmn <- list(, unique(floor(time(co2))))
co2.m <- matrix(co2, 12, dimnames = dmn)

co2.dt <- pracma::detrend(co2.m, tt = 'linear')
co2.dt <- ts(as.numeric(co2.dt), start = c(1959,1), frequency=12)

# first diff
co2.dt.dif <- diff(co2.dt,lag = 12)

# Second diff
co2.dt.dif2 <- diff(co2.dt.dif,lag = 1)

With the data prepared, I ran the following

results <- arima0(co2.dt.dif2, order = c(2,0,0), method = "ML")
resultspredict <- predict(results, n.ahead = 36)

I would like to plot the model and the prediction. I am hoping there is a way to do this in base R. I would also like to be able to plot the predictions as well.

Answer Source

Session 2: A better idea for time series modelling

In previous session, we have seen the difficulty in prediction if we separate de-trending and modelling of de-trended series. Now, we try to combine those two stages in one go.

The seasonal pattern of co2 is really strong, so we need a seasonal differencing anyway:

co2dt <- diff(co2, lag = 12)
par(mfrow = c(1,2)); ts.plot(co2dt); acf(co2dt)

enter image description here

After this seasonal differencing, co2dt does not look stationary. So we need further a non-seasonal differencing.

co2dt.dif <- diff(co2dt)
par(mfrow = c(1,2)); ts.plot(co2dt.dif); acf(co2dt.dif)

enter image description here

The negative spikes within season and between season suggest that a MA process is needed for both. I will not work with co2dt.dif; we can work with co2 directly:

fit <- arima0(co2, order = c(0,1,1), seasonal = c(0,1,1))

enter image description here

Now the residuals are perfectly uncorrelated! So we have an ARIMA(0,1,1)(0,1,1)[12] model for co2 series.

As usual, fitted values are obtained by subtracting residuals from data:

co2fitted <- co2 - fit$residuals

Predictions are made by a single call to predict:

co2pred <- predict(fit, n.ahead = 36, = FALSE)

Let's plot them together:

ts.plot(co2, co2fitted, co2pred, col = 1:3)

enter image description here

Oh, this is just gorgeous!

Recommended from our users: Dynamic Network Monitoring from WhatsUp Gold from IPSwitch. Free Download