Метод декомпозиции временных рядов с ускорением градиента
TLDR:ThymeBoost сочетаеттрадиционный процесс декомпозиции с повышением градиента, чтобы обеспечить гибкую структуру временных рядов для анализа тенденций/сезонности/экзогенных факторов и прогнозирования. прочь.
Весь код живет здесь: ThymeBoost Github
Вот новейшая статья из серии: Конкуренция M4 с ThymeBoost.
Мотивирующий пример
Традиционная декомпозиция временных рядов обычно включает в себя последовательность шагов:
- Приблизительный тренд/уровень
- Отклонение от приблизительной сезонности
- Удаление тренда и десезонализация для приблизительного соответствия другим факторам (экзогенным)
- Добавьте тренд, сезонность и экзогенность для результатов.
В зависимости от определенных алгоритмов на каждом этапе нам может потребоваться сначала изменить порядок и десезонализировать. В любом случае, на первый взгляд, мы можем убедить себя, что этот процесс… хороший.
На самом деле, он может неплохо работать в простой симуляции. Давайте создадим базовый временной ряд тренда + сезонности, чтобы показать, как работает этот процесс.
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
sns.set_style("darkgrid")
#Here we will just create a simple series
np.random.seed(100)
trend = np.linspace(1, 50, 100) + 50
seasonality = ((np.cos(np.arange(1, 101))*10))
y = trend + seasonality
#let's plot it
plt.plot(y)
plt.show()
А теперь, чтобы сделать некоторую базовую декомпозицию.
Следуя описанному выше процессу, мы сначала подгоним простой линейный тренд, чтобы аппроксимировать компонент тренда:
import statsmodels.api as sm
#A simple input matrix for a deterministic trend component
X = np.array([range(1, len(y) + 1),
np.ones(len(y))]).T
mod = sm.OLS(y, X)
res = mod.fit()
fitted_trend = res.predict(X)
plt.plot(trend)
plt.plot(fitted_trend, linestyle='dashed')
plt.show()
Эта линия тренда на высоте! Мы вычтем подобранный тренд из наших фактических значений и используем этот сигнал без тренда для аппроксимации сезонной составляющей:
detrended_series = y - fitted_trend #Set the seasonal period seasonal_period = 25 #Get the average of every 25 data points to use as seasonality avg_season = [np.mean(detrended_series[i::seasonal_period], axis=0) for i in range(seasonal_period)] avg_season = np.array(avg_season) fitted_seasonality = np.resize(avg_season, len(y)) plt.plot(fitted_seasonality) plt.plot(seasonality) plt.show()
Сезонная составляющая тоже хороша! Хотя мы видим небольшие отклонения, в целом две линии почти идентичны. Вы могли заметить, что в сюжетах нет легенды. Это связано с тем, что наш сигнал изначально не имеет шума, поэтому разложение и ряд практически одинаковы. Процесс разложения довольно хорошо себя ведет, не требуя легенды… пока.
Давайте соберем все это вместе и посмотрим на наши полные результаты:
fitted = fitted_trend + fitted_seasonality plt.plot(fitted, linestyle='dashed') plt.plot(y, linestyle='dashed') plt.show()
Окончательные результаты разложения весьма убедительны. Этот процесс дает нам полный контроль над различными компонентами. Допустим, нам нужно относиться к сезонности иначе, чем к тренду, возможно, придавая некоторым образцам больший или меньший вес:
Можем.
Кроме того, если мы добавим экзогенный фактор, мы сможем аппроксимировать этот компонент любой моделью, которую захотим. Если у нас есть много данных и сложных функций, мы могли бы просто вычесть тренд и сезонный компонент, а затем передать эти остатки в XGBoost.
Но для этого примера давайте добавим простой эффект «праздника» к нашему текущему ряду и сопоставим его с линейной регрессией. Во-первых, мы восстановим наш смоделированный ряд, но теперь добавим этот экзогенный фактор.
np.random.seed(100) trend = np.linspace(1, 50, 100) + 50 seasonality = ((np.cos(np.arange(1, 101))*10)) exogenous = np.random.randint(low=0, high=2, size=len(trend)) y = trend + seasonality + exogenous * 20 plt.plot(y) plt.show()
Мы просто создали случайный ряд из 1 и 0 в качестве нашего «внешнего» и дали ему эффект добавления 20 к нашему ряду, когда 1. Теперь, если мы повторим тот же процесс, что и раньше, мы заметим, что он не держится так, как раньше. хорошо.
import statsmodels.api as sm
#A simple input matrix for a deterministic trend component
X = np.array([range(1, len(y) + 1),
np.ones(len(y))]).T
mod = sm.OLS(y, X)
res = mod.fit()
fitted_trend = res.predict(X)
plt.plot(trend, label='actual trend')
plt.plot(fitted_trend, linestyle='dashed', label='fitted trend')
plt.legend()
plt.show()
Ясно, что эта тенденция выключена. Мы аппроксимируем правильный наклон , но у нас есть неправильная точка пересечения. Это весьма тревожно, потому что мы подаем ряды без тренда, чтобы аппроксимировать сезонность, что еще больше усугубляет ошибку. Кроме того, мы не допускаем каких-либо корректировок тренда.
У нас есть один шанс, чтобы все исправить.
В любом случае, давайте перейдем к аппроксимации сезонной составляющей:
#Set the seasonal period seasonal_period = 25 #Get the average of every 25 data points to use as seasonality avg_season = [np.mean(detrended_series[i::seasonal_period], axis=0) for i in range(seasonal_period)] avg_season = np.array(avg_season) fitted_seasonality = np.resize(avg_season, len(y)) plt.plot(fitted_seasonality, label='fitted seasonality') plt.plot(seasonality, label='actual seasonality') plt.legend() plt.show()
Как и ожидалось, на сезонную составляющую не только негативно влияет неправильный тренд, но и, похоже, она поглощает сигнал от экзогенного фактора. Давайте посмотрим, так ли это, построив остатки:
residuals = y - fitted_trend - fitted_seasonality plt.plot(residuals) plt.show()
Напомним, что экзогенный фактор просто прибавлял 20, когда равнялся 1. В идеальном мире эти остатки должны быть случайным шумом около 0 вместе с пиками 20. Без этого сомнительно, что наше экзогенное приближение будет правильным, но давайте попробуем. в любом случае.
Чтобы аппроксимировать экзогенный компонент, мы будем соответствовать остаткам от удаления тренда и десезонализации.
X = np.column_stack([exogenous, np.ones(len(y))]) mod = sm.OLS(residuals, X) res = mod.fit() print(res.summary())
Давайте посмотрим на сводку от OLS:
"""
OLS Regression Results
==============================================================================
Dep. Variable: y R-squared: 0.768
Model: OLS Adj. R-squared: 0.766
Method: Least Squares F-statistic: 325.1
Date: Wed, 10 Nov 2021 Prob (F-statistic): 6.82e-33
Time: 09:46:43 Log-Likelihood: -289.29
No. Observations: 100 AIC: 582.6
Df Residuals: 98 BIC: 587.8
Df Model: 1
Covariance Type: nonrobust
==============================================================================
coef std err t P>|t|
------------------------------------------------------------------------------
x1 15.9072 0.882 18.032 0.000 17.658
const -7.9536 0.624 -12.750 0.000
====================================================================
Ожидаемое влияние нашего экзогенного фактора сократилось с 20 до 16! Исходя из всего, что мы видели до этого момента, это не должно быть сюрпризом. Любая незначительная ошибка в начале процесса приводит к ошибкам компаундирования позже, и не существует механизма для корректировки.
Собираем все вместе:
fitted_exogenous = res.predict(X) fitted = fitted_trend + fitted_seasonality + fitted_exogenous plt.plot(y, label='actuals') plt.plot(fitted, label='fitted') plt.legend() plt.show()
Учитывая все, что пошло не так во время этого процесса, этот результат выглядит хорошо. Помните, что это смоделированный «нетронутый» временной ряд, поэтому тот факт, что мы смогли довольно легко прервать процесс, должен настораживать. Кажется, что процесс может быть слишком жестким. Поскольку сезонная составляющая поглощает некоторый сигнал экзогенной составляющей, следует ли нам поменять порядок аппроксимации? Что, если наш тренд более сложный и нам придется использовать более сложный метод, такой как двойное экспоненциальное сглаживание? Не будет ли это также потенциально съедать сигнал сезонности И экзогенный сигнал?
Мы можем просто решить эту проблему, НЕвыполняя этот процесс вообще. Мы можем перейти к методу, который аппроксимирует все компоненты одновременно, например, к модели Facebook Prophet или, возможно, к какой-либо другой установке, подобной GAM. Но это происходит за счет потери некоторого контроля над отдельными компонентами. Мы больше не можем подгонять разные компоненты разными методами; процесс упростился, но упростились и наши потенциальные прогнозы.
Вопрос в том:
Как мы можем получить преимущества одной процедуры настройки, сохранив при этом полный контроль над отдельными компонентами?
Ответ:
Примените повышение градиента к нашей декомпозиции временных рядов.
Естественный следующий шаг
Простая корректировка предыдущего потока процесса поможет нам эмулировать процедуру одновременной подгонки. После того, как мы аппроксимируем последний шаг, мы будем использовать остатки и перестроить всю процедуру.
Это повышение продолжается до тех пор, пока некоторая глобальная функция стоимости не будет минимизирована. Кроме того, мы можем использовать еще больше концепций, используемых в повышении градиента, таких как скорость обучения. Но теперь мы можем применять скорость обучения к каждому отдельному компоненту, устраняя любые проблемы, связанные с сезонностью, поглощающей экзогенные сигналы. Кроме того, мы можем разрешить компоненту тренда искать точки изменения с помощью бинарной сегментации и, благодаря бустингу, получить очень сложную модель. Сложность, которая теперь может быть достигнута простым усилением, просто безумна (я полагаю, что это техническое слово для этого). Эта недавно достигнутая сложность приводит к большему количеству форм, которые мы можем прогнозировать, что означает более высокий потенциал для большей точности. Мы можем быстро построить ARIMA с 5 точками изменения, сезонностью, аппроксимированной с помощью базисных функций Фурье, и усиленными деревьями для экзогенного прогнозирования.
Но сначала давайте просто спрогнозируем наш простой ряд.
ТимьянBoost
ThymeBoost — это пакет разработчика для Python, который содержит все ранее изложенные функции (и даже больше). Но предупреждаю, прежде чем продолжить, он все еще находится на ранней стадии разработки, поэтому можно столкнуться с ошибками. Пожалуйста используйте на свой страх и риск. Я рекомендую ознакомиться с примерами README на GitHub, а также следовать здесь.
Во-первых, давайте установим из pip, вам обязательно понадобятся самые последние статистические модели и учебные пакеты sci-kit:
pip install ThymeBoostДалее давайте перестроим нашу смоделированную серию здесь (обратите внимание, ThymeBoost требует, чтобы экзогенный был двумерным массивом/DataFrame):
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
sns.set_style('darkgrid')
np.random.seed(100)
trend = np.linspace(1, 50, 100) + 50
seasonality = ((np.cos(np.arange(1, 101))*10))
exogenous = np.random.randint(low=0, high=2, size=len(trend))
y = trend + seasonality + exogenous * 20
#reshape for 2d to pass to ThymeBoost
exogenous = exogenous.reshape(-1, 1)
Импортируйте ThymeBoost и создайте класс:
from ThymeBoost import ThymeBoost as tb boosted_model = tb.ThymeBoost(verbose=1)
Подробный аргумент просто указывает, следует ли печатать журналы процесса повышения. Все, что осталось сделать, это подогнать нашу серию точно так же, как и раньше, но на этот раз с добавлением магии ThymeBoost.
output = boosted_model.fit(y,
trend_estimator='linear',
seasonal_estimator='classic',
exogenous_estimator='ols',
seasonal_period=25,
global_cost='maicc',
fit_type='global',
exogenous=exogenous)
Во-первых, давайте рассмотрим некоторые из этих аргументов:
trend_estimator: как вы хотите аппроксимировать тренд, здесь «линейный» — это простая регрессия, как и раньше.
seasonal_estimator: я уверен, вы уже догадались, метод аппроксимации сезонности, «классический», использует простое среднее значение каждого сезона_периода.
exogenous_estimator: метод, используемый для расчета экзогенного компонента.
seasonal_period: ожидаемый сезонный период.
global_cost: Наконец, интересный аргумент! Это управляет процедурой усиления, maicc — это модифицированный «критерий информации Akaike», где модификации используются для включения раундов повышения в формулировку.
fit_type: здесь мы обозначаем, хотим ли мы, чтобы ThymeBoost искал точки изменения или нет, «глобальный» будет подгонять наш trend_estimator ко всем данным. Если мы пропустим «local», мы выполним бинарную сегментацию, чтобы найти потенциальную точку изменения.
Результаты, сохраненные для «вывода», представляют собой Pandas DataFrame с подогнанными значениями и каждым отдельным компонентом. Если у вас есть опыт работы с Prophet, то это должно напомнить вам об этом.
Давайте нарисуем подобранный тренд, сезонность и экзогенный компонент и посмотрим, как это было.
Приблизительный тренд:
plt.plot(trend, label='actual trend') plt.plot(output['trend'], label='fitted trend') plt.legend() plt.show()
Приблизительная сезонность:
plt.plot(seasonality, label='actual seasonality') plt.plot(output['seasonality'], label='fitted seasonality') plt.legend() plt.show()
Приблизительный экзогенный
plt.plot(exogenous * 20, label='actual exogenous') plt.plot(output['exogenous'], label='fitted exogenous') plt.legend() plt.show()
Это значительное улучшение по сравнению с традиционным процессом. Мы смогли точно аппроксимировать каждый отдельный компонент. ThymeBoost также предоставляет некоторые вспомогательные функции для более быстрого построения графиков:
boosted_model.plot_components(output)
В пакете используются стандартные методы прогнозирования соответствия. Давайте используем ThymeBoost для создания прогнозов и их построения с помощью метода plot_results:
#create a future exogenous input
forecast_horizon = 20
np.random.seed(100)
future_exogenous = np.random.randint(low=0, high=2, size=forecast_horizon)
#use predict method and pass fitted output, forecast horizon, and future exogenous
predicted_output = boosted_model.predict(output,
forecast_horizon=forecast_horizon,
future_exogenous=future_exogenous)
boosted_model.plot_results(output, predicted_output)
Заключение
ThymeBoost устраняет многие основные проблемы с традиционной декомпозицией временных рядов и позволяет нам сохранять полный контроль над каждым компонентом. Процедура даже дает нам доступ к новым интересным функциям благодаря повышающим итерациям. Например, распространенной проблемой временных рядов является множественная сезонность. Thymeboost может решить эту проблему вполне естественным образом, изменяя сезонный период во время каждого раунда повышения с помощью простого изменения нашего вызова fit, который теперь передает список в Season_period.
От:
output = boosted_model.fit(y,
trend_estimator='linear',
seasonal_estimator='classic',
exogenous_estimator='ols',
seasonal_period=25,
global_cost='maicc',
fit_type='global',
exogenous=exogenous)
To:
output = boosted_model.fit(y,
trend_estimator='linear',
seasonal_estimator='classic',
exogenous_estimator='ols',
seasonal_period=[25, 5],
global_cost='maicc',
fit_type='global',
exogenous=exogenous)
ThymeBoost будет переключаться между двумя сезонными периодами, указанными в списке, каждый раз внося коррективы в общую сезонную составляющую. Но дело не в сезонности. Мы также можем сделать это с каждым аргументом, переданным для соответствия. Это означает, что мы можем начать с простого линейного тренда, за которым следует ARIMA, а затем двойное экспоненциальное сглаживание в последнем раунде повышения. Вы имеете полный контроль над всей процедурой.
Как упоминалось ранее, этот пакет все еще находится на ранней стадии разработки. Огромное количество возможных конфигураций делает отладку сложной процедурой, поэтому используйте ее на свой страх и риск. Но, пожалуйста, поэкспериментируйте и открывайте любые проблемы, с которыми вы столкнетесь на GitHub!
