Эта статья является переводом оригинальной статьи «Forecasting at Scale», написанной Шоном Дж. Тейлором и Бенджамином Летэмом из компании Facebook. Перевод выполнен для русскоязычной аудитории с сохранением оригинального смысла и структуры. Название изменено на «Математический аппарат Facebook Prophet», поскольку в нем подробно изложены математические основы расчета прогнозов в рамках данного фреймворка.

Оригинальная статья: https://peerj.com/preprints/3190.pdf

1. Введение

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

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

Первые два аспекта масштабируемости, которые мы рассматриваем, заключаются в том, что методы бизнес-прогнозирования должны: 

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

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

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

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

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

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

2. Особенности временных рядов в бизнес-задачах

Существует большое разнообразие задач бизнес-прогнозирования, однако многие из них обладают общими характеристиками. На рисунке 2 представлен типичный временной ряд данных Facebook, относящийся к событиям (Facebook Events). Пользователи Facebook могут использовать платформу Events для создания страниц мероприятий, приглашения других пользователей и взаимодействия с событиями различными способами.

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

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

Этот временной ряд наглядно демонстрирует трудности, связанные с получением обоснованных прогнозов с использованием полностью автоматизированных методов. На рис. 3 представлены прогнозы, построенные с помощью нескольких автоматических методов из пакета forecast в R, описанных в работе Hyndman et al. (2007). Прогнозы были получены по трем отрезкам времени, для обучении модели использовалась только часть временного ряда до первой точки соответствующего отрезка.

Методы, представленные на рисунке, включают:

• auto.arima – автоматический выбор наилучшей модели из семейства ARIMA;
• ets – подбор и выбор лучшей модели экспоненциального сглаживания (Hyndman et al. 2002);
• snaive – сезонная наивная модель (seasonal naive), которая использует случайное блуждание и делает постоянные прогнозы с учетом недельной сезонности;
• tbats – модель TBATS, учитывающая как недельную, так и годовую сезонности (De Livera et al. 2011).

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

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

3. Модель прогнозирования Prophet

