2- Эконометрические и статистические модели во временных рядах

Обновление: эта статья является частью серии, в которой я исследовал временные ряды. Ознакомьтесь с полной серией: Часть 1, Часть 2, Часть 3.

Вступление

У меня есть среднемесячные температуры во Флориде с 1743 года, о которых я говорил в предыдущей статье. В этой статье я сначала объясню в основном некоторые модели и их определение, а также покажу вам основные операции во временных рядах с использованием таких моделей, как ARIMA, SARIMA с этим набором данных.

Не объясняй мне, покажи код

Вы можете получить доступ к репозиторию GitHub здесь.

AR (p): авторегрессия

Оценка производится путем линейной комбинации наблюдений за предыдущие временные шаги.
Она подходит для одномерных временных рядов без тренда или сезонности.
p: временная задержка. Если p = 1, это означает, что модель была создана с предыдущим временным шагом.

MA (q): скользящая средняя.

Оценка выполняется с помощью линейной комбинации ошибок, полученных на предыдущих временных шагах.
Она подходит для одномерных временных рядов без тренда или сезонности.
q: число временной задержки.

ARMA (p, q) = AR (p) + MA (q) Авторегрессионная скользящая средняя

Сочетает методы AR и MA.
Оценка выполняется с помощью линейной комбинации исторических значений и прошлых ошибок.
Она подходит для одномерных временных рядов без тренда или сезонности.
p и q - время- номера задержки. p для модели AR q для модели MA.

ARIMA (p, d, q): (интегрированное скользящее среднее с авторегрессией)

Модели ARIMA используются для анализа данных временных рядов и составления прогнозов на будущее. Обычное употребление - арима (p, d, q). p используется для авторегрессионной (ar (p)) части, а d указывает, во сколько раз различаются ряды данных. (новый ряд, полученный вычитанием каждых данных из предыдущих данных) Наконец, q используется для части скользящего среднего (ma (q)). В этом случае Arima (1,0,0) и ar (1) совпадают.
arima (0,0,1) и ma (1) совпадают.

Если мы хотим резюмировать шаг за шагом:

Оценка производится путем линейной комбинации наблюдений и ошибок, которые отличались от предыдущих временных шагов.
Он подходит для одномерных, трендовых, но не сезонных данных.
p: номер задержки фактического значения (степень авторегрессии)
- Если p = 2, то это в моделях yt-1 и yt-2.
d: количество различных операций (степень различия, I)
q: номер задержки ошибки (степень скользящего среднего)
Что значит быть замеченным?
- Например, это вычитание значений сегодняшнего и предыдущего дня друг из друга.
Что делает процесс обнаружения, почему нас замечают?
- Чтобы стабилизировать серию.
Что значит стабилизировать ряд?
- Исключить ситуацию, когда статистические свойства ряда меняются со временем.
Как можно понять, процесс определения разницы здесь создал возможность делать прогнозы в сериях, которые сейчас в тренде.

SARIMA (p, d, q) (P, D, Q) m: (Сезонная авторегрессионная интегрированная скользящая средняя)

Это метод, особенно используемый при анализе временных рядов и используемый для объяснения и прогнозирования данных. В дополнение к модели Арима может быть на 3 параметра больше, чем модель очистки для использования в сезонных данных. Эти параметры представляют собой автокорреляцию, разность и т. Д. Сезонных данных. ценности.

Если мы хотим резюмировать шаг за шагом:

ARIMA + - это сезонность
Его можно использовать в одномерных рядах, которые включают тренд и сезонность.
p, d, q - параметры из ARIMA. Элементы тренда *. ARIMA смогла смоделировать тенденцию.
p: номер задержки фактического значения (степень авторегрессии)
- Если P = 2, то это в моделях yt-1 и yt-2.
d: количество разных транзакций (степень различия)
q: номер задержки ошибки (степень скользящего среднего)
Если q = 2, то это в моделях et-1 и et-2.
P, D, Q числа сезонных лагов. Элементы сезона.
m количество временных шагов за один сезонный период. Он выражает структуру сезонности.

0. Подготовка данных

Прежде всего, мы разделяем наши данные как поезд / тест для нашей модели Arima. Как и в предыдущей статье, при определении данных, найденных до 1995 года, в виде набора поездов.

train = florida[:"1994-12-01"]
len(train)
test = florida["1995-01-01":]
len(test)

1. Модель Арима

В качестве первого шага мы настраиваем базовую модель Arima и смотрим на нашу тестовую ошибку, не применяя никаких операций настройки.

arima_model = ARIMA(train, order = (1,0,1)).fit(disp=0)
arima_model.summary()

y_pred= arima_model.forecast(225)[0]
mean_absolute_error(test, y_pred)
# 4.39981

