3  Modelos ARIMA

3.1 Introducción

En este capítulo se introduce la familia de modelos ARIMA para series de tiempo. Se motiva su uso cuando la regresión clásica no captura toda la dinámica (p. ej., autocorrelación en residuos). Se presentan los modelos autorregresivos (AR), de media móvil (MA) y mixtos (ARMA), así como nociones de causalidad, invertibilidad y problemas de redundancia de parámetros. Se utilizan notación de operador de rezago y resultados sobre soluciones estacionarias, junto con ejemplos y código en R.

3.2 Modelos Autorregresivos de Media Móvil

3.2.1 Motivación y contexto

La regresión clásica es estática (solo entradas contemporáneas). En series de tiempo se desea permitir dependencia respecto a valores pasados de la propia serie y de entradas externas, lo cual abre la puerta a la predicción.

3.2.2 Introducción a modelos autorregresivos (AR)

3.2.2.1 Ejemplo motivador

Modelo generador:

\[ x_t = x_{t-1} - 0.90\,x_{t-2} + w_t,\quad w_t \sim \text{WN}(0,\sigma_w^2=1). \]

Un pronóstico ingenuo:

\[ x_{n+1}^{\,n} = x_n - 0.90\,x_{n-1}. \]

3.2.2.2 Definición (Modelo \(\mathbf{AR}(p)\))

Un modelo autorregresivo de orden \(p\) es

\[ \begin{equation*} x_t=\phi_1 x_{t-1}+\phi_2 x_{t-2}+\cdots+\phi_p x_{t-p}+w_t, \end{equation*} \]

donde \(x_t\) es estacionaria, \(w_t\sim wn(0,\sigma_w^2)\) y \(\phi_p\neq 0\). Si \(\mu=\mathrm{E}(x_t)\neq 0\), equivalentemente

\[ \begin{equation*} x_t=\alpha+\phi_1 x_{t-1}+\cdots+\phi_p x_{t-p}+w_t,\quad \alpha=\mu(1-\phi_1-\cdots-\phi_p). \end{equation*} \]

3.2.2.3 Notación de operador de rezago

Con el operador \(B\), el modelo se escribe

\[ \begin{equation*} (1-\phi_1 B-\cdots-\phi_p B^p)\,x_t=w_t, \end{equation*} \]

o concisamente

\[ \begin{equation*} \phi(B)\,x_t=w_t, \end{equation*} \]

donde

\[ \begin{equation*} \phi(B)=1-\phi_1 B-\cdots-\phi_p B^p. \end{equation*} \]

3.2.2.4 Ejemplo (AR(1))

Para \(x_t=\phi x_{t-1}+w_t\) e iterando hacia atrás:

\[ x_t=\phi^k x_{t-k}+\sum_{j=0}^{k-1}\phi^j w_{t-j}. \]

Si \(|\phi|<1\) y \(\sup_t \mathrm{var}(x_t)<\infty\), la solución estacionaria como proceso lineal es

\[ \begin{equation*} x_t=\sum_{j=0}^{\infty}\phi^j w_{t-j}. \end{equation*} \]

Con esto, \(\mathrm{E}(x_t)=0\) y la autocovarianza

\[ \begin{align*} \gamma(h) &=\mathrm{cov}(x_{t+h},x_t) =\sigma_w^2 \sum_{j=0}^{\infty}\phi^{h+j}\phi^j =\frac{\sigma_w^2\,\phi^{h}}{1-\phi^2},\quad h\ge 0. \end{align*} \]

Así, la ACF es

\[ \begin{equation*} \rho(h)=\frac{\gamma(h)}{\gamma(0)}=\phi^{h},\quad h\ge 0, \end{equation*} \]

y satisface

\[ \begin{equation*} \rho(h)=\phi\,\rho(h-1),\quad h=1,2,\ldots \end{equation*} \]

3.2.2.5 Ejemplo (Trayectorias simuladas AR(1))

Comparación de \(\phi=0.9\) (trayectoria suave) vs \(\phi=-0.9\) (trayectoria “rugosa”):

par(mfrow=c(2,1))
plot(arima.sim(list(order=c(1,0,0), ar=.9), n=100), ylab="x",
     main=(expression(AR(1)~~~phi==+.9)))
plot(arima.sim(list(order=c(1,0,0), ar=-.9), n=100), ylab="x",
     main=(expression(AR(1)~~~phi==-.9)))

3.2.2.6 Ejemplo (Explosión y causalidad)

Para \(|\phi|>1\), se obtiene al iterar hacia adelante:

\[ \begin{align*} x_t&=\phi^{-k}x_{t+k}-\sum_{j=1}^{k-1}\phi^{-j}w_{t+j} \end{align*} \]

y, en el límite,

\[ \begin{equation*} x_t=-\sum_{j=1}^{\infty}\phi^{-j}w_{t+j}, \end{equation*} \]

que es estacionaria pero dependiente del futuro (no causal) y por tanto inútil para pronóstico.

3.2.2.7 Ejemplo (Toda explosión tiene un “causal” equivalente)

Si \(x_t=\phi x_{t-1}+w_t\) con \(|\phi|>1\), existe un proceso causal equivalente

\[ y_t=\phi^{-1}y_{t-1}+v_t,\quad v_t\sim \text{iid }N\!\left(0,\sigma_w^2\phi^{-2}\right), \]

que es estocásticamente igual a \(x_t\) (misma familia finito-dimensional), ya que ambas comparten la misma estructura de covarianza:

\[\begin{aligned} \gamma_x(h) & =\operatorname{cov}\left(x_{t+h}, x_t\right)=\operatorname{cov}\left(-\sum_{j=1}^{\infty} \phi^{-j} w_{t+h+j},-\sum_{k=1}^{\infty} \phi^{-k} w_{t+k}\right) \\ & =\sigma_w^2 \phi^{-2} \phi^{-h} /\left(1-\phi^{-2}\right) \end{aligned}\]

3.2.2.8 Coincidencia de coeficientes y operadores

Si \(\phi(B)x_t=w_t\) y \(x_t=\psi(B)w_t=\sum_{j=0}^\infty \psi_j B^j w_t\), entonces

\[ \phi(B)\psi(B)w_t=w_t \;\Rightarrow\; (1-\phi B)\left(1+\psi_1 B+\psi_2 B^2+\cdots\right)=1. \]

Igualando coeficientes se obtiene \(\psi_j=\phi^j\). En general, tratar \(B\) como \(z\) complejo:

\[ \phi^{-1}(z)=\frac{1}{1-\phi z}=\sum_{j=0}^{\infty}\phi^j z^j,\quad |z|\le 1. \]

3.2.3 Introducción a modelos de media móvil (MA)

3.2.3.1 Definición (Modelo \(\mathbf{MA}(q)\))

\[ \begin{equation*} x_t=w_t+\theta_1 w_{t-1}+\cdots+\theta_q w_{t-q}, \end{equation*} \]

con \(w_t\sim wn(0,\sigma_w^2)\) y \(\theta_q\neq 0\). Equivalentemente,

\[ \begin{equation*} x_t=\theta(B)w_t, \end{equation*} \]

donde el operador MA es

\[ \begin{equation*} \theta(B)=1+\theta_1 B+\cdots+\theta_q B^q. \end{equation*} \]

Los procesos MA son estacionarios para cualquier \((\theta_1,\ldots,\theta_q)\).

3.2.3.2 Ejemplo (MA(1))

Para \(x_t=w_t+\theta w_{t-1}\), \(\mathrm{E}(x_t)=0\) y

\[ \gamma(h)= \begin{cases} (1+\theta^2)\sigma_w^2, & h=0,\\ \theta\,\sigma_w^2, & h=1,\\ 0, & h>1, \end{cases} \quad \rho(h)= \begin{cases} \dfrac{\theta}{1+\theta^2}, & h=1,\\ 0, & h>1. \end{cases} \]

par(mfrow = c(2,1))
plot(arima.sim(list(order=c(0,0,1), ma=.9), n=100), ylab="x",
     main=(expression(MA(1)~~~theta==+.5)))
plot(arima.sim(list(order=c(0,0,1), ma=-.9), n=100), ylab="x",
     main=(expression(MA(1)~~~theta==-.5)))

3.2.3.3 Ejemplo (No unicidad e invertibilidad)

Para MA(1), \(\rho(h)\) es igual para \(\theta\) y \(1/\theta\); distintos pares \((\sigma_w^2,\theta)\) pueden inducir la misma \(\gamma(h)\): por ejemplo tomando \(\theta=5\) y \(\sigma^2_w=1\) o \(\theta=1/5\) y \(\sigma^2_w=25\). Se elige la representación invertible: aquella con \(|\theta|<1\), que permite

\[ w_t=\sum_{j=0}^{\infty}(-\theta)^j x_{t-j}. \]

