5. ARIMA: integración, ajuste y predicción con ventana móvil

Nivel 9: Series Temporales - ARIMA y GARCH

El modelo autorregresivo integrado de media móvil, o ARIMA, es el modelo estadístico estándar para analizar y predecir series temporales. Sus siglas describen exactamente sus tres piezas:

• AR: autorregresión. Usa la relación de dependencia entre una observación y cierto número de observaciones retardadas.
• I: integrado. Aplica diferenciación a las observaciones originales (restar a cada observación la del paso anterior) para hacer la serie estacionaria.
• MA: media móvil. Usa la dependencia entre una observación y los errores residuales de un modelo de media móvil aplicado a observaciones retardadas.

Sus tres parámetros:

• p: número de observaciones retardadas incluidas, también llamado orden de retardo.
• d: número de veces que se diferencian las observaciones originales, o grado de diferenciación.
• q: tamaño de la ventana de media móvil, u orden de media móvil.

Cuando d es cero tenemos un ARMA(p,q), que combina ambos efectos:

# ARMA(p,q):
# x[t] = alpha_1*x[t-1] + ... + alpha_p*x[t-p]
#        + w[t] + beta_1*w[t-1] + ... + beta_q*w[t-q]

Un ARMA suele necesitar menos parámetros que un AR(p) o un MA(q) por separado para el mismo ajuste. Y una advertencia que será el puente hacia los dos últimos modulos: un ARMA no tiene en cuenta la agrupación de volatilidad, un fenómeno empírico clave de muchas series financieras.

Ajuste basico (usa la clase moderna, no la obsoleta):

from statsmodels.tsa.arima.model import ARIMA

model = ARIMA(Re, order=(5, 1, 0))
model_fit = model.fit()
print(model_fit.summary())

residuals = pd.DataFrame(model_fit.resid)
residuals.plot()
residuals.plot(kind='kde')
print(residuals.describe())

Una media de residuos cercana a cero indica ausencia de sesgo; una media positiva o negativa indica un sesgo que puede corregirse directamente antes de usar las previsiones.

Ojo con un detalle que el código del curso deja pasar: Re ya son retornos logarítmicos, es decir, la serie ya está diferenciada una vez. Poner d = 1 encima significa sobrediferenciar, lo que introduce una estructura MA artificial en los residuos e infla la varianza. Sobre retornos, d = 0 es casi siempre lo correcto.

Previsión ARIMA con ventana móvil. Igual que con el AR, queremos reconstruir el modelo cada día. Aquí si se reajusta desde cero en cada paso:

from tqdm import tqdm

X = Re.values
size = int(len(X) * 0.66)
train, test = X[0:size], X[size:len(X)]
history = [x for x in train]
predictions = []

for t in tqdm(range(len(test)), desc='Procesando'):
    model = ARIMA(history, order=(5, 1, 0))
    model_fit = model.fit()
    output = model_fit.forecast()
    yhat = output[0]
    predictions.append(yhat)
    history.append(test[t])       # solo despues de predecir

rmse = np.sqrt(mean_squared_error(test, predictions))
print('Test RMSE: %.3f' % rmse)

Esto es honesto pero costoso: son tantos ajustes de máxima verosimilitud como días de prueba.

Busqueda en rejilla de hiperparametros. En lugar de fijar (5,1,0) a ojo, evaluamos combinaciones y nos quedamos con la de menor RMSE fuera de muestra:

def evaluate_arima_model(X, arima_order):
    train_size = int(len(X) * 0.66)
    train, test = X[0:train_size], X[train_size:]
    history = [x for x in train]
    predictions = list()
    for t in range(len(test)):
        model_fit = ARIMA(history, order=arima_order).fit()
        predictions.append(model_fit.forecast()[0])
        history.append(test[t])
    return np.sqrt(mean_squared_error(test, predictions))

def evaluate_models(dataset, p_values, d_values, q_values):
    dataset = dataset.astype('float32')
    best_score, best_cfg = float('inf'), None
    for p in p_values:
        for d in d_values:
            for q in q_values:
                order = (p, d, q)
                try:
                    rmse = evaluate_arima_model(dataset, order)
                    if rmse < best_score:
                        best_score, best_cfg = rmse, order
                    print('ARIMA%s RMSE=%.3f' % (order, rmse))
                except:
                    continue
    print('Mejor ARIMA%s RMSE=%.3f' % (best_cfg, best_score))

evaluate_models(Re.values, [0, 1, 2], [0, 1], [0, 1])

El except: continue desnudo es cómodo (muchas combinaciones no convergen) pero peligroso: silencia también errores de programación. En producción, captura la excepción y regístrala.