‘USA Total GDP’ time series data

Forecasting using ARIMA & ETS models

Author

Tom Maccor

The objective is to model time-series data, select the model that best fits, and forecast 5 years ahead.

In this instance, we will use a dataset that contains the Total GDP for the USA per year, starting in 1960 and finishing in 2017:

Rows: 58
Columns: 9
Key: Country [1]
$ Country    <fct> "United States", "United States", "United States", "United …
$ Code       <fct> USA, USA, USA, USA, USA, USA, USA, USA, USA, USA, USA, USA,…
$ Year       <dbl> 1960, 1961, 1962, 1963, 1964, 1965, 1966, 1967, 1968, 1969,…
$ GDP        <dbl> 543300, 563300, 605100, 638600, 685800, 743700, 815000, 861…
$ Growth     <dbl> NA, 2.3000000, 6.1000000, 4.4000000, 5.8000000, 6.4000000, …
$ CPI        <dbl> 13.56306, 13.70828, 13.87261, 14.04459, 14.22421, 14.44969,…
$ Imports    <dbl> 4.196576, 4.029824, 4.131549, 4.087065, 4.097404, 4.235579,…
$ Exports    <dbl> 4.969630, 4.899698, 4.809122, 4.870028, 5.103529, 4.988571,…
$ Population <dbl> 180671000, 183691000, 186538000, 189242000, 191889000, 1943…

This time series displays a trend and no seasonality; thus, it is non-stationary data. In order to be able to model it using an ARIMA model, we need to make this time-series stationary.

Understanding the ARIMA models

An autoregressive integrated moving average (ARIMA) model is a form of regression analysis that gauges the strength of one dependent variable relative to other changing variables. The model’s goal is to predict future securities or financial market moves by examining the differences between values in the series instead of through actual values.

An ARIMA model can be understood by outlining each of its components as follows:

Autoregression (AR): refers to a model that shows a changing variable that regresses on its own lagged, or prior, values.
Integrated (I): represents the differencing of raw observations to allow the time series to become stationary (i.e., data values are replaced by the difference between the data values and the previous values).
Moving average (MA): incorporates the dependency between an observation and a residual error from a moving average model applied to lagged observations.

Time series analysis & data preparation

First we have a look at the variance, and we can see it does not change (increase or decrease) over time. Therefore, there is no need to perform any transformations to the data (such as a Box-Cox logarithm being applied).

We proceed to difference the data once, and plot to see what is the result:

As we can see, the differenced time-series is still non-stationary -as evidenced by the pattern of the ACF plot (which is not decreasing with time as expected).

If we do a 2nd difference, the time-series now has the same appearance:

BUT…. the ACF and PACF plots now tell us that the time-series is now stationary -evidenced by the

USA_GDP %>%  gg_tsdisplay(difference(difference(Exports)), plot_type='partial')

MODEL FITTING

We can now proceed to fit several ARIMA models to our time-series:

USA_GDP_fit <- USA_GDP %>%
  model(arima012 = ARIMA(GDP ~ pdq(0,1,2)),
        stepwise = ARIMA(GDP),
        search = ARIMA(GDP, stepwise=FALSE),
        )

USA_GDP_fit

# A mable: 1 x 4
# Key:     Country [1]
  Country                      arima012       stepwise         search
  <fct>                         <model>        <model>        <model>
1 United States <ARIMA(0,1,2) w/ drift> <ARIMA(0,2,2)> <ARIMA(0,2,2)>

glance(USA_GDP_fit) |> arrange(AICc) |> select(.model:BIC)

# A tibble: 3 × 6
  .model         sigma2 log_lik   AIC  AICc   BIC
  <chr>           <dbl>   <dbl> <dbl> <dbl> <dbl>
1 stepwise 26149937130.   -750. 1507. 1507. 1513.
2 search   26149937130.   -750. 1507. 1507. 1513.
3 arima012 31919699906.   -769. 1546. 1547. 1554.

Upon inspection of the models’ performance (by comparing their AICc score -Aikaiken’s Information Criterion), the model that seems to perform the best is the ‘search’ model.

We then proceed to check the residuals of this model:

USA_GDP_fit |>
  select(search) |>
  gg_tsresiduals()

The ACF plot of the residuals from the ARIMA ‘search’ model (0,2,2) model shows that all autocorrelations are within the threshold limits, indicating that the residuals are behaving like white noise.

FORECASTING

So we can now proceed to produce a forecast for the next 5 years:

### Forecasts  ###
USA_GDP_fit %>%
  forecast(h=5) |>
  filter(.model=='search') |>
  autoplot(USA_GDP %>% filter(Year >= 2000)) +
  labs(title = "USA GDP",  y = "Year")

ARIMA vs. ETS model

We can compare our best ARIMA model, with an ETS (Error, Trend and Seasonality, also known as exponential smoothing) model.

Below is the forecast that we obtain:

### Compare results with using a model with ETS transformation
USA_GDP |>
  slice(-n()) |>
  stretch_tsibble(.init = 10) |>
  model(
    ETS(GDP),
    ARIMA(GDP)
  ) |>
  forecast(h = 1) |>
  accuracy(USA_GDP) |>
  select(.model, RMSE:MAPE)

# A tibble: 2 × 5
  .model        RMSE     MAE   MPE  MAPE
  <chr>        <dbl>   <dbl> <dbl> <dbl>
1 ARIMA(GDP) 222698. 136765. 0.493  1.76
2 ETS(GDP)   190831. 126608. 0.589  1.73

model_comparison <- USA_GDP %>%
  model (ETS = ETS(GDP),
        stepwise = ARIMA(GDP),
        search = ARIMA(GDP, stepwise=FALSE),
        )

glance(model_comparison) |> arrange(AICc) |> select(.model:BIC)

# A tibble: 3 × 6
  .model     sigma2 log_lik   AIC  AICc   BIC
  <chr>       <dbl>   <dbl> <dbl> <dbl> <dbl>
1 stepwise 2.61e+10   -750. 1507. 1507. 1513.
2 search   2.61e+10   -750. 1507. 1507. 1513.
3 ETS      6.78e- 4   -789. 1588. 1589. 1598.

USA_GDP |>
  model(ETS(GDP)) |>
  forecast(h = "5 years") |>
  autoplot(USA_GDP |> filter(Year >= 2000)) +
  labs(title = "USA GDP",  y = "Year")

If we take the AICc as a measure of the model’s fit, ARIMA has a better performance that this ETS model, even though the RMSE (root square Mean Error) is lower for the ETS model.

AICc and RMSE are inter-related but they represent different objectives in choosing the best model.

RMSE/MAPE are measures of error and disregards the “complexity” of the model. Optimizing for RMSE/MAPE can give you accurate results, but could lead to overly complex model that captures too much noise in the data, otherwise known as overfitting.

This is where AIC/AICc and their relative BIC comes in. They take the error term and add a penalty related to the number of predictors used in the model, such that more complex models are less favorable and allows you to strike a balance between a complex but accurate model, vs a simpler but still reasonably accurate model.

It ultimately comes down to the purpose of your model. If having the most accurate prediction matters then you might simply look at RMSE/MAPE -but if you need a model that is more interpretable/explainable, then you might want to consider AICc which better balances complexity and accuracy.

As we can see when we compare both forecasts in the plots below, the prediction intervals in the ARIMA model have a lot less variance.