Vía operadores: si \(\theta(B)=1+\theta B\), entonces para \(|\theta|<1\),

\[ \pi(B)=\theta^{-1}(B)=\sum_{j=0}^{\infty}(-\theta)^j B^j. \]

3.2.4 Modelos ARMA

3.2.4.1 Definición (Proceso \(\mathbf{ARMA}(p,q)\))

Una serie \({x_t}\) es ARMA\((p,q)\) si es estacionaria y

\[ \begin{equation*} x_t=\phi_1 x_{t-1}+\cdots+\phi_p x_{t-p}+w_t+\theta_1 w_{t-1}+\cdots+\theta_q w_{t-q}, \end{equation*} \]

con \(\phi_p\neq 0\), \(\theta_q\neq 0\), \(\sigma_w^2>0\). Si \(\mu\neq 0\),

\[ \begin{equation*} x_t=\alpha+\phi_1 x_{t-1}+\cdots+\phi_p x_{t-p}+w_t+\theta_1 w_{t-1}+\cdots+\theta_q w_{t-q}, \end{equation*} \]

y, en forma concisa,

\[ \begin{equation*} \phi(B)\,x_t=\theta(B)\,w_t. \end{equation*} \]

3.2.4.2 Ejemplo (Redundancia de parámetros)

Partiendo de \(x_t=w_t\) y multiplicando por \(\eta(B)=1-0.5B\):

\[ \begin{equation*} x_t=0.5\,x_{t-1}-0.5\,w_{t-1}+w_t, \end{equation*} \]

que parece un ARMA(1,1) pero es ruido blanco: hay sobre–parametrización.

Demostración empírica de sobre–parametrización al ajustar ARMA(1,1) a datos iid:

set.seed(8675309) 
x = rnorm(150, mean=5) # iid N(5,1)
arima(x, order=c(1,0,1)) # estimación

Call:
arima(x = x, order = c(1, 0, 1))

Coefficients:
          ar1     ma1  intercept
      -0.9595  0.9527     5.0462
s.e.   0.1688  0.1750     0.0727

sigma^2 estimated as 0.7986:  log likelihood = -195.98,  aic = 399.96

# Coeficientes típicos: ar1 ~ -0.96, ma1 ~ 0.95, ambos "significativos"

3.2.5 Polinomios AR y MA

3.2.5.1 Definición (Polinomios)

\[ \begin{equation*} \phi(z)=1-\phi_1 z-\cdots-\phi_p z^p,\quad \phi_p\neq 0; \end{equation*} \]

\[ \begin{equation*} \theta(z)=1+\theta_1 z+\cdots+\theta_q z^q,\quad \theta_q\neq 0. \end{equation*} \]

Forma mínima: \(\phi(z)\) y \(\theta(z)\) sin factores comunes (evita redundancia).

3.2.5.2 Definición (Causalidad)

ARMA\((p,q)\) es causal si puede escribirse como proceso lineal unilateral

\[ \begin{equation*} x_t=\sum_{j=0}^{\infty}\psi_j w_{t-j}=\psi(B)w_t,\quad \sum_{j=0}^{\infty}|\psi_j|<\infty,\ \psi_0=1. \end{equation*} \]

3.2.5.3 Propiedad (Causalidad)

ARMA\((p,q)\) es causal si y solo si \(\phi(z)\neq 0\) para \(|z|\le 1\). Los pesos \(\psi_j\) satisfacen

\[ \psi(z)=\sum_{j=0}^{\infty}\psi_j z^j=\frac{\theta(z)}{\phi(z)},\quad |z|\le 1. \]

3.2.5.4 Definición (Invertibilidad)

ARMA\((p,q)\) es invertible si existe

\[ \begin{equation*} \pi(B)x_t=\sum_{j=0}^{\infty}\pi_j x_{t-j}=w_t,\quad \sum_{j=0}^{\infty}|\pi_j|<\infty,\ \pi_0=1. \end{equation*} \]

3.2.5.5 Propiedad 3.2 (Invertibilidad)

ARMA\((p,q)\) es invertible si y solo si \(\theta(z)\neq 0\) para \(|z|\le 1\), y

\[ \pi(z)=\sum_{j=0}^{\infty}\pi_j z^j=\frac{\phi(z)}{\theta(z)},\quad |z|\le 1. \]

3.2.6 Ejemplos de aplicación de las condiciones

3.2.6.1 Ejemplo (Redundancia, causalidad, invertibilidad)

Dado

\[ x_t=0.4\,x_{t-1}+0.45\,x_{t-2}+w_t+w_{t-1}+0.25\,w_{t-2}, \]

u operador:

\[ (1-0.4B-0.45B^2)\,x_t=(1+B+0.25B^2)\,w_t. \]

Factorizando:

\[ \phi(B)=(1+0.5B)(1-0.9B),\quad \theta(B)=(1+0.5B)^2. \]

Cancelando el factor común \((1+0.5B)\):

\[ \begin{equation*} (1-0.9B)\,x_t=(1+0.5B)\,w_t \;\Rightarrow\; x_t=0.9\,x_{t-1}+0.5\,w_{t-1}+w_t, \end{equation*} \]

que es ARMA(1,1), causal (raíz de \(\phi(z)=1-0.9z\) es \(z=10/9>1\)) e invertible (raíz de \(\theta(z)=1+0.5z\) es \(z=-2\)).

Pesos \(\psi_j\) (linealización vía \(\phi(z)\psi(z)=\theta(z)\)):

\[ (1-0.9 z)\left(1+\psi_1 z+\psi_2 z^2+\cdots\right)=1+0.5 z \;\Rightarrow\; \psi_1=1.4,\ \psi_j=0.9\,\psi_{j-1}\ (j>1), \]

de modo que \(\psi_j=1.4(0.9)^{j-1},, j\ge 1\) y

\[ x_t=w_t+1.4\sum_{j=1}^{\infty} 0.9^{\,j-1} w_{t-j}. \]

En R:

ARMAtoMA(ar=.9, ma=.5, 10) # primeros 10 psi-weights

 [1] 1.4000000 1.2600000 1.1340000 1.0206000 0.9185400 0.8266860 0.7440174
 [8] 0.6696157 0.6026541 0.5423887

Pesos \(\pi_j\) (invertibilidad vía \(\theta(z)\pi(z)=\phi(z)\)):

\[ (1+0.5z)\left(1+\pi_1 z+\pi_2 z^2+\cdots\right)=1-0.9z \;\Rightarrow\; \pi_j=(-1)^j\,1.4\,(0.5)^{\,j-1},\ j\ge 1, \]

y como \(w_t=\sum_{j=0}^{\infty}\pi_j x_{t-j}\),

\[ x_t=1.4\sum_{j=1}^{\infty}(-0.5)^{\,j-1} x_{t-j}+w_t. \]

ARMAtoMA(ar=-.5, ma=-.9, 10) # primeros 10 pi-weights

 [1] -1.400000000  0.700000000 -0.350000000  0.175000000 -0.087500000
 [6]  0.043750000 -0.021875000  0.010937500 -0.005468750  0.002734375

3.2.6.2 Ejemplo 3.9 (Condiciones de causalidad para AR(2))

Para \((1-\phi_1 B-\phi_2 B^2)\,x_t=w_t\) es causal si las raíces de \(\phi(z)=1-\phi_1 z-\phi_2 z^2\) están fuera del círculo unidad. Equivalente en parámetros: \[ \begin{equation*} \phi_1+\phi_2<1,\quad \phi_2-\phi_1<1,\quad |\phi_2|<1. \end{equation*} \] (Define una región triangular de causalidad en el plano \((\phi_1,\phi_2)\).)

# (Opcional) Visualización de la región causal AR(2) con ggplot2:
# Grid de parámetros y prueba de desigualdades
grid <- expand_grid(phi1 = seq(-2, 2, by = 0.01),
                    phi2 = seq(-1, 1, by = 0.01)) |>
  mutate(causal = (phi1 + phi2 < 1) & (phi2 - phi1 < 1) & (abs(phi2) < 1))

ggplot(grid, aes(phi1, phi2, fill = causal)) +
  geom_raster() +
  scale_fill_manual(values = c("TRUE"="#4daf4a","FALSE"="#f0f0f0")) +
  labs(x=expression(phi[1]), y=expression(phi[2]),
       title="Región de causalidad para AR(2)") +
  theme_minimal()

3.3 Ecuaciones en diferencias

3.3.1 Ecuaciones en diferencias de primer orden

En esta sección se estudian las ecuaciones en diferencias y su relación con los procesos ARMA y sus funciones de autocorrelación. Se presentan las soluciones para ecuaciones homogéneas de primer y segundo orden, su generalización al caso de orden \(p\), y se muestran aplicaciones al cálculo de funciones de autocorrelación y de los pesos \(\psi\) en modelos ARMA.

Consideremos la secuencia \(u_{0}, u_{1}, u_{2}, \ldots\) que satisface:

