Зачем маркетологу прогнозирование на sktime? Часть 1: Теоретическая

Предположим, вам задают простой вопрос: "Сколько трафика на сайте будет вечером?" Можно, конечно прикинуть на глазок или ответить, примерно столько же, сколько и вчера (такой ответ называется наивным прогнозом).

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

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

Но, для начала немного разберемся с терминологией. В разных библиотеках машинного обучения, терминология может немного различаться. Например, в sklearn (scikit-learn) в идеологии которого сделан sktime популярно определение оценщика (оценивателя), который отвечает за настройку модели на данных, включая выбор гирерпараметров и их обучение. Модель — это результат обучения оценщика.

Так же нам могут понадобиться:

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

Таким образом ответ на вопрос: "Сколько трафика на сайте будет вечером?" с точки зрения машинного обучения сводится к задаче прогнозирования (предсказания) значений временного ряда (time series forecasting). И состоит из следующий этапов:

Сбор\Загрузка данных
Предобработка данных
Разделение на обучающую и тестовую выборку
Выбор модели и подбор гиперпараметров
Оценка модели (расчет ошибок, кросс-валидация и т.д.)
Выполнение прогноза и визуализация
Размещение модели в продакшен (это маркетологу не нужно)

Для первого урока мы свои данные собирать не будем, а воспользуемся уже готовым датасетом по авиаперевозкам:

from sktime.datasets import load_airline y = load_airline() # (целевая переменная, то что будем предсказывать)

Взглянем на данные:

y.plot();

В датасете 144 наблюдения с января 1949 года по декабрь 1960 (за 12 лет)

Одномерные временные ряды в sktime должны быть представлены в Pandas Series, что очень удобно. Dataframe используется для многомерных рядов или для признаков. Если вы плохо понимаете о чем идет речь, то ничего страшного, со временем разберетесь, просто выполняйте полный код в конце статьи и смотрите, как работает.

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

from sktime.split import temporal_train_test_split y_train, y_test = temporal_train_test_split(y=y, test_size=12)

В данных данных предсказание обычно делают на 12 месяцев, поэтому я указал test_size=12

Теперь нужно указать fh (future horizon) т.е. на какой период будем делать прогноз. В sktime это делается немного сложно:

from sktime.forecasting.base import ForecastingHorizon fh = ForecastingHorizon(y_test.index, is_relative=False)

Или можно по-простому, по-деревенски ) Тоже будет работать.

fh = range(1, 12+1)

И остался самый сложный момент, это указать сезонность. В данном датасете она тоже равна 12, потому что, прогноз обычно делают, на один сезон. И сезонность в месячных данных обычно равна 12. Но предположим, у вас какие-то другие данные и\или вы не знаете какая сезонность:

from sktime.param_est.seasonality import SeasonalityACF # Автоматическое определение сезонности sp_est = SeasonalityACF() sp_est.fit(y.diff().dropna()) sp = sp_est.get_fitted_params().get("sp") print(sp)

Результат выполнения кода: 12

Автоматическое определение сезонности! Просто? Элегантно? Ох, я месяцами мучался с этим определением сезонности, жаль я этого способа не знал раньше. Важно: перед определением сезонности этим способом всегда надо сделать вычисление разности .diff() и избавится от первого NaN в начале.

Все, самая сложная часть кода пройдена, теперь осталось только откинуться в кресле и позволить библиотеке сделать всю работу за нас:

from sktime.forecasting.ets import AutoETS from sktime.utils.plotting import plot_series from sktime.performance_metrics.forecasting import mean_squared_error # Создаем объект оценщика forecaster = AutoETS(auto=True, sp=sp) # Обучаем модель forecaster.fit(y_train, fh=fh) # Делаем прогноз y_pred = forecaster.predict() # Считаем ошибку rmse = mean_squared_error(y_test, y_pred, square_root=True) print(f"RMSE: {rmse}") plot_series(y_train, y_test, y_pred, labels=["Тренировочные данные", "Тестовые данные", "Прогноз"]);

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

MAE (Mean Absolute Error) — Средняя абсолютная ошибка. Показывает среднее абсолютное отклонение предсказаний от фактических значений.

RMSE (Root Mean Squared Error) — Корень из среднеквадратичной ошибки. Оценка среднего отклонения без квадратных единиц отклонения.

MASE (Mean Absolute Scaled Error) — Средняя абсолютная масштабированная ошибка. Применяется для сравнения моделей на разных временных рядах, масштабируя ошибку.

RMSSE (Root Mean Squared Scaled Error) — Корень из среднеквадратичной масштабированной ошибки. Тоже самое, но из MASE.

MAPE (Mean Absolute Percentage Error) — Средняя абсолютная процентная ошибка. Отражает среднюю относительную ошибку предсказаний. Удобно, что измеряется в процентах. Из минусов, что на маленьких значениях получается очень большой.

SMAPE (Symmetric Mean Absolute Percentage Error) — Симметричная MAPE.

R² (Coefficient of Determination) — Коэффициент детерминации, специально придуман для сравнения с наивным прогнозом, где 1 — идеальная модель, 0 — случайные предсказания. Если значение отрицательное, значит модель проигрывает наивному прогнозу.