Теперь мы опишем модель прогнозирования временных рядов, разработанную с учетом общих характеристик бизнес-временных рядов, показанных на рис. 2. Важно отметить, что она также имеет интуитивно понятные параметры, которые можно настраивать, не зная деталей базовой модели. Это необходимо аналитику для эффективной настройки модели, как показано на рис. 1. Наша реализация доступна в виде программного обеспечения с открытым исходным кодом на Python и R под названием Prophet — (https://facebook.github.io/prophet/)

Мы используем декомпозируемую модель временных рядов (Harvey & Peters 1990) с тремя основными компонентами: трендом, сезонностью и праздниками. Эти компоненты объединяются в следующем уравнении:$$y(t) = g(t) + s(t) + h(t) + \varepsilon_t$$

Здесь [math]g(t)[/math] — функция тренда, которая моделирует непериодические изменения значений временного ряда, [math]s(t)[/math] представляет периодические изменения (например, еженедельную и ежегодную сезонность), а [math]h(t)[/math] представляет влияние праздников, которые происходят по потенциально нерегулярному графику в течение одного или нескольких дней. Член ошибки [math]\varepsilon_t[/math] обозначает любые индивидуальные изменения, не учтённые моделью; позже мы выдвинем параметрическое предположение о том, что [math]\varepsilon_t[/math] имеет нормальное распределение.

Эта спецификация аналогична обобщенной аддитивной модели (GAM) (Hastie & Tibshirani, 1987), классу регрессионных моделей с потенциально нелинейными сглаживателями, применяемыми к регрессорам. Здесь мы используем только время в качестве регрессора, но, возможно, несколько линейных и нелинейных функций времени в качестве компонентов. Моделирование сезонности как аддитивного компонента — это тот же подход, что и при экспоненциальном сглаживании (Гарднер, 1985). Мультипликативная сезонность, где сезонный эффект является фактором, умножающим [math]g(t)[/math], может быть достигнута с помощью логарифмического преобразования.

Выбор в пользу GAM имеет то преимущество, что она легко декомпозируется и при необходимости позволяет добавлять новые компоненты, например, когда обнаруживается новый источник сезонности. Кроме того, GAM обучаются очень быстро — либо с помощью бэкфиттинга (обратной подгонки), либо с использованием метода L-BFGS (Byrd et al. 1995) (мы предпочитаем последний), что позволяет пользователю интерактивно изменять параметры модели.

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

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

1. Гибкость: Мы можем легко учитывать сезонность с несколькими периодами и позволить аналитику делать различные предположения о трендах.
2. В отличие от моделей ARIMA, измерения не обязательно должны быть стационарными (равномерно распределены), и нам не нужно интерполировать отсутствующие значения, например, при удалении выбросов.
3. Обучение модели происходит очень быстро, что позволяет аналитику интерактивно исследовать множество вариантов спецификаций модели, например, в приложении Shiny (Chang et al. 2015).
4. Прогнозная модель обладает легко интерпретируемыми параметрами, которые аналитик может изменять, чтобы накладывать определённые предположения на прогноз. Более того, аналитики обычно имеют опыт работы с регрессией и без труда могут расширить модель, добавив новые компоненты.

Автоматическое прогнозирование имеет долгую историю и включает в себя множество методов, адаптированных к конкретным типам временных рядов (Tashman & Leach, 1991, De Gooijer & Hyndman, 2006). Наш подход обусловлен как характером временных рядов, которые мы прогнозируем в Facebook (кусочные тренды, множественная сезонность, плавающие праздники), так и сложностями, связанными с масштабируемым прогнозированием.

3.1 Модель тренда

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

3.1.1 Нелинейный, насыщенный рост

Для прогнозирования роста основным компонентом процесса генерации данных является модель, предсказывающая, как выросла численность населения и как она будет расти в будущем. Моделирование роста в Facebook в целом похоже на рост населения в природных экосистемах (например, Хатчинсон, 1978), где наблюдается нелинейный рост, который насыщается при максимальной пропускной способности. Например, количество пользователей Facebook в определенной области может быть рассчитано как количество людей, имеющих доступ к сети Интернет. Такого рода рост обычно моделируется с использованием логистической модели роста, которая в своей самой базовой форме представляет собой $$g(t) = \frac{C}{1 + \exp(-k(t — m))}, \tag{2}$$

где [math]C[/math] — верхний порог (пропускная способность), [math]k[/math] — темп роста, а [math]m[/math] — параметр смещения, позволяющий «сдвигать» функцию вдоль оси времени.

Есть два важных аспекта роста Facebook, которые не отражены в уравнении (2).
Во–первых, пропускная способность непостоянна — с увеличением числа людей в мире, имеющих доступ к Интернету, увеличивается и предел роста. Таким образом, мы заменяем фиксированную пропускную способность [math]C[/math] на изменяющуюся во времени пропускную способность [math]C(t[/math] ). Во-вторых, темпы роста не являются постоянными. Выход новых продуктов может существенно изменить темпы роста в регионе, поэтому модель должна учитывать изменяющиеся темпы, чтобы соответствовать историческим данным.

Мы учитываем изменения тренда в модели роста, четко определяя точки изменений, в которых допускается изменение темпов роста. Мы включаем изменения тренда в модель роста, явно определяя точки изменения (change points), в которых темп роста может измениться. Предположим, что существует [math]S[/math] точек изменения в моменты времени [math]s_j[/math], [math]j = 1,…,S[/math].

Мы определяем вектор корректировки скорости роста [math]\boldsymbol{\delta} \in \mathbb{R}^S[/math], где [math]\delta_j[/math] — изменение скорости которое происходит в момент времени [math]s_j[/math]. Тогда скорость в любой момент [math]t[/math] равна базовой скорости [math]k[/math] + все корректировки до этой точки: [math]k + \sum_{j : t > s_j} \boldsymbol{\delta}_j[/math]. Это можно выразить более лаконично, определив вектор [math]\mathbf{a}(t) \in \{0, 1\}^S[/math], таким образом, что $$a_j(t) = \begin{cases} 1, & \text{if } t \geq s_j; \\ 0, & \text{в остальных случаях}. \end{cases}$$ Выходит, что скорость в момент времени [math]t[/math] равна [math]k+\boldsymbol{a}(t)^T \boldsymbol{\delta}[/math]. При изменении скорости [math]k[/math] необходимо также настроить параметр смещения [math]m[/math], чтобы соединить концы сегментов. Правильная настройка в точке изменения [math]j[/math] легко вычисляется как $$\gamma_j = \left(s_j — m — \sum_{l<j} \gamma_l\right) \left(1 — \frac{k + \sum_{l<j} \delta_l}{k + \sum_{l \leq j} \delta_l}\right)$$. Тогда кусочно-логистическая модель роста принимает следующий вид:$$g(t) = \frac{C(t)}{1 + \exp\left(-(k + \mathbf{a}(t)^\top \boldsymbol{\delta})(t — (m + \mathbf{a}(t)^\top \boldsymbol{\gamma}))\right)}, \tag{3}$$Важным параметром в нашей модели является [math]C(t)[/math], или ожидаемая пропускная способность системы в любой момент времени. Аналитики часто имеют представление о размерах рынка и могут настроить их соответствующим образом. Могут также существовать внешние источники данных, которые могут обеспечить оперативность, такие как демографические прогнозы Всемирного банка.

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

3.1.2 Линейный тренд с точками изменения

Для задач прогнозирования, в которых нет роста с насыщением, кусочно-постоянный темп роста обеспечивает экономичную и часто полезную модель. Здесь модель тренда представляет собой $$g(t) = \left(k + \mathbf{a}(t)^\top \boldsymbol{\delta}\right)t + \left(m + \mathbf{a}(t)^\top \boldsymbol{\gamma}\right), \tag{4}$$ где, как и прежде, [math]k[/math] — скорость роста, [math]\boldsymbol{\delta}[/math] — корректировка скорости, [math]m[/math] — параметр смещения, а [math]\gamma_j[/math] устанавливается равной — [math]s_j \delta_j[/math], чтобы функция была непрерывной.

3.1.3 Автоматический отбор точек изменения

Точки изменения [math]s_j​[/math] могут быть заданы аналитиком, используя известные даты запусков продуктов или другие события, изменяющие рост, либо могут быть выбраны автоматически из заданного набора значений. Автоматический выбор может быть выполнен совершенно естественным образом с помощью формул, приведенных в (3) и (4), выбрав априорное распределение Лапласа для [math]\boldsymbol{\delta}[/math]

Мы указываем большое количество точек изменения (например, по одной в месяц на протяжении нескольких лет) и используем [math]\delta_j \sim \text{Laplace}(0, \tau)[/math]. Параметр [math]\tau[/math] напрямую управляет гибкостью модели в изменении её скорости. Важно отметить, что применение распределение Лапласа для корректировок [math]\boldsymbol{\delta}[/math] не влияет на основную скорость роста [math]\tau[/math]. Таким образом, при стремлении [math]\tau[/math] к нулю модель сводится к стандартному (не кусочному) логистическому или линейному росту.

3.1.4 Неопределенность прогноза тренда

Когда происходит экстраполяции модели за пределы исторических данных для получения прогноза, тренд будет иметь постоянную скорость.Мы оцениваем неопределенность прогнозируемого тренда, расширяя генеративную модель на будущее. Генеративная модель тренда заключается в том, что существует [math]S[/math] точек изменения на протяжении истории из [math]T[/math] точек, каждая из которых имеет скорость изменения [math]\delta_j \sim \text{Laplace}(0, \tau)[/math]. Мы моделируем будущие изменения скорости, которые повторяют прошлые, заменяя [math]\tau[/math] на дисперсию, полученную на основе данных.В полностью байесовской структуре это можно было бы сделать с помощью иерархического априора на [math]\tau[/math], чтобы получить его апостериор, в противном случае мы можем использовать оценку максимального правдоподобия параметра шкалы скорости: [math]\lambda = \frac{1}{S} \sum_{j=1}^S |\delta_j|[/math]. Точки будущих изменений выбираются случайным образом так, чтобы средняя частота точек изменения соответствовала средней частоте точек изменения в исторических данных: $$\forall j > T, \quad \begin{cases} \delta_j = 0 \quad \text{с вероятностью } \frac{T — S}{T}, \\ \delta_j \sim \text{Laplace}(0, \lambda) \quad \text{с вероятностью } \frac{S}{T}. \end{cases}$$ Таким образом, мы измеряем неопределенность прогнозируемого тренда, предполагая, что в будущем будет наблюдаться та же средняя частота и величина изменений скорости, которые наблюдались в прошлом. После того, как [math]\lambda[/math] был выведен из данных, мы используем эту генеративную модель для моделирования возможных будущих трендов и используем смоделированные тренды для вычисления доверительных интервалов.

Предположение о том, что тренд будет продолжать меняться с той же частотой и скоростью изменений, что и на протяжении всех исторических данных , является достаточно сильным, поэтому мы не ожидаем высокой точности от доверительных интервалов. Тем не менее, они являются полезным показателем уровня неопределенности и, в частности, показателем переобучения. С увеличением [math]\tau[/math] модель становится более гибкой, обучаясь на исторических данных, поэтому ошибки обучения уменьшаются. Однако при прогнозировании будущего, эта гибкость приведет к широким доверительным интервалам.

3.2 Сезонность

Временные ряды в бизнес-задачах часто имеют многопериодную сезонность в результате человеческого поведения, которое они отражают. Например 5-ти дневная рабочая неделя может оказывать влияние на временной ряд, повторяющееся каждую неделю, в то время как график отпусков и школьных каникул может оказывать влияние, повторяющееся каждый год. Чтобы моделировать и прогнозировать данные эффекты, мы должны задать модели сезонности, которые являются периодическими функциями времени [math]t[/math].

Мы используем ряды Фурье для создания гибкой модели периодических эффектов (Harvey & Shephard, 1993). Пусть [math]P[/math] — это постоянное значение периода который мы ожидаем у временного ряда (например, [math]P=365,25[/math] для годовых данных или [math]P=7[/math] для недельных данных, когда переменная времени масштабируется в днях). Мы можем аппроксимировать произвольные гладкие сезонные эффекты с помощью стандартного ряда Фурье²: $$s(t) = \sum_{n=1}^N \left( a_n \cos\left(\frac{2\pi n t}{P}\right) + b_n \sin\left(\frac{2\pi n t}{P}\right) \right).$$ Подгонка сезонности требует оценки [math]2N[/math] параметров [math]\boldsymbol{\beta} = [a_1, b_1, \dots, a_N, b_N]^\top[/math]. Это достигается путем построения матрицы векторов сезонности для каждого значения [math]t[/math] в наших исторических и будущих данных. Например, для годовой сезонности и [math]N = 10[/math]: $$X(t) = \left[ \cos\left(\frac{2\pi (1)t}{365.25}\right), \dots, \sin\left(\frac{2\pi (10)t}{365.25}\right) \right]. \tag{5}$$ Тогда сезонная компонента будет иметь вид $$s(t) = X(t) \boldsymbol{\beta}. \tag{6}$$В нашей генеративной модели мы принимаем [math]\boldsymbol{\beta} \sim \text{Normal}(0, \sigma^2)[/math], чтобы сгладить сезонность.

При сокращении ряда до [math]N[/math] применяется фильтр нижних частот для учета сезонности, поэтому увеличение [math]N[/math] позволяет подобрать сезонные модели, которые меняются быстрее, хотя и с повышенным риском переобучения. Для годовой и недельной сезонности мы нашли значения [math]N = 10[/math] и [math]N = 3[/math] соответственно, которые хорошо подходят для большинства задач. Выбор этих параметров можно автоматизировать с помощью процедуры выбора модели, такой как AIC.

3.3 Праздники и события

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

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

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

Включение этого списка праздников в модель упрощается если предположить, что эффекты от праздников независимы друг от друга. Для каждого праздника [math]i[/math] вводится [math]D_i[/math] набор прошлых и будущих дат этого праздника. Мы добавляем функцию-индикатор, проверяющую, является ли момент времени [math]t[/math] событием [math]i[/math], и назначаем каждому событию [math]i[/math] параметр [math]k_i[/math], который является соответствующим изменением прогноза. Это достигается тем же путем, что и в части с сезонной составляющей, путем создания матрицы регрессоров $$Z(t) = \left[ 1(t \in D_1), \dots, 1(t \in D_L) \right]$$ и ввода ее в модель следующим образом: $$h(t) = Z(t) \boldsymbol{\kappa}. \tag{7}$$ Как и в случае с сезонностью, мы используем априорное распределение [math]\boldsymbol{\kappa} \sim \text{Normal}(0, \nu^2)[/math]. Часто важно учитывать эффекты создаваемые днями, хронологически расположенными рядом с определенными праздниками, например, выходных в честь Дня благодарения. Для этого мы добавляем дополнительные параметры для дней, окружающих праздник, по сути, рассматривая каждый из дней в этом окне как сам праздник.

3.4 Обучение модели

Когда сезонность и индикаторы праздников для каждого наблюдения объединены в матрицу [math]X[/math], а индикаторы точек изменений [math]a(t)[/math] – в матрицу [math]A[/math] — модель может быть выражена в нескольких строках кода Stan. Для оценивания параметров обучаемой модели мы используем L-BFGS в Stan, но также можем выполнить полное апостериорный вывод (full posterior inference), чтобы учесть неопределенность параметров модели в неопределенности прогноза.

model {
// Priors
k ∼ normal(0, 5);
m ∼ normal(0, 5);
epsilon ∼ normal(0, 0.5);
delta ∼ double_exponential(0, tau);
beta ∼ normal(0, sigma);

// Logistic likelihood
y ∼ normal(C ./ (1 + exp(-(k + A * delta) .* (t - (m + A * gamma)))) +
          X * beta, epsilon);
// Linear likelihood
y ∼ normal((k + A * delta) .* t + (m + A * gamma) + X * beta, sigma);
}

// Модель временного ряда с логистическим и линейным ростом

// Параметры модели:
k      // Скорость изменения тренда (наклон)
m      // Начальное значение тренда
epsilon // Стандартное отклонение шума (логистическая модель)
sigma  // Стандартное отклонение шума (линейная модель)
delta  // Влияние точек изменений (changepoints) на тренд (~Laplace)
beta   // Коэффициенты сезонных и праздничных эффектов (~Normal)
tau    // Гиперпараметр распределения Laplace для delta
gamma  // Влияние точек изменений на положение тренда

// Логистическая модель:
y ~ normal(C / (1 + exp(-(k + A * delta) * (t - (m + A * gamma)))) + X * beta, epsilon);
// Используется для моделирования ограниченного роста (например, насыщение рынка).

// Линейная модель:
y ~ normal((k + A * delta) * t + (m + A * gamma) + X * beta, sigma);
// Используется для моделирования неограниченного линейного тренда.

Рис. 4 показывает прогноз модели Prophet для временного ряда событий Facebook, представленного на рис. 3. Эти прогнозы были сделаны на тех же трех датах, что и на рис. 3, используя только данные, доступные на момент предшествующей первой прогнозируемой даты. Прогноз Prophet способен учитывать как недельную, так и годовую сезонность. В отличие от базовых моделей, представленных на рис. 3, он не чрезмерно реагирует на падение активности в праздничный период в первый год. В первом прогнозе модель Prophet слегка переобучилась на годовой сезонности, так как использовала только один год данных. В третьем прогнозе модель еще не «поняла», что тренд изменился.

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

3.5 Моделирование в соответствии с принципом «‎analyst-in-the-loop»‎

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

Размеры рынка: Аналитики могут располагать внешними данными об общем объеме рынка и могут применять эти знания непосредственно, определяя мощность.

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

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

Параметры сглаживания: Регулируя [math]\tau[/math], аналитик может выбрать модель с более глобальным или более локально сглаженным подходом. Параметры сглаживания сезонности и праздников [math](\sigma, \nu)[/math]позволяют аналитику указать модели, какая доля исторической сезонной вариации ожидается в будущем.

С помощью хороших средств визуализации аналитики могут использовать эти параметры для улучшения обучения модели. Когда обучение модели строится на исторических данных, быстро становится очевидно, если автоматический выбор точек изменений пропустил какие-либо из них. Параметр [math]\tau[/math] представляет собой единый регулятор, который можно настроить для увеличения или уменьшения гибкости тренда, а [math]\sigma[/math] — регулятор, позволяющий усилить или ослабить влияние компонента сезонности. Визуализация открывает множество других возможностей для плодотворной корректировки прогноза: это может быть выбор между линейным трендом и логистическим ростом, определение временных масштабов сезонности или выявление аномальных периодов, которые следует исключить из модели. Данные корректировки можно осуществлять без глубоких статистических знаний, что позволяет аналитикам применять свои инсайты и знания в предметной области.

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

Наш подход к моделированию в соответствии принципом «analyst-in-the-loop» — является альтернативным, поскольку он стремится объединить преимущества статистических так и экспертных прогнозов, концентрируя усилия аналитика на улучшении модели в случае необходимости, а не на непосредственном формировании прогнозов с использованием какой-либо неформализованной процедуры. Мы обнаружили, что наш подход очень похож на цикл «преобразование-визуализация-моделирование», предложенный Wickham & Grolemund (2016), в котором знания эксперта из предметной области воплощаются в улучшенной модели после некоторого количества итераций.

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

Реализация подхода «analyst-in-the-loop» на крупномасштабном уровне зависит, в значительной степени зависит от автоматической оценки качества прогнозов и хороших инструментов визуализации. Теперь мы расскажем, как можно автоматизировать оценку прогнозов, чтобы определить наиболее релевантные прогнозы для аналитиков.

4. Автоматизация оценки качества прогнозов

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

4.1 Использование прогнозов базовых моделей

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

4.2 Оценка качества прогноза

Прогнозы составляются на определенный период, который мы обозначаем как H. Горизонт – это количество дней в будущем, на которое мы хотим сделать прогноз — обычно это 30, 90, 180 или 365 дней в наших задачах. Таким образом, для любого прогноза, основанного на ежедневных наблюдениях, мы получаем H прогнозов, каждая из которых будет связана с некоторой ошибкой. Нам нужно задать прогнозируемую целевую функцию, чтобы сравнивать методы и отслеживать их эффективность. Кроме того, понимание того, насколько наша процедура прогнозирования подвержена ошибкам, может позволить потребителям прогнозов в бизнес-среде определить, стоит ли ей вообще доверять.

Пусть [math] \hat{y}(t | T)[/math] представляет прогноз для момента времени [math]t[/math] , полученный с учетом исторической информации до момента времени [math]T[/math], а [math]d(y, y’)[/math] — это метрика расстояния, такая как средняя абсолютная ошибка, [math]d(y, y’) = |y-y’|[/math]. Выбор функции расстояния должен зависеть от конкретной задачи. De Gooijer & Hyndman (2006) сделали обзор нескольких таких метрик — на практике мы предпочитаем среднюю абсолютную процентную ошибку (MAPE) из-за её интерпретируемости. Мы определяем эмпирическую точность прогноза на [math]h \in (0, H)[/math] периодов вперед относительно времени [math]T[/math] следующим образом: $$\phi(T, h) = d(\hat{y}(T + h | T), y(T + h)).$$

Чтобы оценить качество и определить, как она изменяется с [math]h[/math], часто задают параметрическую модель для члена ошибки и оценивают её параметры на основе данных. Например, если мы используем модель [math]AR(1) \), \( y(t) = \alpha + \beta y(t — 1) + \nu(t)[/math], мы предполагаем, что [math]\nu(t) \sim \text{Normal}(0, \sigma_{\nu}^{2})[/math] и сосредотачиваемся на оценке дисперсии [math]\sigma_{\nu}^{2}[/math] по данным. Тогда можно вычислить математические ожидания, используя любую функцию расстояния через моделирование или аналитическое выражение для математического ожидания суммы ошибок. Однако такие подходы дают правильные оценки ошибок только при условии, что мы правильно задали модель процесса, — а это условие в реальности редко выполняется.

Мы предпочитаем непараметрический подход к оценке ожидаемых ошибок, который применим к разным моделям. Этот подход аналогичен применению кросс-валидации для оценки ошибки на тестовых данных для моделей, предсказывающих на [math]i.i.d.[/math] данных. Имея набор исторических прогнозов, мы строим модель ожидаемой ошибки при различных горизонтах прогноза [math]h[/math].$$\begin{equation}\xi(h) = \mathbb{E}[\phi(T, h)].\tag{8}\end{equation}$$Эта модель должна быть гибкой, но может выдвигать некоторые простые предположения. Во-первых, функция должна быть локально гладкой по [math]h[/math], поскольку мы ожидаем, что ошибки, которые мы допускаем на последовательные дни, будут относительно схожими. Во-вторых, мы можем наложить предположение о том, что функция должна быть слабо возрастающей по [math]h[/math], хотя это не обязательно для всех моделей прогнозирования.

На практике мы используем локальную регрессию (Cleveland & Devlin, 1988) или изотоническую регрессию (Dykstra, 1981) в качестве гибких непараметрических моделей кривых ошибок. Чтобы сгенерировать исторические ошибки прогноза и подогнать эту модель, мы используем процедуру, которую называем смоделированные исторические прогнозы.

4.3 Смоделированные исторические прогнозы

Мы хотели бы обучить модель ожидаемой ошибки, описанную в (8), чтобы выполнить отбор и оценку модели. К сожалению, использование таких методов, как кросс-валидация, затруднено, поскольку наблюдения не являются взаимозаменяемыми – мы не можем просто случайным образом разделить данные.

Мы применяем симулированные исторические прогнозы (SHFs — Simulated Historical Forecasts), создавая [math]K [/math] прогнозов в различных точках отсечения в истории, выбранных таким образом, чтобы горизонты находились в пределах истории и можно было оценить суммарную ошибку. Эта процедура основана на классических методах оценки прогнозов с «скользящим началом» (Tashman 2000), но использует лишь небольшую последовательность дат отсечения вместо одного прогноза для каждой исторической даты. Основное преимущество использования меньшего количества симулированных дат (при оценке со скользящим началом получается один прогноз на дату) заключается в экономии вычислительных ресурсов и снижении корреляции между измерениями точности.

SHF имитируют ошибки, которые мы допустили бы, если бы использовали этот метод прогнозирования в те моменты в прошлом. Прогнозы, представленные на рисунках 3 и 4, являются примерами SHF. Этот метод обладает преимуществом простоты, его легко объяснить аналитикам и лицам, принимающим решения, и он относительно не вызывает споров при получении представления об ошибках прогнозирования. Существует два основных момента, на которые следует обратить внимание при использовании методологии SHF для оценки и сравнения подходов к прогнозированию.

Во-первых, чем больше симулированных прогнозов мы делаем, тем более коррелированными оказываются их оценки ошибки. В крайнем случае, если мы делаем симулированный прогноз для каждого дня в истории, изменения между прогнозами, основанными на дополнительном дне информации, будут незначительными, и ошибки от одного дня к другому будут практически идентичными. С другой стороны, если мы делаем очень мало симулированных прогнозов, то у нас будет меньше наблюдений за историческими ошибками прогнозирования, на которых можно основывать выбор модели. Как эвристическое правило, для горизонта прогноза [math]H[/math] мы обычно делаем симулированный прогноз каждые [math]H/2[/math] периодов. Хотя коррелированные оценки не вводят смещения в оценку точности модели, они дают менее полезную информацию и замедляют процесс оценки прогноза.

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

На рисунке 7 представлены наши оценки функции [math]\xi(h)[/math], то есть ожидаемой средней абсолютной процентной ошибки на всём горизонте прогноза, полученные с помощью [math]LOESS[/math] для временных рядов, показанных на рисунках 3 и 4. Оценка была сделана, исходя из девяти симулированных дат прогноза, по одной на квартал, начиная после первого года. Модель Prophet показывает меньшую ошибку прогноза на всех горизонтах. Прогнозы Prophet строились с настройками по умолчанию; при этом более тонкая настройка параметров может ещё больше повысить точность.

При визуализации прогнозов мы предпочитаем использовать точки, а не линии, чтобы отображать исторические данные, поскольку они представляют собой точные измерения, которые никогда не интерполируются. Затем поверх точек мы накладываем линии, отражающие прогнозы. Для SHF полезно визуализировать ошибки, которые модель допустила на разных горизонтах, как в виде временного ряда (как на Рисунке 3), так и в агрегированном виде по всем SHF (как на Рисунке 7).

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

4.4 Выявление больших ошибок прогноза

Когда мы получаем большое количество прогнозов в ручную проверить представляется невозможным то, важно иметь возможность автоматически выявлять проблемные прогнозы. Автоматическая идентификация некачественных прогнозов позволяет аналитикам наиболее эффективно использовать своё ограниченное время и применять свой опыт для устранения выявленных проблем. Существует несколько способов использования SHF для обнаружения потенциальных проблем с прогнозами:

• Если прогноз имеет значительные ошибки по сравнению с базовыми моделями, это может указывать на неправильную спецификацию модели. Аналитики могут скорректировать модель тренда или сезонность в зависимости от необходимости.
• Если для всех методов в определённую дату наблюдаются большие ошибки, это может свидетельствовать о наличии выбросов. Аналитики могут выявить и удалить выбросы.
• Если ошибка SHF для какого-либо метода резко увеличивается от одного момента отсечения к следующему, это может указывать на изменение процесса генерации данных. В таком случае добавление точек изменения или моделирование различных фаз отдельно может помочь устранить проблему.

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

5. Вывод

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

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

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

6. Благодарности

Мы благодарим Дэна Мерла за содействие в разработке Prophet, а также за предложения и идеи на всех этапах её создания. Мы также выражаем благодарность Дирку Эдделбюттелю, Даниэлю Каплану, Робу Хайндману, Алексу Гилгуру и Ладе Адамич за полезные рецензии этой работы. Особенно благодарим Роба Хайндмана за ценные идеи, связавшие нашу работу с экспертными прогнозами.

1. Модели ARIMA способны включать сезонные ковариации, но добавление этих ковариат приводит к чрезвычайно длительному времени подгонки и требует опыта моделирования, которого нет у многих новичков в прогнозировании. ↩︎
2. Термин свободного члена (intercept) можно опустить, поскольку одновременно подбирается компонент тренда. ↩︎