В визуализации я визуализирую после 1985 года для удобного и четкого изображения.

train["1985":].plot(legend=True, label = 'TRAIN')
test.plot(legend=True, label = 'TEST', figsize = (6,4))
pd.Series(y_pred, index=test.index).plot(legend=True, label = 'Prediction')
plt.title("Train, Test and Predicted Test")
plt.show()

Как видите, наша догадка не кажется удачной. Итак, приступим к настройке.

Есть два обычно используемых метода для оптимизации этих статистических моделей. Во-первых, определение оценки модели на основе диаграмм ACF и PACF. Во-вторых, определение оценки модели на основе статистики AIC.

1.1 Определение класса модели по графикам ACF и PACF

def acf_pacf(y, lags=30):
plt.figure(figsize=(12, 7))
layout = (2, 2)
ts_ax = plt.subplot2grid(layout, (0, 0), colspan=2)
acf_ax = plt.subplot2grid(layout, (1, 0))
pacf_ax = plt.subplot2grid(layout, (1, 1))
y.plot(ax=ts_ax)
# Durağanlık testi (HO: Series is not Stationary. H1: Series is Stationary.)
p_value = sm.tsa.stattools.adfuller(y)[1]
ts_ax.set_title(‘Time Series Analysis Plots\n Dickey-Fuller: p={0:.5f}’.format(p_value))
smt.graphics.plot_acf(y, lags=lags, ax=acf_ax)
smt.graphics.plot_pacf(y, lags=lags, ax=pacf_ax)
plt.tight_layout()
plt.show()
acf_pacf(florida)

Здесь мы пытаемся дать лучшее значение для значений p и q на графике.

Здесь необходимо посмотреть на корреляции между графиком ACF и значениями предыдущих периодов во временном ряду. Здесь мы смотрим на транзакции со значением лага 30.

Итак, мы смотрим на yt-1, yt-2,…. расширяется до yt-30.

Автокорреляция не будет значительной в точках, оставшихся в голубой области на графике, у нас много значимых, но мало бессмысленных автокорреляций.

Если мы хотим кратко подвести итог:

- Это означает, что если ширина ACF «УМЕНЬШАЕТСЯ» в соответствии с задержками, а PACF p - «CUT OFF »после шаблона задержки AR (p).

- Если после задержки для ширины ACF q установлено значение «ВЫКЛЮЧЕНО», а для ширины PACF «УМЕНЬШАЕТСЯ» по задержке это означает, что это паттерн MA (q).

- Если ширина ACF и PACF уменьшается в соответствии с лагами, это означает, что это модель ARMA.

df_diff = florida.diff()
df_diff.dropna(inplace=True)
acf_pacf(df_diff)

Мы легко можем увидеть сезонность в автокорреляции. Мы не можем увидеть оценку с помощью этой визуальной техники, которая могла бы помочь определить оценку.

2 Определение рейтинга модели согласно статистике AIC и BIC

Теперь мы пытаемся найти лучшие параметры, пробуя все комбинации в диапазоне, который мы определили, чтобы найти значения, которые могут дать нам лучший результат.

# Generation of combinations of p and q

p = d = q = range(0, 4)
pdq = list(itertools.product(p, d, q))

def arima_optimizer_aic(train, orders):
    best_aic, best_params = float("inf"), None

    for order in orders:
        try:
            arma_model_result = ARIMA(train, order).fit(disp=0)
            aic = arma_model_result.aic
            if aic < best_aic:
                best_aic, best_params = aic, order
            print('ARIMA%s AIC=%.2f' % (order, aic))
        except:
            continue
    print('Best ARIMA%s AIC=%.2f' % (best_params, best_aic))
    return best_params

best_params_aic = arima_optimizer_aic(train, pdq)

Добавляя лучшие значения, которые мы нашли, в нашу модель ARIMA, мы создаем нашу настроенную модель.

arima_model = ARIMA(train, best_params_aic).fit(disp=0)
y_pred = arima_model.forecast(225)[0]
mean_absolute_error(test, y_pred)
# 1.4470

Когда мы визуализируем нашу настроенную модель, мы можем сказать, что у нас есть очень успешный результат модели по сравнению с предыдущим.

2. Модель SARIMA

Прежде всего, мы выполнили этот шаг, так как нам нужно будет создать проверку в нашем наборе данных.

from statsmodels.tsa.statespace.sarimax import SARIMAX

train = florida[:"1994-12-01"]
len(train)
test = florida["1995-01-01":]
len(test)
val = train["1991-01-01":]
len(val)

2.1 Базовая модель

Мы сделаем это с помощью функции SARIMAX.

model = SARIMAX(train, order=(1,0,1), seasonal_order=(0,0,0,12))
sarima_model = model.fit(disp=0)