Я использовал довольно мощную AutoETS модель портированную из пакета R. Расшифровывается как E = Error (ошибка), T = Trend (Тренд), S = Seasonal (Сезонность), а параметр auto=True (по умолчанию False) отвечает за автоматический подбор гиперпараметров. Данная модель является продолжением модели Хольта-Винтерса и входит в класс экспоненциального сглаживания. В библиотеке statasmodels такой модели нет, но есть аналогичные.

Результат получился не супер пупер, но зато все полностью 100% автоматизировано, включая определение сезонности. Главное, модель работает очень быстро.

Осталось сделать прогноз на следующие 3 года вперед.

y_pred = forecaster.fit_predict(y, fh=range(1, 36+1)) plot_series(y, y_pred, labels=["Исходные данные", "Прогноз"]);

Т.к. за 61-63 года данных нет, ошибку посчитать мы не можем.

Сразу можем использовать метод .fit_predict() для обучения модели на полном датасете y и выполнении прогноза.

Осталась небольшая факультативная часть, которая маркетологу вообще не нужна (можете смело пропускать), это саммари модели:

forecaster.summary()

ETS Results Dep. Variable: Number of airline passengers No. Observations: 132 Model: ETS(MAM) Log Likelihood -465.990 Date: Thu, 26 Sep 2024 AIC 967.981 Time: 12:51:50 BIC 1019.871 Sample: 01-31-1949 HQIC 989.066 - 12-31-1959 Scale 0.001 Covariance Type: approx coef std err z P>|z| [0.025 0.975] smoothing_level 0.7731 0.087 8.907 0.000 0.603 0.943 smoothing_trend 7.731e-05 nan nan nan nan nan smoothing_seasonal 2.269e-05 nan nan nan nan nan initial_level 111.0577 1279.351 0.087 0.931 -2396.424 2618.539 initial_trend 1.9264 22.213 0.087 0.931 -41.609 45.462 initial_seasonal.0 0.9865 11.365 0.087 0.931 -21.288 23.261 initial_seasonal.1 0.8783 10.118 0.087 0.931 -18.953 20.710 initial_seasonal.2 1.0111 11.648 0.087 0.931 -21.819 23.841 initial_seasonal.3 1.1652 13.423 0.087 0.931 -25.144 27.475 initial_seasonal.4 1.3404 15.442 0.087 0.931 -28.925 31.606 initial_seasonal.5 1.3502 15.554 0.087 0.931 -29.136 31.836 initial_seasonal.6 1.2214 14.071 0.087 0.931 -26.357 28.800 initial_seasonal.7 1.0787 12.426 0.087 0.931 -23.277 25.434 initial_seasonal.8 1.0812 12.456 0.087 0.931 -23.331 25.494 initial_seasonal.9 1.1280 12.995 0.087 0.931 -24.342 26.598 initial_seasonal.10 0.9846 11.343 0.087 0.931 -21.247 23.216 initial_seasonal.11 1.0000 11.520 0.087 0.931 -21.579 23.579 Ljung-Box (Q): 34.86 Jarque-Bera (JB): 3.03 Prob(Q): 0.07 Prob(JB): 0.22 Heteroskedasticity (H): 0.51 Skew: -0.37 Prob(H) (two-sided): 0.03 Kurtosis: 3.01 Warnings: [1] Covariance matrix calculated using numerical (complex-step) differentiation.

Тут нас интересует маленький кусочек: Model: ETS(MAM), где:

M — мультипликативная ошибка,
A — аддитивный тренд,
M — мультипликативная сезонность.

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

Таким образом финальный код полностью:

from sktime.datasets import load_airline from sktime.forecasting.ets import AutoETS from sktime.forecasting.base import ForecastingHorizon from sktime.param_est.seasonality import SeasonalityACF from sktime.utils.plotting import plot_series from sktime.split import temporal_train_test_split from sktime.performance_metrics.forecasting import mean_squared_error y = load_airline() y_train, y_test = temporal_train_test_split(y=y, test_size=12) # fh = range(1, 12+1) # Одно и тоже fh = ForecastingHorizon(y_test.index, is_relative=False) # Оцениваем сезонность sp_est = SeasonalityACF() sp_est.fit(y.diff().dropna()) sp = sp_est.get_fitted_params().get("sp") forecaster = AutoETS(auto=True, sp=sp) forecaster.fit(y_train, fh=fh) y_pred = forecaster.predict() rmse = mean_squared_error(y_test, y_pred, square_root=True) print(f"RMSE: {rmse}") plot_series(y_train, y_test, y_pred, labels=["Тренировочные данные", "Тестовые данные", "Прогноз"]); y_pred = forecaster.fit_predict(y, fh=range(1, 36+1)) plot_series(y, y_pred, labels=["Исходные данные", "Прогноз"]);

Можете просто скопировать и выполнить.

Выводы:

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

Так же, в следующих уроках я хочу рассмотреть:

Прогнозирование на реальных данных
Разделение временного ряда на компоненты
Автоматический выбор модели для прогноза
Подбор гиперпараметров
Кросс-валидацию
Бектестинг

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

UPD: Вторая часть