\[ u_{n} - \alpha u_{n-1} = 0, \quad \alpha \neq 0, \quad n=1,2,\ldots \]

Ejemplo: la ACF de un AR(1) cumple:

\[ \rho(h) - \phi \rho(h-1) = 0, \quad h=1,2,\ldots \]

La solución es:

\[ u_{n} = \alpha^n u_0 \]

con condición inicial \(u_{0}=c\), se tiene \(u_{n} = \alpha^n c\).

En notación operador:

\[ (1 - \alpha B)u_n = 0 \]

y el polinomio asociado es \(\alpha(z) = 1 - \alpha z\), con raíz \(z_0 = 1/\alpha\). La solución:

\[ u_{n} = \alpha^n c = (z_0^{-1})^n c \]

3.3.2 Ecuaciones en diferencias de segundo orden

Ahora supongamos:

\[ u_n - \alpha_1 u_{n-1} - \alpha_2 u_{n-2} = 0, \quad \alpha_2 \neq 0, \quad n=2,3,\ldots \]

El polinomio asociado es:

\[ \alpha(z) = 1 - \alpha_1 z - \alpha_2 z^2 \]

con raíces \(z_1, z_2\).

• Caso \(z_1 \neq z_2\) (raíces distintas):

\[ u_n = c_1 z_1^{-n} + c_2 z_2^{-n} \]

• Caso \(z_1 = z_2 = z_0\) (raíz repetida):

\[ u_n = z_0^{-n}(c_1 + c_2 n) \]

3.3.3 Generalización a orden \(p\)

La ecuación en diferencias:

\[ u_n - \alpha_1 u_{n-1} - \cdots - \alpha_p u_{n-p} = 0, \quad \alpha_p \neq 0 \]

con polinomio asociado \(\alpha(z) = 1 - \alpha_1 z - \cdots - \alpha_p z^p\).

La solución general es:

\[ u_n = z_1^{-n} P_1(n) + z_2^{-n} P_2(n) + \cdots + z_r^{-n} P_r(n) \]

donde \(z_j\) son raíces con multiplicidad \(m_j\) y \(P_j(n)\) son polinomios de grado \(m_j - 1\).

3.3.3.1 Ejemplo: ACF de un AR(2)

Sea:

\[ x_t = \phi_1 x_{t-1} + \phi_2 x_{t-2} + w_t \]

Se obtiene:

\[ \gamma(h) = \phi_1 \gamma(h-1) + \phi_2 \gamma(h-2), \quad h=1,2,\ldots \]

Dividiendo por \(\gamma(0)\):

\[ \rho(h) - \phi_1 \rho(h-1) - \phi_2 \rho(h-2) = 0, \quad h=1,2,\ldots \]

La solución depende de las raíces de \(\phi(z) = 1 - \phi_1 z - \phi_2 z^2\):

• Raíces reales distintas: \(\rho(h) = c_1 z_1^{-h} + c_2 z_2^{-h}\)

• Raíz doble: \(\rho(h) = z_0^{-h}(c_1 + c_2 h)\)

• Raíces complejas conjugadas: \(\rho(h) = a |z_1|^{-h} \cos(h\theta + b)\)

Supongamos que las raíces son \(z_{1}\) y \(z_{2} = \bar{z}_{1}\), es decir, un par conjugado complejo. En este caso, los coeficientes de la solución cumplen \(c_{2} = \bar{c}_{1}\), ya que la función de autocorrelación \(\rho(h)\) debe ser un número real.

La solución general para \(\rho(h)\) es:

\[ \rho(h) = c_{1} z_{1}^{-h} + \bar{c}_{1} \bar{z}_{1}^{-h}. \]

Recordemos que cualquier número complejo \(z = x + iy\) se puede escribir en forma polar como:

\[ z = |z|(\cos \theta + i \sin \theta), \]

donde:

• \(|z| = \sqrt{x^{2} + y^{2}}\) es el módulo o distancia desde el origen hasta el punto \((x,y)\) en el plano complejo.
• \(\theta = \arg(z) = \arctan\left(\tfrac{y}{x}\right)\) es el argumento o ángulo medido desde el eje real positivo hasta el vector que representa a \(z\).

De forma compacta, esta representación se escribe como:

\[ z = |z| e^{i\theta}. \]

Aplicando esto a \(z_{1}\):

\[ z_{1} = |z_{1}| e^{i\theta}. \]

De aquí se sigue que:

\[ z_{1}^{-h} = |z_{1}|^{-h} e^{-i h \theta}. \]

Como \(c_{1}\) también puede escribirse en forma polar, la expresión completa de \(\rho(h)\) se puede reorganizar usando la identidad de Euler:

\[ e^{i\alpha} + e^{-i\alpha} = 2\cos(\alpha). \]

De este modo, la solución toma la forma:

\[ \rho(h) = a |z_{1}|^{-h} \cos(h\theta + b), \]

donde:

• \(a\) es una constante que depende de la magnitud de \(c_{1}\),
• \(b\) es un ángulo de fase determinado por las condiciones iniciales.

3.3.3.2 Ejemplo: AR(2) con raíces complejas

Modelo:

\[ x_t = 1.5 x_{t-1} - 0.75 x_{t-2} + w_t \]

Polinomio: \(\phi(z) = 1 - 1.5z + 0.75z^2\).

Las raíces: \(1 \pm i/\sqrt{3}\), con \(\theta = \tan^{-1}(1/\sqrt{3}) = 2\pi/12\).

En R:

z <- c(1,-1.5,.75)
a <- polyroot(z)[1]
a

[1] 1+0.5773503i

arg <- Arg(a)/(2*pi)
1/arg

[1] 12

Simulación:

set.seed(8675309)
ar2 <- arima.sim(list(order=c(2,0,0), ar=c(1.5,-.75)), n=144)
ar2 %>% 
  tsibble::as_tsibble() %>% 
  ggplot(aes(index, value)) +
  geom_line() +
  labs(x="Tiempo", y="Serie simulada")

Registered S3 method overwritten by 'tsibble':
  method               from 
  as_tibble.grouped_df dplyr

ACF en R:

acf_vals <- ARMAacf(ar=c(1.5,-.75), ma=0, lag.max=50)
tibble(lag=0:50, acf=acf_vals) %>%
  ggplot(aes(lag, acf)) +
  geom_col() +
  geom_hline(yintercept=0, color="red")

3.3.3.3 Ejemplo: Pesos \(\psi\) en un ARMA

Para un ARMA(\(p,q\)):

\[ \phi(B)x_t = \theta(B)w_t \]

con expansión:

\[ x_t = \sum_{j=0}^\infty \psi_j w_{t-j} \]

Las condiciones son:

\[ \psi_j - \sum_{k=1}^p \phi_k \psi_{j-k} = 0, \quad j \geq \max(p,q+1) \]

y

\[ \psi_j - \sum_{k=1}^j \phi_k \psi_{j-k} = \theta_j, \quad 0 \leq j < \max(p,q+1) \]

Ejemplo:

\[ x_t = 0.9 x_{t-1} + 0.5 w_{t-1} + w_t \]

Se obtiene \(\psi_0 = 1\), \(\psi_1 = 1.4\), y en general:

\[ \psi_j = 1.4 (0.9)^{j-1}, \quad j \geq 1 \]

En R:

ARMAtoMA(ar=0.9, ma=0.5, 50) %>% 
  tibble(lag=1:50, psi=.) %>%
  ggplot(aes(lag, psi)) +
  geom_line() +
  labs(x="Rezago", y="Peso psi")

3.4 Autocorrelación de procesos MA y ARMA

En esta sección se estudian las funciones de autocorrelación (ACF) y autocorrelación parcial (PACF) para procesos de series de tiempo. Se muestran sus propiedades en modelos MA, AR y ARMA, así como su utilidad para la identificación del orden de dependencia en estos procesos. Además, se presentan ejemplos numéricos y código en R para ilustrar los cálculos.

3.4.1 ACF de un proceso MA(q)

Un proceso \(x_t = \theta(B) w_t\), con \(\theta(B) = 1 + \theta_1 B + \cdots + \theta_q B^q\), tiene media cero y función de autocovarianza:

\[ \gamma(h) = \begin{cases} \sigma_w^2 \sum_{j=0}^{q-h} \theta_j \theta_{j+h}, & 0 \leq h \leq q \\ 0, & h > q \end{cases} \]

Al dividir por \(\gamma(0)\), se obtiene la ACF:

\[ \rho(h) = \begin{cases} \dfrac{\sum_{j=0}^{q-h} \theta_j \theta_{j+h}}{1+\theta_1^2+\cdots+\theta_q^2}, & 1 \leq h \leq q \\ 0, & h > q \end{cases} \]

3.4.2 ACF de un proceso ARMA(p, q)

Para un ARMA causal \(\phi(B) x_t = \theta(B) w_t\), se escribe:

