ARMA Model

An ARMA Model (Autoregressive Moving Average) is a time-series model that represents a stationary sequence using past values (AR part) and past forecast errors (MA part). It is important for capturing temporal dependence, enabling forecasting and uncertainty estimation in ML pipelines. Example: modeling hourly demand as ARMA to predict the next hour’s load and quantify prediction error.

Imagine tracking the temperature each day. Today’s temperature is usually related to recent days, plus a bit of random “weather noise.” An ARMA model is a way to describe and predict this kind of time-ordered data by combining two simple ideas: autoregressive (AR), meaning the value now depends on a few past values, and moving average (MA), meaning it also depends on recent random shocks or errors. In statistics and machine learning, ARMA models are used for forecasting things like sales, demand, or sensor readings when the past helps explain the present.

import numpy as np  # numerical arrays for time series values
import pandas as pd  # time-indexed Series container
from statsmodels.tsa.arima.model import ARIMA  # ARMA is ARIMA with d=0
from sklearn.metrics import mean_absolute_error  # simple forecast accuracy metric

# Create a realistic monthly time series (e.g., demand index) with trend + seasonality + noise
rng = np.random.default_rng(7)  # reproducible randomness
idx = pd.date_range("2018-01-01", periods=96, freq="MS")  # 8 years of monthly starts
t = np.arange(len(idx))  # time step index
y = 120 + 0.15 * t + 6 * np.sin(2 * np.pi * t / 12) + rng.normal(0, 2.5, len(t))  # data-generating process
ts = pd.Series(y, index=idx, name="demand_index")  # store as a labeled time series

# Split into train/test to evaluate out-of-sample forecasting
train, test = ts.iloc[:-12], ts.iloc[-12:]  # last year as holdout

# Fit an ARMA(p, q) by using ARIMA(p, 0, q) (no differencing)
# *** KEY LINE ***
arma = ARIMA(train, order=(2, 0, 1), trend="c").fit()  # ARMA(2,1) with an intercept

# Produce forecasts for the test horizon
pred = arma.get_forecast(steps=len(test))  # forecast object with mean + intervals
yhat = pred.predicted_mean  # point forecasts aligned to future dates
ci = pred.conf_int(alpha=0.05)  # 95% prediction intervals

# Quick, practical evaluation for ML workflows
mae = mean_absolute_error(test, yhat)  # average absolute forecast error
print(arma.summary().tables[1])  # parameter estimates for AR and MA terms
print(f"Test MAE: {mae:.3f}")  # scalar forecast accuracy
print(pd.concat([test.rename("actual"), yhat.rename("forecast"), ci], axis=1).head())  # inspect outputs

• The ARIMA call with order=(2, 0, 1) is the ARMA fit; d=0 means no differencing is applied.
• get_forecast returns both point forecasts and uncertainty; conf_int provides prediction intervals useful for risk-aware decisions.
• Common mistake: fitting ARMA to a series with strong trend/nonstationarity; in practice you’d often difference (ARIMA) or detrend before fitting.
• In an AI/ML pipeline, the forecast (and optionally interval width) can be used as features or as a baseline model to beat with more complex approaches.

set.seed(7)  # make randomness reproducible

# Create a realistic monthly time series (e.g., demand index) with trend + seasonality + noise
n <- 96  # 8 years of monthly data
t <- 0:(n - 1)  # time step index
y <- 120 + 0.15 * t + 6 * sin(2 * pi * t / 12) + rnorm(n, mean = 0, sd = 2.5)  # data-generating process
ts_y <- ts(y, start = c(2018, 1), frequency = 12)  # store as a monthly ts object

# Split into train/test (last 12 months as holdout)
train <- window(ts_y, end = c(2024, 12))  # first 84 points (2018-01 to 2024-12)
test  <- window(ts_y, start = c(2025, 1))  # last 12 points (2025-01 to 2025-12)

# Fit an ARMA(p, q) using arima() with no differencing (d = 0)
# *** KEY LINE ***
fit <- arima(train, order = c(2, 0, 1), include.mean = TRUE)  # ARMA(2,1) with intercept

# Forecast the next 12 months and extract uncertainty
fc <- predict(fit, n.ahead = length(test))  # returns mean forecasts and standard errors
yhat <- ts(fc$pred, start = c(2025, 1), frequency = 12)  # align forecasts to test period
se <- fc$se  # forecast standard errors for intervals
lower <- yhat - 1.96 * se  # 95% lower bound (approx.)
upper <- yhat + 1.96 * se  # 95% upper bound (approx.)

# Simple out-of-sample evaluation
mae <- mean(abs(as.numeric(test) - as.numeric(yhat)))  # average absolute forecast error
print(fit)  # parameter estimates for AR and MA terms
cat("Test MAE:", round(mae, 3), "\n")  # scalar forecast accuracy
print(cbind(actual = as.numeric(test), forecast = as.numeric(yhat),
            lower95 = as.numeric(lower), upper95 = as.numeric(upper))[1:5, ])  # inspect outputs

• The arima call with order = c(2, 0, 1) is the ARMA fit; setting the middle term to 0 avoids differencing.
• predict provides forecast means and standard errors; the interval here uses a normal approximation (± 1.96 * se).
• Common mistake: interpreting se as observation noise only; it reflects increasing forecast uncertainty with horizon under the fitted model.
• In ML workflows, the MAE on a holdout period is a quick baseline metric before trying richer models (e.g., SARIMA, state space, or ML regressors with lag features).