Другое введение в анализ временных рядов

Ключевые концепции для понимания временных рядов и выбора лучшей прогностической модели

Введение. Что такое временной ряд?

Временной ряд — это упорядоченная последовательность мер величины во времени, где с каждой мерой связано определенное время. Наряду с этим постом мы собираемся использовать набор данных из задачи Kaggle Прогнозирование спроса на товары в магазине. Среднесуточные продажи этого набора данных имеют следующую закономерность.

df.set_index('date').plot()
plt.title('Store Item Demand Forecasting')

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

1. Статистический анализ — понимание поведения временного ряда

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

1.1 Характеристика временного ряда — стационарный анализ

Временной ряд Y(t) называется стационарным, если его статистические свойства не меняются во времени. Эта концепция важна, потому что будущие значения будут иметь такое же поведение. Предыдущее не означает, что данные будут предсказуемыми и наоборот (например, временной ряд с трендом или сезонностью не является стационарным), но в худшем сценарии мы можем оценить ожидаемое значение с определенной дисперсией при любой заданный момент времени.

Следовательно,набор данных Y(t) для прогнозирования спроса на товары в магазине не является стационарным. Хотя иногда временной ряд не является стационарным, можно выполнить некоторые преобразования, чтобы преобразовать его в стационарный, например, путем устранения эффекта тренда с дифференциацией n-го порядка. Для ряда Y (t) дифференцирование n-го порядка определяется как

Y’(t) = Y(t) — Y(t-1) — … — Y(t-n)

Делая дифференциацию 1-го порядка в последнем году нашего набора данных, мы получаем следующее

df.sales=df.sales.diff()
df.plot()
plt.title('Rolling diff Demand')

Несмотря на это, иногда данные остаются нестационарными после дифференцирования n-го порядка, потому что рост не является линейным. Например, если мы увидим весь временной ряд Y(t), а не только последний год, мы увидим, что он не является стационарным.

В таких случаях мы можем попытаться заранее выполнить логарифмическое или экспоненциальноепреобразование Y(t). С помощью этого преобразования мы наконец можем увидеть, что наш временной ряд Y (t) является стационарным.

Y’(t) = log(Y(t))

Y’(t) = exp(Y(t))

df.sales=np.log(df.sales).diff()
df.set_index('date').plot()
plt.title('Rolling diff Demand')

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

1.2 Автокорреляция и частичная автокорреляция (ACF PACF)

Чтобы предсказать будущие значения, мы не только хотим понять временную эволюцию статистических свойств Y(t), но также понять, насколько связаны будущие значения с прошлыми значениями. Для этого мы используем ACF и PACF.

ACF — это корреляция временного ряда Y(t) с самим собой в предыдущий момент времени ряда. Вычисляя корреляцию с k прошлыми значениями ряда, мы можем проверить, какие из k предыдущих значений коррелируют с текущим значением

ACF = корреляция ( Y (t), Y (t-k))

from statsmodels.graphics.tsaplots import plot_acf
plot_acf(df.sales)

Аналогичным образом PACF выполняет корреляцию Y(t) с прошлым значением Y(t-k), но удаляет влияние промежуточных значений: Y(t-1), Y(t-2), … Y(t-k+1):

PACF = корреляция (Y (t), Y (t-k) | Y (t-1), Y (t-2), … Y (t-k + 1))

from statsmodels.graphics.tsaplots import plot_pacf
plot_pacf(df.sales)

Основное различие между ними состоит в том, что АКФ с Y(t-k) учитывает не только свою зависимость от Y(t), но и все промежуточное звено между Y(t) и Y(t-k) косвенно. В то время как PACF с Y(t-k) учитывает только взаимосвязь с Y(t) и Y(t-k).

1.3 Шум и остатки

Иногда бывает так, что дисперсия данных очень высока, чтобы увидеть какую-либо закономерность. Такой шаблон называется шумной серией. Нет немедленного решения этой проблемы. Одна из возможностей состоит в том, чтобы агрегировать данные на более высоком уровне данных (например, на уровне дня вместо уровня часа). Другим распространенным решением для обработки данных является вычисление скользящего среднего ряда.

df.rolling(window=7).mean().plot()

В нашем случае это не требуется, потому что мы работаем с ежедневными данными и хотим прогнозировать ежедневные продажи. Если бы наши данные были почасовыми продажами, было бы интересно объединить их в дневную сумму.

1.4 Сезонность

Обычно значение Y(t) не только связано с ближайшими предыдущими значениями, но и может иметь длительный период повторения. Такое поведение называется сезонностью. Более того, обычно данные имеют не только один сезонный вклад, но и многие из них одновременно (например, годовая сезонность, месячная сезонность и т. д.).

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

Предполагая, что сезонность имеет постоянный период, мы можем вычислить эти периоды каждого вклада сезонности, используя разложение ряда Фурье. Во-первых, удобно убрать эффект трендаиз данных. Если данные имеют линейный тренд, мы можем сделать это с помощью линейной регрессии.

trend = LinearRegression(
    ).fit(df.index.values.reshape([-1,1]),df.sales
    ).predict(df.index.values.reshape([-1,1]))
df.sales = df.sales-trend