\[ x_t = \sum_{j=0}^\infty \psi_j w_{t-j} \]

con \(\mathrm{E}(x_t)=0\) y

\[ \gamma(h) = \sigma_w^2 \sum_{j=0}^\infty \psi_j \psi_{j+h}, \quad h \geq 0 \]

La ecuación de recurrencia para \(\gamma(h)\) es:

\[ \gamma(h) - \phi_1 \gamma(h-1) - \cdots - \phi_p \gamma(h-p) = 0, \quad h \geq \max(p, q+1) \]

con condiciones iniciales:

\[ \gamma(h) - \sum_{j=1}^p \phi_j \gamma(h-j) = \sigma_w^2 \sum_{j=h}^q \theta_j \psi_{j-h}, \quad 0 \leq h < \max(p,q+1) \]

3.4.3 Ejemplos de ACF

3.4.3.1 Ejemplo: ACF de un AR(p)

Se cumple:

\[ \rho(h) - \phi_1 \rho(h-1) - \cdots - \phi_p \rho(h-p) = 0, \quad h \geq p \]

La solución general es:

\[ \rho(h) = z_1^{-h} P_1(h) + \cdots + z_r^{-h} P_r(h), \quad h \geq p \]

donde \(z_i\) son raíces de \(\phi(z)\) y \(P_j(h)\) son polinomios.

3.4.3.2 Ejemplo: ACF de un ARMA (1, 1)

Considere el proceso ARMA(1,1) \(x_{t}=\phi\, x_{t-1}+\theta\, w_{t-1}+w_{t}\), donde \(|\phi|<1\). La función de autocovarianza satisface

\[ \gamma(h)-\phi\, \gamma(h-1)=0, \quad h=2,3,\ldots \]

y por lo tanto la solución general es

\[ \begin{equation*} \gamma(h)=c\, \phi^{h}, \quad h=1,2,\ldots \end{equation*} \]

Para obtener las condiciones iniciales usamos:

\[ \gamma(0)=\phi\, \gamma(1)+\sigma_{w}^{2}\big[1+\theta \phi+\theta^{2}\big] \quad \text{y} \quad \gamma(1)=\phi\, \gamma(0)+\sigma_{w}^{2}\theta . \]

Al resolver para \(\gamma(0)\) y \(\gamma(1)\), obtenemos:

\[ \gamma(0)=\sigma_{w}^{2}\,\frac{1+2\theta \phi+\theta^{2}}{1-\phi^{2}} \quad \text{y} \quad \gamma(1)=\sigma_{w}^{2}\,\frac{(1+\theta \phi)(\phi+\theta)}{1-\phi^{2}} . \]

Para resolver \(c\), note que, \(\gamma(1)=c\,\phi\), es decir, \(c=\gamma(1)/\phi\). Por lo tanto, la solución específica para \(h \ge 1\) es

\[ \gamma(h)=\frac{\gamma(1)}{\phi}\,\phi^{h} =\sigma_{w}^{2}\,\frac{(1+\theta \phi)(\phi+\theta)}{1-\phi^{2}}\,\phi^{\,h-1}. \]

Finalmente, al dividir por \(\gamma(0)\) se obtiene la ACF

\[ \begin{equation*} \rho(h)=\frac{(1+\theta \phi)(\phi+\theta)}{1+2\theta \phi+\theta^{2}}\,\phi^{\,h-1}, \quad h \ge 1 \end{equation*} \]

Obsérvese que el patrón general de \(\rho(h)\) frente a \(h\) no es distinto del de un AR(1). Por lo tanto, es poco probable que podamos distinguir entre un \(\operatorname{ARMA}(1,1)\) y un \(\operatorname{AR}(1)\) basándonos únicamente en una ACF estimada a partir de una muestra. Esta consideración nos llevará a la función de autocorrelación parcial.

3.4.4 Función de autocorrelación parcial (PACF)

Anteriormente vimos que, para un modelo MA(\(q\)), la función de autocorrelación (ACF) se anula para rezagos mayores que \(q\). Además, como \(\theta_q \neq 0\), la ACF no es cero en el rezago \(q\). Esto significa que la ACF aporta bastante información sobre el orden de dependencia en procesos de media móvil.

Sin embargo, si el proceso es ARMA o AR, la ACF por sí sola no permite identificar claramente el orden. Por ello, se introduce la función de autocorrelación parcial (PACF), que se comporta para los modelos AR de forma similar a la ACF en los modelos MA.

3.4.4.1 Definición intuitiva de correlación parcial

Si \(X\), \(Y\) y \(Z\) son variables aleatorias, la correlación parcial entre \(X\) y \(Y\) dado \(Z\) se obtiene mediante:

1. Regresar \(X\) sobre \(Z\) para obtener \(\hat{X}\).
2. Regresar \(Y\) sobre \(Z\) para obtener \(\hat{Y}\).
3. Calcular la correlación entre los residuos:

\[ \rho_{XY|Z} = \operatorname{corr}\{X - \hat{X}, \, Y - \hat{Y}\}. \]

La idea es medir la asociación de \(X\) y \(Y\) una vez eliminada la influencia lineal de \(Z\). En el caso normal multivariado, esta definición coincide con:

\[ \rho_{XY|Z} = \operatorname{corr}(X,Y \mid Z). \]

3.4.4.2 Ejemplo

Consideremos un modelo causal AR(1):

\[ x_t = \phi x_{t-1} + w_t. \]

El cálculo de la covarianza en el rezago 2 es:

\[ \begin{aligned} \gamma_x(2) &= \operatorname{cov}(x_t, x_{t-2}) \\ &= \operatorname{cov}(\phi x_{t-1} + w_t, \, x_{t-2}) \\ &= \operatorname{cov}(\phi^2 x_{t-2} + \phi w_{t-1} + w_t, \, x_{t-2}) \\ &= \phi^2 \gamma_x(0). \end{aligned} \]

Este resultado se debe a la causalidad: \(x_{t-2}\) depende de \({w_{t-2}, w_{t-3},\ldots}\), que son no correlacionados con \(w_t\) y \(w_{t-1}\). Así, \(x_t\) y \(x_{t-2}\) están correlacionados a través de \(x_{t-1}\).

Para romper esta cadena de dependencia, eliminamos el efecto de \(x_{t-1}\):

\[ \operatorname{cov}\big(x_t - \phi x_{t-1}, \, x_{t-2} - \phi x_{t-1}\big) = \operatorname{cov}(w_t, x_{t-2} - \phi x_{t-1}) = 0. \]

De aquí surge la PACF: la correlación entre \(x_s\) y \(x_t\) una vez eliminada la influencia de las variables intermedias.

Para una serie estacionaria de media cero:

• Sea \(\hat{x}_{t+h}\) la regresión de \(x_{t+h}\) en \({x_{t+h-1}, x_{t+h-2}, \ldots, x_{t+1}}\):

\[ \hat{x}_{t+h} = \beta_1 x_{t+h-1} + \beta_2 x_{t+h-2} + \cdots + \beta_{h-1} x_{t+1}. \]

• Sea \(\hat{x}_t\) la regresión de \(x_t\) en \({x_{t+1}, \ldots, x_{t+h-1}}\):

\[ \hat{x}_t = \beta_1 x_{t+1} + \beta_2 x_{t+2} + \cdots + \beta_{h-1} x_{t+h-1}. \]

Por estacionariedad, los coeficientes \(\beta_1, \ldots, \beta_{h-1}\) son iguales en ambas expresiones.

Formalmente, la PACF mide la correlación entre \(x_{t+h}\) y \(x_t\) eliminando el efecto lineal de los valores intermedios. Formalmente:

• Para \(h=1\): \(\phi_{11} = \rho(1)\)
• Para \(h \geq 2\):

\[ \phi_{hh} = \operatorname{corr}(x_{t+h} - \hat{x}_{t+h}, x_t - \hat{x}_t) \]

donde \(\hat{x}_{t+h}\) y \(\hat{x}_t\) son regresiones sobre los valores intermedios.

3.4.4.3 Ejemplo: PACF de un AR(1)

Consideremos el proceso AR(1):

\[ x_t = \phi x_{t-1} + w_t, \quad |\phi| < 1. \]

Por definición, el primer valor de la PACF es simplemente:

\[ \phi_{11} = \rho(1) = \phi. \]

Cálculo de \(\phi_{22}\)

1. Regresión de \(x_{t+2}\) sobre \(x_{t+1}\) Supongamos \(\hat{x}_{t+2} = \beta x_{t+1}\). Para encontrar \(\beta\), minimizamos el error cuadrático esperado:

   \[ \mathrm{E}\left[(x_{t+2} - \hat{x}_{t+2})^2\right] = \mathrm{E}\left[(x_{t+2} - \beta x_{t+1})^2\right] = \gamma(0) - 2\beta \gamma(1) + \beta^2 \gamma(0). \]

   Derivando respecto a \(\beta\) y resolviendo, se obtiene:

   \[ \beta = \frac{\gamma(1)}{\gamma(0)} = \rho(1) = \phi. \]