2.2 Ошибка проверки

pred = sarima_model.get_prediction(start = pd.to_datetime('1991-01-01'),dynamic=False)
pred_ci = pred.conf_int()

y_pred = pred.predicted_mean
mean_absolute_error(val, y_pred)
# 1.9192

2.3 Визуализация проверочного прогноза

Наша ошибка проверки выглядит довольно хорошо по сравнению с базовой моделью.

2.4 Ошибка теста

y_pred_test = sarima_model.get_forecast(steps=225)
pred_ci = y_pred_test.conf_int()
y_pred = y_pred_test.predicted_mean
mean_absolute_error(test, y_pred)
# 17.0051

Мы видим, что наша реальная ошибка теперь относительно высока. Если мы думаем об этом как о фактической оценке, низкая ошибка проверки, которую мы видели на предыдущем шаге, может указывать на переобучение.

2.5 Визуализация прогнозов

Мы видим, насколько плох график.

2.6 Настройка модели

p = d = q = range(0, 2)
pdq = list(itertools.product(p, d, q))
seasonal_pdq = [(x[0], x[1], x[2], 12) for x in list(itertools.product(p, d, q))]

def sarima_optimizer_aic(train, pdq, seasonal_pdq):
    best_aic, best_order, best_seasonal_order = float("inf"), float("inf"), None
    for param in pdq:
        for param_seasonal in seasonal_pdq:
            try:
                sarimax_model = SARIMAX(train, order=param, seasonal_order=param_seasonal)
                results = sarimax_model.fit(disp=0)
                aic = results.aic
                if aic < best_aic:
                    best_aic, best_order, best_seasonal_order = aic, param, param_seasonal
                print('SARIMA{}x{}12 - AIC:{}'.format(param, param_seasonal, aic))
            except:
                continue
    print('SARIMA{}x{}12 - AIC:{}'.format(best_order, best_seasonal_order, best_aic))
    return best_order, best_seasonal_order

best_order, best_seasonal_order = sarima_optimizer_aic(train, pdq, seasonal_pdq)

Здесь мы постарались сохранить небольшой диапазон (2) из-за избытка переменных. Вы можете увеличить дальность. Позже мы написали 12, потому что знаем сезонный период.

2.7 Окончательная модель и ее ошибка тестирования

model = SARIMAX(train, order=best_order, seasonal_order=best_seasonal_order)
sarima_final_model = model.fit(disp=0)
###############################
### Final Model Test Error ###
###############################

y_pred_test = sarima_final_model.get_forecast(steps=225)
pred_ci = y_pred_test.conf_int()
y_pred = y_pred_test.predicted_mean
mean_absolute_error(test, y_pred)
# 1.013578872841977

2.8 Визуализация окончательной модели

2.9 Оптимизатор SARIMA на основе MAE

Мы оптимизировали модель в соответствии с AIC, но мы также могли сделать это в соответствии с MAE. Для этого я хотел показать это, открыв здесь дополнительный раздел.

def fit_model_sarima(train, val, pdq, seasonal_pdq):
    sarima_model = SARIMAX(train, order=pdq, seasonal_order=seasonal_pdq).fit(disp=0)
    y_pred_val = sarima_model.get_forecast(steps=48)
    y_pred = y_pred_val.predicted_mean
    return mean_absolute_error(val, y_pred)

fit_model_sarima(train, val, (0, 1, 0), (0, 0, 0, 12))

p = d = q = range(0, 2)
pdq = list(itertools.product(p, d, q))
seasonal_pdq = [(x[0], x[1], x[2], 12) for x in list(itertools.product(p, d, q))]


def sarima_optimizer_mae(train, val, pdq, seasonal_pdq):
    best_mae, best_order, best_seasonal_order = float("inf"), float("inf"), None
    for param in pdq:
        for param_seasonal in seasonal_pdq:
            try:
                mae = fit_model_sarima(train, val, param, param_seasonal)
                if mae < best_mae:
                    best_mae, best_order, best_seasonal_order = mae, param, param_seasonal
                print('SARIMA{}x{}12 - MAE:{}'.format(param, param_seasonal, mae))
            except:
                continue
    print('SARIMA{}x{}12 - MAE:{}'.format(best_order, best_seasonal_order, best_mae))
    return best_order, best_seasonal_order


best_order, best_seasonal_order = sarima_optimizer_mae(train, val, pdq, seasonal_pdq)

model = SARIMAX(train, order=best_order, seasonal_order=best_seasonal_order)
sarima_final_model = model.fit(disp=0)

y_pred_test = sarima_final_model.get_forecast(steps=225)
pred_ci = y_pred_test.conf_int()
y_pred = y_pred_test.predicted_mean
mean_absolute_error(test, y_pred)
# 0.92