Имея ряд Y(t) без тренда, мы можем провести анализ разложения Фурье. При анализе мы можем увидеть периодичность сезонности.

from scipy.fftpack import fft
plt.plot(fft(df.sales.values))

1.5 Анализ результатов

При анализе мы увидели следующие свойства ряда продаж:

• Стационарное поведение при логарифмическом росте.
• Периодичность сезонности 8,8 мес.
• Последний месяц сильно коррелирует с текущим значением.
• Рецидив 7 дней.

2. Прогнозные модели — прогнозируйте эволюцию временного ряда.

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

2.1 Авторегрессионные модели

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

Y(t) = Y(t-1)*a_t-1 + Y(t-2)*a_t-2 +… Y(t-p)*a_t-p

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

Y(t) = Y(t-1)*a_t-1 + Y(t-2)*a_t-2 +… Y(t-p)*a_t-p +

+ ошибка(t-1)*b_t-1 + ошибка(t-2)*b_t-2 + … + ошибка(t-q)*b_t-q

Где ошибка (t) = Y (t) — Y (t-1). Эта модель требует, чтобы временной ряд был стационарным, поэтому она автоматически выполняет дифференциацию ряда n-го порядка, чтобы при необходимости преобразовать ряд в стационарный (важно отметить, что она выполняет только дифференцирование, не экспоненциальное логарифмическое преобразование).

Y(t) = Y(t) — Y(t-1) — Y(t-2) — … — Y(t-d)

Значения p, d и q из предыдущих формул являются параметрами, необходимыми для настройки ARIMA. p и q могут быть получены из последних значений частичной автокорреляции и корреляции соответственно с высокой корреляцией. Значение d — это просто порядок дифференцирования, при котором ряд остается стационарным. Однако ARIMA не справляется с сезонностью. Вот почему существует расширение алгоритма под названием SARIMA, которое может учитывать сезонность.

import statsmodels.api as sm
model = sm.tsa.statespace.SARIMAX(df_train.sales,
    trend='t', freq='D', order=(30,1,9)
).fit(disp=False, solver='powell')
model.predict(365)

2.2 Регрессионные модели

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

Решением для обработки нелинейной авторегрессии является использование модели нелинейной регрессии, способной суммировать значения из окна времени, чтобы делать прогнозы на будущее. Одной из таких моделей является сеть LSTM (конкретная рекуррентная нейронная сеть). Вы можете добавить сезонный эффект к данным, преобразовав временной ряд в табличный набор данных, подобный этому.

weekday_df = pd.get_dummies(
    df.date.dt.dayofweek.to_numpy(), prefix='weekday')
df = pd.merge(weekday_df, df, how='inner', on = 'date')

Чтобы взять ближайшие прошлые значения в качестве авторегрессионных моделей, мы должны построить стек движущихся окон ближайших значений. Модель LSTM создается путем объединения значений окна в массив, но другие модели регрессии (например, RandomForest) создаются по-другому.

x = np.zeros((len(df)-seq_size, seq_size, df.shape[-1]-1))
y = np.zeros((len(df)-seq_size, 1))
    
for i in range(len(df)-seq_size):
    x[i,:,:]=df[i:(i+seq_size),0:-1]
    y[i,:]=df[i+seq_size, -1:]

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

from keras.models import Sequential
from keras.layers import LSTM, Flatten,Dropout,Dense
inputshape = (X_train.shape[1], X_train.shape[2])
model = Sequential()
model.add(LSTM(64, activation='tanh', return_sequences=False))
model.add(Dense(500))
model.add(Dense(1))
model.compile(optimizer='adam', loss = 'mse', metrics = ['accuracy'])
model.fit(X_train, y_train, verbose=1, epochs = 10, batch_size = 1,shuffle = False)

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

2.3 Модели ряда Фурье

Авторегрессионные модели хорошо предсказывают временные ряды, где будущие значения связаны с близкими предыдущими значениями (скользящее среднее) с определенной периодичностью, как в предыдущем примере. Однако эти модели не могут предсказать временные ряды, которые не имеют вклада скользящего среднего и с несколькими сезонностями одновременно (ежедневно, ежемесячно, триместрально, еженедельно,…).

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

Для нашего текущего набора данных мы собираемся использовать сезонность по умолчанию (годовую и еженедельную), а также добавим добавим сезонное поведение 8,8 месяцев, которое мы видели в анализе Фурье.

from prophet import Prophet
model = Prophet(
    seasonality_mode='multiplicative',
    ).add_seasonality(name='monthly', period=260, fourier_order=5
                     ).fit(df_train)
forecast = model.predict(model.make_future_dataframe(periods=365))

Как мы видим в коде, мы устанавливаем сезонность как «мультипликативную». В авторегрессионной модели мы увидели эффект сезонности как дополнительный вклад во временной ряд.

Y(t) = Y_тренд(t) + Y_сезон(t)

Однако иногда влияние сезонности меняет тренд не на добавку, а на коэффициент. Пророк способен соответствовать этим двум видам поведения

Y(t) = Y_тренд(t) * Y_сезон(t)

С Prophet мы также можем видеть компонентную декомпозицию ряда, поскольку мы можем видеть различные компоненты временного ряда.

3. Выводы

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

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

Надеюсь, вам понравилось :)