2. Regresión de \(x_t\) sobre \(x_{t+1}\) Supongamos \(\hat{x}_t = \beta x_{t+1}\). Nuevamente, al minimizar:

   \[ \mathrm{E}\left[(x_t - \hat{x}_t)^2\right] = \mathrm{E}\left[(x_t - \beta x_{t+1})^2\right] = \gamma(0) - 2\beta \gamma(1) + \beta^2 \gamma(0), \] Con estos valores:

\[ \begin{aligned} \phi_{22} &= \operatorname{corr}(x_{t+2} - \hat{x}_{t+2}, \, x_t - \hat{x}_t) \\ &= \operatorname{corr}(x_{t+2} - \phi x_{t+1}, \, x_t - \phi x_{t+1}) \\ &= \operatorname{corr}(w_{t+2}, \, x_t - \phi x_{t+1}) \\ &= 0, \end{aligned} \]

debido a la causalidad: \(w_{t+2}\) es no correlacionado con las variables del pasado.

3.4.4.4 Ejemplo: PACF de un AR(p)

Para \(h > p\), \(\phi_{hh} = 0\) por causalidad. La PACF corta en el rezago \(p\).

Código en R:

ACF <- ARMAacf(ar=c(1.5, -0.75), ma=0, 24)[-1]
PACF <- ARMAacf(ar=c(1.5, -0.75), ma=0, 24, pacf=TRUE)

par(mfrow=c(1,2))
plot(ACF, type="h", xlab="Rezago", ylim=c(-0.8,1)); abline(h=0)
plot(PACF, type="h", xlab="Rezago", ylim=c(-0.8,1)); abline(h=0)

3.4.4.5 Ejemplo 3.17: PACF de un MA(q)

Para un MA(1), \(x_t = w_t + \theta w_{t-1}\):

\[ \phi_{hh} = -\frac{(-\theta)^h (1-\theta^2)}{1-\theta^{2(h+1)}}, \quad h \geq 1 \]

En general, la PACF de un MA no corta, a diferencia de un AR.

3.4.5 Ejemplo: Análisis preliminar de la serie Recruitment

Con 453 meses de datos (1950–1987), la ACF muestra ciclos de 12 meses y la PACF valores altos en \(h=1,2\), sugiriendo un AR(2):

\[ x_t = \phi_0 + \phi_1 x_{t-1} + \phi_2 x_{t-2} + w_t \]

Estimaciones:

• \(\hat{\phi}_0 = 6.74 , (1.11)\)
• \(\hat{\phi}_1 = 1.35 , (0.04)\)
• \(\hat{\phi}_2 = -0.46 , (0.04)\)
• \(\hat{\sigma}_w^2 = 89.72\)

Código en R:

acf2(rec, 48) # ACF y PACF

     [,1]  [,2]  [,3]  [,4] [,5]  [,6]  [,7] [,8] [,9] [,10] [,11] [,12] [,13]
ACF  0.92  0.78  0.63  0.48 0.36  0.26  0.18 0.13 0.09  0.07  0.06  0.02 -0.04
PACF 0.92 -0.44 -0.05 -0.02 0.07 -0.03 -0.03 0.04 0.05 -0.02 -0.05 -0.14 -0.15
     [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24] [,25]
ACF  -0.12 -0.19 -0.24 -0.27 -0.27 -0.24 -0.19 -0.11 -0.03  0.03  0.06  0.06
PACF -0.05  0.05  0.01  0.01  0.02  0.09  0.11  0.03 -0.03 -0.01 -0.07 -0.12
     [,26] [,27] [,28] [,29] [,30] [,31] [,32] [,33] [,34] [,35] [,36] [,37]
ACF   0.02 -0.02 -0.06 -0.09 -0.12 -0.13 -0.11 -0.05  0.02  0.08  0.12  0.10
PACF -0.03  0.05 -0.08 -0.04 -0.03  0.06  0.05  0.15  0.09 -0.04 -0.10 -0.09
     [,38] [,39] [,40] [,41] [,42] [,43] [,44] [,45] [,46] [,47] [,48]
ACF   0.06  0.01 -0.02 -0.03 -0.03 -0.02  0.01  0.06  0.12  0.17  0.20
PACF -0.02  0.05  0.08 -0.02 -0.01 -0.02  0.05  0.01  0.05  0.08 -0.04

regr <- ar.ols(rec, order=2, demean=FALSE, intercept=TRUE)
regr$asy.se.coef # errores estándar

$x.mean
[1] 1.110599

$ar
[1] 0.04178901 0.04187942

3.4.6 Resumen final

Modelo	ACF	PACF
AR(p)	Decae	Corta en p
MA(q)	Corta en q	Decae
ARMA(p,q)	Decae	Decae

3.5 Pronóstico

Este apartado presenta los fundamentos del pronóstico en series de tiempo estacionarias con parámetros conocidos. Se define el mejor predictor en media cuadrática como una esperanza condicional y, dentro de los predictores lineales, se desarrollan los mejores predictores lineales (BLP). Se introducen ecuaciones de predicción en forma escalar y matricial, el algoritmo de Durbin–Levinson, el algoritmo de innovaciones, así como la formulación de pronósticos ARMA (incluyendo versiones truncadas), pronósticos de largo plazo y backcasting.

3.5.1 Objetivo de pronóstico y predictor MMCE

Dado \(x_{1:n}={x_1,\dots,x_n}\) y una serie estacionaria \(x_t\) con parámetros conocidos, el objetivo es predecir \(x_{n+m}\), \(m=1,2,\ldots\). El predictor que minimiza el error cuadrático medio (MMCE) es

\[ x_{n+m}^{\,n}=\mathrm{E}\!\left(x_{n+m}\mid x_{1:n}\right), \]

pues minimiza

\[ \mathrm{E}\!\left[x_{n+m}-g(x_{1:n})\right]^2. \]

3.5.2 Mejores predictores lineales (BLP)

Se restringe la atención a predictores lineales:

\[ x_{n+m}^{\,n}=\alpha_0+\sum_{k=1}^n \alpha_k x_k, \]

donde los coeficientes dependen de \(n\) y \(m\). Los BLP son aquellos que minimizan (en el conjunto lineal) el error cuadrático. Bajo gaussianidad, MMCE y BLP coinciden.

Propiedad (Mejor predicción lineal para procesos estacionarios). Dado \(x_{1},\ldots,x_n\), el BLP \(x_{n+m}^{,n}=\alpha_0+\sum_{k=1}^n \alpha_k x_k\) satisface

\[ \mathrm{E}\!\left[\left(x_{n+m}-x_{n+m}^{\,n}\right)x_k\right]=0,\quad k=0,1,\ldots,n, \]

donde \(x_0=1\), para resolver \({\alpha_0,\alpha_1,\ldots,\alpha_n}\). Si \(\mathrm{E}(x_t)=\mu\), el primer momento implica \(\mathrm{E}(x_{n+m}^{,n})=\mu\) y

\[ \alpha_0=\mu\Big(1-\sum_{k=1}^n \alpha_k\Big), \]

por lo que, sin pérdida de generalidad (mientras no estimemos), se toma \(\mu=0\) y \(\alpha_0=0\).

3.5.3 Predicción un paso adelante y forma matricial

El BLP un paso adelante tiene forma

\[ x_{n+1}^{\,n}=\phi_{n1}x_n+\phi_{n2}x_{n-1}+\cdots+\phi_{nn}x_1. \]

Los coeficientes satisfacen

\[ \sum_{j=1}^n \phi_{nj}\,\gamma(k-j)=\gamma(k),\quad k=1,\ldots,n, \]

que en notación matricial son

\[ \Gamma_n\,\phi_n=\gamma_n, \]

con \(\Gamma_n={\gamma(k-j)}_{j,k=1}^n\), \(\phi_n=(\phi_{n1},\ldots,\phi_{nn})'\) y \(\gamma_n=(\gamma(1),\ldots,\gamma(n))'\). Si \(\Gamma_n\) es no singular,

\[ \phi_n=\Gamma_n^{-1}\gamma_n. \]

El error cuadrático medio (ECM) un paso adelante es

\[ P_{n+1}^{\,n}=\mathrm{E}\!\left(x_{n+1}-x_{n+1}^{\,n}\right)^2=\gamma(0)-\gamma_n'\Gamma_n^{-1}\gamma_n. \]

3.5.3.1 Ejemplo: Predicción para un AR(2)

Para \(x_t=\phi_1x_{t-1}+\phi_2x_{t-2}+w_t\) causal:

• Con un dato \(x_1\):

\[ x_2^{\,1}=\phi_{11}x_1=\frac{\gamma(1)}{\gamma(0)}x_1=\rho(1)x_1. \]

• Con \(x_1,x_2\):

\[ x_3^{\,2}=\phi_{21}x_2+\phi_{22}x_1, \]

y por el modelo, \(x_3^{,2}=\phi_1x_2+\phi_2x_1\), de modo que \(\phi_{21}=\phi_1,\ \phi_{22}=\phi_2\). En general, para \(n\ge2\):

\[ x_{n+1}^{\,n}=\phi_1x_n+\phi_2x_{n-1}. \]

Generalización AR(\(p\)): para \(n\ge p\),

\[ x_{n+1}^{\,n}=\phi_1x_n+\phi_2x_{n-1}+\cdots+\phi_p x_{n-p+1}. \]

3.5.4 Algoritmo de Durbin–Levinson

El sistema \(\Gamma_n \phi_n = \gamma_n\) es costoso de resolver para grandes \(n\). Para evitar invertir matrices grandes, se usan soluciones recursivas:

• Condiciones iniciales:

\[ \phi_{00}=0,\quad P_{1}^{\,0}=\gamma(0). \]

• Para \(n\ge1\):

\[ \phi_{nn}=\frac{\rho(n)-\sum_{k=1}^{n-1}\phi_{n-1,k}\rho(n-k)}{1-\sum_{k=1}^{n-1}\phi_{n-1,k}\rho(k)},\qquad P_{n+1}^{\,n}=P_n^{\,n-1}\left(1-\phi_{nn}^2\right), \]

y, para \(k=1,\ldots,n-1\),

\[ \phi_{nk}=\phi_{n-1,k}-\phi_{nn}\,\phi_{n-1,n-k}. \]

Error de predicción un paso adelante (forma producto):

\[ P_{n+1}^{\,n}=\gamma(0)\prod_{j=1}^n\left[1-\phi_{jj}^2\right]. \]

3.5.5 PACF vía Durbin–Levinson

Como consecuencia de Durbin–Levinson:

Propiedad La PACF se obtiene iterativamente como \(\phi_{nn}\), \(n=1,2,\ldots\).

Para un AR(\(p\)), tomando \(n=p\):

\[ x_{p+1}^{\,p}=\phi_{p1}x_p+\cdots+\phi_{pp}x_1=\phi_1x_p+\cdots+\phi_p x_1, \]

lo que muestra que \(\phi_{pp}=\phi_p\).

3.5.5.1 Ejemplo: PACF de un AR(2)

Usando \(\rho(h)-\phi_1\rho(h-1)-\phi_2\rho(h-2)=0\) (para \(h\ge1\)):

\[ \phi_{11}=\rho(1)=\frac{\phi_1}{1-\phi_2},\qquad \phi_{22}=\frac{\rho(2)-\rho(1)^2}{1-\rho(1)^2}=\phi_2,\qquad \phi_{21}=\rho(1)(1-\phi_2)=\phi_1,\qquad \phi_{33}=0. \]

3.5.6 Predicción a \(m\) pasos

Dado \(x_{1:n}\), el BLP a \(m\) pasos es

\[ x_{n+m}^{\,n}=\phi_{n1}^{(m)}x_n+\phi_{n2}^{(m)}x_{n-1}+\cdots+\phi_{nn}^{(m)}x_1, \]

con

\[ \sum_{j=1}^n \phi_{nj}^{(m)}\gamma(k-j)=\gamma(m+k-1),\quad k=1,\ldots,n, \]

o en forma matricial

\[ \Gamma_n\,\phi_n^{(m)}=\gamma_n^{(m)},\quad \gamma_n^{(m)}=(\gamma(m),\ldots,\gamma(m+n-1))'. \]

El ECM es

\[ P_{n+m}^{\,n}=\gamma(0)-\gamma_n^{(m)'}\Gamma_n^{-1}\gamma_n^{(m)}. \]

3.5.7 Algoritmo de innovaciones

Para un proceso estacionario de media cero:

• Inicialización:

\[ x_1^{\,0}=0,\quad P_1^{\,0}=\gamma(0). \]

• Para \(t=1,2,\ldots\):

\[ x_{t+1}^{\,t}=\sum_{j=1}^t \theta_{tj}\big(x_{t+1-j}-x_{t+1-j}^{\,t-j}\big), \]

\[ P_{t+1}^{\,t}=\gamma(0)-\sum_{j=0}^{t-1}\theta_{t,t-j}^2\,P_{j+1}^{\,j}, \]

con, para \(j=0,1,\ldots,t-1\),

\[ \theta_{t,t-j}=\frac{\gamma(t-j)-\sum_{k=0}^{j-1}\theta_{j,j-k}\theta_{t,t-k}P_{k+1}^{\,k}}{P_{j+1}^{\,j}}. \]

El predictor a \(m\) pasos y su ECM:

\[ x_{n+m}^{\,n}=\sum_{j=m}^{n+m-1}\theta_{n+m-1,j}\big(x_{n+m-j}-x_{n+m-j}^{\,n+m-j-1}\big), \]

\[ P_{n+m}^{\,n}=\gamma(0)-\sum_{j=m}^{n+m-1}\theta_{n+m-1,j}^2\,P_{n+m-j}^{\,n+m-j-1}. \]

Nota: las innovaciones \(x_t-x_t^{t-1}\) y \(x_s-x_s^{s-1}\) son no-correlacionadas para \(t\neq s\).

3.5.7.1 Ejemplo: Predicción para un MA(1)

Para \(x_t=w_t+\theta w_{t-1}\), con \(\gamma(0)=(1+\theta^2)\sigma_w^2\), \(\gamma(1)=\theta\sigma_w^2\), \(\gamma(h)=0\) si \(h>1\):

\[ \theta_{n1}=\theta\sigma_w^2/P_n^{\,n-1},\qquad \theta_{nj}=0\ (j\ge2), \]

\[ P_1^{\,0}=(1+\theta^2)\sigma_w^2,\qquad P_{n+1}^{\,n}=(1+\theta^2-\theta\theta_{n1})\sigma_w^2, \]

y

\[ x_{n+1}^{\,n}=\theta\,(x_n-x_n^{\,n-1})\,\sigma_w^2/ P_n^{\,n-1}. \]

3.5.8 Pronósticos ARMA: forma causal/invertible e historia infinita

Asumiendo \(x_t\) causal e invertible ARMA(\(p,q\)), \(\phi(B)x_t=\theta(B)w_t\), con \(w_t\sim\text{iid }N(0,\sigma_w^2)\). Definimos

\[ \tilde{x}_{n+m}=\mathrm{E}\!\left(x_{n+m}\mid x_n,x_{n-1},\ldots\right). \]

Formas causal e invertible:

\[ x_{n+m}=\sum_{j=0}^\infty \psi_j w_{n+m-j},\ \ \psi_0=1;\qquad w_{n+m}=\sum_{j=0}^\infty \pi_j x_{n+m-j},\ \ \pi_0=1. \]

Entonces

\[ \tilde{x}_{n+m}=\sum_{j=0}^\infty \psi_j \tilde{w}_{n+m-j}=\sum_{j=m}^\infty \psi_j w_{n+m-j}, \]

y

\[ \tilde{x}_{n+m}=-\sum_{j=1}^{m-1}\pi_j \tilde{x}_{n+m-j}-\sum_{j=m}^{\infty}\pi_j x_{n+m-j}. \]

El ECM:

\[ P_{n+m}^{\,n}=\sigma_w^2\sum_{j=0}^{m-1}\psi_j^2, \]

y la covarianza de errores de distinto horizonte (\(k\ge1\)):

\[ \mathrm{E}\!\left[(x_{n+m}-\tilde{x}_{n+m})(x_{n+m+k}-\tilde{x}_{n+m+k})\right]=\sigma_w^2\sum_{j=0}^{m-1}\psi_j\psi_{j+k}. \]

3.5.8.1 Ejemplo: Pronósticos de largo plazo

Si \(\mu_x\neq 0\), entonces

\[ \tilde{x}_{n+m}=\mu_x+\sum_{j=m}^\infty \psi_j w_{n+m-j}\ \ \Rightarrow\ \ \tilde{x}_{n+m}\to \mu_x, \]

exponencialmente en media cuadrática cuando \(m\to\infty\). Además,

\[ P_{n+m}^{\,n}\to \sigma_w^2\sum_{j=0}^\infty \psi_j^2=\gamma_x(0)=\sigma_x^2. \]

3.5.9 Pronóstico truncado (historia finita)

Como sólo observamos \(x_1,\ldots,x_n\), truncamos:

\[ \tilde{x}_{n+m}^{\,n}=-\sum_{j=1}^{m-1}\pi_j \tilde{x}_{n+m-j}^{\,n}-\sum_{j=m}^{n+m-1}\pi_j x_{n+m-j}. \]

Para AR(\(p\)) y \(n>p\), el predictor \(x_{n+m}^{n}\) es exacto y el error un paso es \(\sigma_w^2\).

3.5.10 Propiedad (Predicción truncada ARMA)

Para ARMA(\(p,q\)):

\[ \tilde{x}_{n+m}^{\,n}=\phi_1\tilde{x}_{n+m-1}^{\,n}+\cdots+\phi_p\tilde{x}_{n+m-p}^{\,n} +\theta_1\tilde{w}_{n+m-1}^{\,n}+\cdots+\theta_q\tilde{w}_{n+m-q}^{\,n}, \]

donde \(\tilde{x}_t^{,n}=x_t\) si \(1\le t\le n\) y \(\tilde{x}_t^{,n}=0\) si \(t\le0\). Los errores truncados cumplen \(\tilde{w}_t^{,n}=0\) para \(t\le0\) o \(t>n\) y, para \(1\le t\le n\),

\[ \tilde{w}_t^{\,n}=\phi(B)\tilde{x}_t^{\,n}-\theta_1\tilde{w}_{t-1}^{\,n}-\cdots-\theta_q\tilde{w}_{t-q}^{\,n}. \]

3.5.10.1 Ejemplo: Pronóstico ARMA(1,1)

Modelo:

\[ x_{n+1}=\phi x_n+w_{n+1}+\theta w_n. \]

Predicción truncada un paso:

\[ \tilde{x}_{n+1}^{\,n}=\phi x_n+\theta \tilde{w}_n^{\,n}. \]

Para \(m\ge2\):

\[ \tilde{x}_{n+m}^{\,n}=\phi\,\tilde{x}_{n+m-1}^{\,n}. \]

Inicialización de errores truncados (\(\tilde{w}_0^{,n}=0\), \(x_0=0\)):

\[ \tilde{w}_t^{\,n}=x_t-\phi x_{t-1}-\theta \tilde{w}_{t-1}^{\,n},\quad t=1,\ldots,n. \]

ECM aproximado usando \(\psi_j=(\phi+\theta)\phi^{,j-1}\) para \(j\ge1\):

\[ P_{n+m}^{\,n}=\sigma_w^2\left[1+(\phi+\theta)^2\sum_{j=1}^{m-1}\phi^{2(j-1)}\right] =\sigma_w^2\left[1+\frac{(\phi+\theta)^2\left(1-\phi^{2(m-1)}\right)}{1-\phi^2}\right]. \]

Intervalos de predicción \((1-\alpha)\):

\[ x_{n+m}^{\,n}\pm c_{\alpha/2}\sqrt{P_{n+m}^{\,n}}. \]

3.5.10.2 Ejemplo: Pronóstico de la serie Recruitment

Usando estimaciones como valores verdaderos, para horizonte \(m=1,\ldots,24\) y \(n=453\):

\[ x_{n+m}^{\,n}=6.74+1.35\,x_{n+m-1}^{\,n}-0.46\,x_{n+m-2}^{\,n}. \]

Con \(\hat{\sigma}_w^2=89.72\) y recurriendo a \(\psi_0=1\), \(\psi_1=1.35\), \(\psi_j=1.35\psi_{j-1}-0.46\psi_{j-2}\):

\[ P_{n+1}^{\,n}=89.72,\quad P_{n+2}^{\,n}=89.72\,(1+1.35^2),\quad P_{n+3}^{\,n}=89.72\big(1+1.35^2+[1.35^2-0.46]^2\big),\ \ldots \]

regr <- ar.ols(rec, order = 2, demean = FALSE, intercept = TRUE)
fore <- predict(regr, n.ahead = 24)
ts.plot(rec, fore$pred, col = 1:2, xlim = c(1980, 1990), ylab = "Recruitment")
U <- fore$pred + fore$se
L <- fore$pred - fore$se
xx <- c(time(U), rev(time(U)))
yy <- c(L, rev(U))
polygon(xx, yy, border = 8, col = gray(.6, alpha = .2))
lines(fore$pred, type = "p", col = 2)

3.6 Estimación

En esta sección se estudian los métodos de estimación de parámetros en procesos ARMA \((p, q)\) causales e invertibles. Se asume que los órdenes \(p\) y \(q\) son conocidos, y se busca estimar los parámetros \(\phi_{1}, \ldots, \phi_{p}\), \(\theta_{1}, \ldots, \theta_{q}\) y \(\sigma_{w}^{2}\). Se presentan los estimadores por momentos (método de Yule-Walker), los de máxima verosimilitud y los de mínimos cuadrados (condicionales e incondicionales), junto con sus propiedades asintóticas y ejemplos ilustrativos.

3.6.1 Estimación por el método de momentos

3.6.1.1 Modelos AR(p) y ecuaciones de Yule-Walker

Para el modelo

\[ x_{t}=\phi_{1}x_{t-1}+\cdots+\phi_{p}x_{t-p}+w_{t}, \]

las ecuaciones de Yule-Walker se definen como:

\[ \gamma(h)=\phi_{1}\gamma(h-1)+\cdots+\phi_{p}\gamma(h-p), \quad h=1,\ldots,p \]

\[ \sigma_{w}^{2}=\gamma(0)-\phi_{1}\gamma(1)-\cdots-\phi_{p}\gamma(p) \]

En forma matricial:

\[ \Gamma_{p}\phi=\gamma_{p}, \quad \sigma_{w}^{2}=\gamma(0)-\phi^{\prime}\gamma_{p}, \]

donde \(\Gamma_{p}={\gamma(k-j)}_{j,k=1}^{p}\), \(\phi=(\phi_{1},\ldots,\phi_{p})^{\prime}\) y \(\gamma_{p}=(\gamma(1),\ldots,\gamma(p))^{\prime}\).

Los estimadores de Yule-Walker se obtienen reemplazando momentos poblacionales por muestrales:

\[ \hat{\phi}=\hat{\Gamma}_{p}^{-1}\hat{\gamma}_{p}, \quad \hat{\sigma}_{w}^{2}=\hat{\gamma}(0)-\hat{\gamma}_{p}^{\prime}\hat{\Gamma}_{p}^{-1}\hat{\gamma}_{p} \]

En términos de la ACF muestral \(\hat{R}_{p}=\{\hat \rho(k-j)\}_{j,k=1}^p\) y asumiendo \(\hat \rho_p=(\hat \rho(1),\ldots,\rho(p))\):

\[ \hat{\phi}=\hat{R}_{p}^{-1}\hat{\rho}_{p}, \quad \hat{\sigma}_{w}^{2}=\hat{\gamma}(0)\left[1-\hat{\rho}_{p}^{\prime}\hat{R}_{p}^{-1}\hat{\rho}_{p}\right]. \]

Propiedad: Resultados asintóticos

Si \(n\to\infty\):

\[ \sqrt{n}(\hat{\phi}-\phi)\xrightarrow{d}N(0,\sigma_{w}^{2}\Gamma_{p}^{-1}), \quad \hat{\sigma}_{w}^{2}\xrightarrow{p}\sigma_{w}^{2}. \]

Propiedad: Distribución del PACF

Para \(h>p\):

\[ \sqrt{n}\,\hat{\phi}_{hh}\xrightarrow{d}N(0,1). \]

3.6.2 Ejemplos de estimación Yule-Walker

3.6.2.1 Ejemplo: AR(2) simulado

Modelo:

\[ x_{t}=1.5x_{t-1}-0.75x_{t-2}+w_{t}, \quad w_{t}\sim N(0,1). \]

Con \(n=144\), se obtuvo:

\[ \hat{\phi}=\binom{1.463}{-0.723}, \quad \hat{\sigma}_{w}^{2}=1.187. \]

Varianzas asintóticas:

\[ \begin{bmatrix} 0.058^{2} & -0.003 \\ -0.003 & 0.058^{2} \end{bmatrix}. \] Código en R:

# Ilustración en R del ejemplo AR(2) con n = 144
# Modelo: x_t = 1.5 x_{t-1} - 0.75 x_{t-2} + w_t,  w_t ~ N(0,1)

set.seed(8675309)  # misma semilla usada en ejemplos clásicos
n <- 144
x <- arima.sim(n = n, list(order = c(2,0,0), ar = c(1.5, -0.75)), sd = 1)

# 1) Estimación Yule-Walker (coeficientes y varianza del error)
fit_yw <- stats::ar.yw(x, order.max = 2, aic = FALSE, demean = TRUE)
phi_hat <- as.numeric(fit_yw$ar)            # (phi1_hat, phi2_hat)
sigma2_hat <- as.numeric(fit_yw$var.pred)   # estimador de sigma_w^2

cat("Coeficientes YW (phi1_hat, phi2_hat):\n")

Coeficientes YW (phi1_hat, phi2_hat):

print(round(phi_hat, 3))

[1]  1.438 -0.705

cat("\nVarianza del ruido (sigma_w^2_hat):\n")

Varianza del ruido (sigma_w^2_hat):

print(round(sigma2_hat, 3))

[1] 1.327

# 2) Matriz Γ_p estimada con autocovarianzas muestrales
#    Γ_2 = [[γ(0), γ(1)],
#            [γ(1), γ(0)]]
acf_cov <- acf(x, type = "covariance", plot = FALSE)$acf
gamma0_hat <- as.numeric(acf_cov[1])  # γ(0)
gamma1_hat <- as.numeric(acf_cov[2])  # γ(1)

Gamma_hat <- matrix(c(gamma0_hat, gamma1_hat,
                      gamma1_hat, gamma0_hat), nrow = 2, byrow = TRUE)

# 3) Varianzas asintóticas de (phi1_hat, phi2_hat):
#    Var_asint(hat{phi}) ≈ (sigma_w^2_hat / n) * Γ_p^{-1}
asymp_cov <- (sigma2_hat / n) * solve(Gamma_hat)

cat("\nMatriz de varianzas-covarianzas asintótica (~ (sigma_w^2/n) * Γ^{-1}):\n")

Matriz de varianzas-covarianzas asintótica (~ (sigma_w^2/n) * Γ^{-1}):

print(round(asymp_cov, 6))

          [,1]      [,2]
[1,]  0.003566 -0.003007
[2,] -0.003007  0.003566

cat("\nComparación con la forma reportada (≈ diag(0.058^2) y off-diag ≈ -0.003):\n")

Comparación con la forma reportada (≈ diag(0.058^2) y off-diag ≈ -0.003):

cat("0.058^2 =", round(0.058^2, 6), "\n")

0.058^2 = 0.003364 

3.6.2.2 Ejemplo: Serie de reclutamiento

rec.yw = ar.yw(rec, order=2)
rec.yw$x.mean

[1] 62.26278

rec.yw$ar

[1]  1.3315874 -0.4445447

sqrt(diag(rec.yw$asy.var.coef))

[1] 0.04222637 0.04222637

rec.yw$var.pred

[1] 94.79912

Predicción a 24 pasos:

rec.pr = predict(rec.yw, n.ahead=24)
ts.plot(rec, rec.pr$pred, col=1:2)
lines(rec.pr$pred + rec.pr$se, col=4, lty=2)
lines(rec.pr$pred - rec.pr$se, col=4, lty=2)

3.6.2.3 Ejemplo 3.29 MA(1)

Modelo:

\[ x_{t}=w_{t}+\theta w_{t-1}, \quad |\theta|<1. \]

El estimador de \(\theta\) proviene de:

\[ \hat{\rho}(1)=\frac{\hat{\theta}}{1+\hat{\theta}^{2}}. \]

Solución invertible:

\[ \hat{\theta}=\frac{1-\sqrt{1-4\hat{\rho}(1)^{2}}}{2\hat{\rho}(1)}. \]

Simulación en R:

set.seed(2)
true_rho = 0.9 / (1 + 0.9^2) 
ma1 = arima.sim(list(order = c(0,0,1), ma = 0.9), n = 50)
acf(ma1, plot=FALSE)[1]

Autocorrelations of series 'ma1', by lag

    1 
0.507 

3.7 Estimación de máxima verosimilitud y mínimos cuadrados

3.7.1 Caso AR(1)

Sea el modelo causal AR(1):

\[ \begin{equation*} x_{t}=\mu+\phi\left(x_{t-1}-\mu\right)+w_{t} \end{equation*} \]

con \(|\phi|<1\) y \(w_{t}\sim \text{iid } \mathrm{N}(0,\sigma_w^2)\). Dado \(x_{1\:n}\), la verosimilitud buscada es

\[ L\left(\mu, \phi, \sigma_{w}^{2}\right)=f\left(x_{1}, x_{2}, \ldots, x_{n} \mid \mu, \phi, \sigma_{w}^{2}\right). \]

Para AR(1), usando la factorización condicional:

\[ L\left(\mu, \phi, \sigma_{w}^{2}\right)=f\left(x_{1}\right) f\left(x_{2} \mid x_{1}\right) \cdots f\left(x_{n} \mid x_{n-1}\right), \]

y como \(x_t|x_{t-1}\sim \mathrm{N}!\left(\mu+\phi(x_{t-1}-\mu),\sigma_w^2\right)\),

\[ f\left(x_{t} \mid x_{t-1}\right)=f_{w}\left[\left(x_{t}-\mu\right)-\phi\left(x_{t-1}-\mu\right)\right], \]

donde \(f_w\) es la densidad normal de media \(0\) y varianza \(\sigma_w^2\). Entonces:

\[ L\left(\mu, \phi, \sigma_{w}\right)=f\left(x_{1}\right) \prod_{t=2}^{n} f_{w}\left[\left(x_{t}-\mu\right)-\phi\left(x_{t-1}-\mu\right)\right] . \]

Usando la representación causal

\[ x_{1}=\mu+\sum_{j=0}^{\infty} \phi^{j} w_{1-j}, \]

se obtiene \(x_1\sim \mathrm{N}\left(\mu,\ \sigma_w^2/(1-\phi^2)\right)\) y, por tanto, la verosimilitud incondicional es

\[ \begin{equation*} L\left(\mu, \phi, \sigma_{w}^{2}\right)=\left(2 \pi \sigma_{w}^{2}\right)^{-n / 2}\left(1-\phi^{2}\right)^{1 / 2} \exp \left[-\frac{S(\mu, \phi)}{2 \sigma_{w}^{2}}\right] \end{equation*} \]

donde

\[ \begin{equation*} S(\mu, \phi)=\left(1-\phi^{2}\right)\left(x_{1}-\mu\right)^{2}+\sum_{t=2}^{n}\left[\left(x_{t}-\mu\right)-\phi\left(x_{t-1}-\mu\right)\right]^{2}. \end{equation*} \]

Derivando respecto de \(\sigma_w^2\) en el logaritmo de (3.107) y anulando, para valores dados de \((\mu,\phi)\) se maximiza con

\[ \begin{equation*} \hat{\sigma}_{w}^{2}=n^{-1} S(\hat{\mu}, \hat{\phi}). \end{equation*} \]

Sustituyendo y concentrando, el criterio a minimizar en \((\mu,\phi)\) es

\[ \begin{equation*} l(\mu, \phi)=\log \left[n^{-1} S(\mu, \phi)\right]-n^{-1} \log \left(1-\phi^{2}\right). \end{equation*} \]

3.7.2 Verosimilitud condicional y conexión con mínimos cuadrados

Condicionando en \(x_1\):

\[ \begin{align*} L\left(\mu, \phi, \sigma_{w}^{2} \mid x_{1}\right) & =\prod_{t=2}^{n} f_{w}\left[\left(x_{t}-\mu\right)-\phi\left(x_{t-1}-\mu\right)\right] \\ & =\left(2 \pi \sigma_{w}^{2}\right)^{-(n-1) / 2} \exp \left[-\frac{S_{c}(\mu, \phi)}{2 \sigma_{w}^{2}}\right], \end{align*} \]

con

\[ \begin{equation*} S_{c}(\mu, \phi)=\sum_{t=2}^{n}\left[\left(x_{t}-\mu\right)-\phi\left(x_{t-1}-\mu\right)\right]^{2}. \end{equation*} \]

El EMV condicional de \(\sigma_w^2\) es

\[ \begin{equation*} \hat{\sigma}_{w}^{2}=S_{c}(\hat{\mu}, \hat{\phi}) /(n-1), \end{equation*} \]

y los EMV condicionales de \((\mu,\phi)\) minimizan \(S_c(\mu,\phi)\). Sea \(\alpha=\mu(1-\phi)\):

\[ \begin{equation*} S_{c}(\mu, \phi)=\sum_{t=2}^{n}\left[x_{t}-\left(\alpha+\phi x_{t-1}\right)\right]^{2}, \end{equation*} \]

lo que es un problema de regresión lineal. De mínimos cuadrados:

\[ \begin{gather*} \hat{\mu}=\frac{\bar{x}_{(2)}-\hat{\phi} \bar{x}_{(1)}}{1-\hat{\phi}}, \\ \hat{\phi}=\frac{\sum_{t=2}^{n}\left(x_{t}-\bar{x}_{(2)}\right)\left(x_{t-1}-\bar{x}_{(1)}\right)}{\sum_{t=2}^{n}\left(x_{t-1}-\bar{x}_{(1)}\right)^{2}}. \end{gather*} \]

En consecuencia, \(\hat{\mu}\approx \bar{x}\) y \(\hat{\phi}\approx \hat{\rho}(1)\); Yule–Walker y MC condicional son cercanos (difieren en el tratamiento de \(x_1\) y \(x_n\)). Puede ajustarse dividiendo por \((n-3)\) en lugar de \((n-1)\) para equivaler a MC.