Всем привет!

Как и во многих других компаниях, в X5 существует огромное количество данных, зависящих от времени. Такие данные принято называть временными рядами (time-series). Это могут быть данные о продажах в магазинах, об остатках на складах или об удовлетворенности клиентов. Используя эти данные, мы хотим искать инсайты и приносить пользу бизнесу.

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

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

В нашем обзоре мы рассмотрим, как различные модификации метода бутстрапа учитывают структурные особенности и зависимости в данных временных рядов. Особое внимание будет уделено объяснению, почему нельзя применять стандартный подход бутстрапа к временным рядам без учёта их структуры. Затем мы перейдем к обзору методов, которые позволяют эффективно решить эту проблему.

Временной ряд и его свойства

Общие понятия

Временной ряд (time‑series) — это собранный в разные моменты времени статистический материал о значении каких‑либо параметров (в простейшем случае одного) исследуемого процесса. Каждая единица статистического материала называется измерением или точкой данных. Во временном ряде для каждого отсчёта должно быть указано время измерения или номер измерения по порядку:

Y = Y_{t} : t \in T

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

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

Выявление компонентов временного ряда (time-series decomposition) предполагает его разложение на тренд (trend), сезонность (seasonality) и остатки (residuals).

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

Сезонность (seasonality) отражает периодические колебания в данных. Эти колебания могут быть связаны с сезонными факторами, такими как времена года, месяцы или дни недели. Например, продажи мороженого могут иметь сезонность с пиком в летние месяцы и спадом зимой.

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

Пример декомпозиции временного ряда продаж магазинов из открытых данных Facebook Prophet
Пример декомпозиции временного ряда продаж магазинов из открытых данных Facebook Prophet

Временные ряды могут обладать разными свойствами. Остановимся на некоторых из них подробнее:

Строгая стационарность (strict stationarity). Пусть (X_t, ; t \in T) будет стохастическим (случайным) процессом, то есть семейством случайных величин,
индексируемых параметром t, изменяющихся в индексном множестве T.
F_X(x_{t_1+\tau}, \ldots, x_{t_n+\tau})представляет собой функцию совместного распределения {X_t} в моменты времени t_1 + \tau, \ldots, t_n + \tau. Тогда {X_t}называется строго стационарным, если
F_X(x_{t_1+\tau}, \ldots, x_{t_n+\tau}) = F_X(x_{t_1}, \ldots, x_{t_n}), ; для всех  \tau, t_1, \ldots, t_n \in T
и для всех n \in \mathbb{N}_0.

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

Строго стационарный ряд
Строго стационарный ряд

Слабая стационарность (weak stationary) или просто стационарность. Для временного ряда X_tс дискретными моментами времени t, можно сказать, что ряд считается слабо стационарным, если выполняются следующие условия:

  1. Математическое ожидание (\mu)временного ряда не зависит от времени:
    E[X_t] = \muдля всех t

  2. Дисперсия (\sigma^2)временного ряда не зависит от времени:

    Var[X_t] = \sigma^2 \leq \infty
  3. Автоковариационная функция Cov(x_t, x_{t + s}) = \gamma_s зависит только от длины интервала между моментами времени (t_1)и (t_2):

    Cov(x_t, x_{t + s}) = Cov(x_{t + k}, x_{t + k + s}) = \gamma_s, ; k \in T
Слабостационарный ряд

Свойства статистических оценок

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

    \mathbb{E}(\hat{\theta}) = \theta

    Если \hat{\theta} — оценка параметра \theta, то оценка называется несмещённой, если:

    где \mathbb{E}(\hat{\theta}) обозначает математическое ожидание (среднее) оценки \hat{\theta}, а \theta - истинное значение параметра.

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

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

    Если \hat{\theta} — оценка параметра \theta, то оценка называется состоятельной, если:

    \lim_{n \to \infty} \text{P}(|\hat{\theta} - \theta| < \epsilon) = 1, ; \forall \epsilon > 0

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

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

    Если \hat{\theta} — оценка параметра \theta, то оценка \hat{\theta} называется эффективной, если для любой другой несмещенной оценки \tilde{\theta} параметра \theta выполняется неравенство Крамера-Рао:

    Var(\hat{\theta}) \leq \text{Var}(\tilde{\theta})

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

Почему это важно?

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

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

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

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

Bootstrap

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

На входе у нас есть выборка x_1, x_2,  \ldots , x_N \in X , где N— это размер выборки.

Для b = 1, 2,3, \ldots ,B повторяем следующее (B это параметр, указывающий, сколько псевдо выборок подлежат к рассмотрению):

  1. Семплируем с возвращением точки x_1^{*(b)}, x_2^{*(b)},x_3^{*(b)}, \ldots , x_N^{*(b)} из того набора данных, который мы пронаблюдали (эмпирического распределения). Назовем эту выборку бутстрапным набором данных.

  2. Посчитаем значение искомого параметра \bar A_b^{*} = A(x_1^{*(b)}, x_2^{*(b)},x_3^{*(b)}, \ldots , x_N^{*(b)}) и где-то сохраним результат.

  3. Произведем шаг 1. и шаг 2. B раз.

Как итог, мы получаем новую выборку  \bar A_1^{}, \bar A_2^{}, \ldots ,\bar A_b^{*}размераB, состоящую из оценок искомого параметра.

На ее основе мы можем судить о доверительном интервале, в котором находится истинное значение статистики с заданной доверительной вероятностью. Обычно эту вероятность (уровень доверия) берут на уровне 95\%. Если говорить строго, то:

\theta — оцениваемый параметр (например среднее), L — нижняя граница интервала, U — верхняя граница, а 1 - \alpha — уровень доверия (где \alpha — уровень значимости), тогда доверительный интервал для \theta можно записать как:

P(L < \theta < U) = 1 - \alpha

где Pобозначает вероятность, что параметр находится в интервале между случайными величинами L и U.

И как будто на этом можно было бы и закончить, однако, этот метод работает корректно только с выборками, которые независимо и одинаково распределены (i.i.d.).

Из статьи про бутстрап
Из статьи про бутстрап

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

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

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

Постановка задачи

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

Пусть у нас есть временной ряд, полученный из авторегрессионной AR(2) модели вида:

x_t = \phi_1 x_{t-1} + \phi_2 x_{t-2} + \epsilon_t, \; \epsilon_t \sim \mathcal{N}(0, 1), \; \phi_1 = 0.2, \; \phi_2 = 0.4, \; N = 300

Представим, что мы знаем вид модели, из которой мы получили этот ряд. Мы хотим сделать оценку параметров \hat{\phi}_1 и \hat{\phi}_2 модели на основе этого ряда.

Рассматриваемый временной ряд
Рассматриваемый временной ряд

Воспользуемся тем, что мы знаем вид модели и ее параметры, сгенерируем в дополнение к этому ряду еще B - 1 временных рядов, где B = 1000 (параметр равен числу бутстрап итераций). Мы будем рассматривать эту выборку как генеральную совокупность (обозначим ее как True), к которой у нас нет доступа в реальной жизни, но есть доступ в рамках симуляции. Оценки параметров на основе этих данных будут являться эталоном для сравнения разных методов бутстрапа.

Block bootstrap

Блочный бутстрап является наиболее общим методом для улучшения точности бутстрапа для временных рядов.

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

Для того чтобы реализовать любой из методов блочного бутстрапа, нам необходимо выбрать (ожидаемую) длину блока. Оптимальная длина блока зависит от размера выборки и ее корреляционной структуры, и различается для разных методов блочного бутстрапа. Существуют два основных подхода к выбору длины блока: метод перекрестной проверки и \emph{метод подстановки}. В методе перекрестной проверки мы выбираем критерий (например, среднеквадратичную ошибку) и минимизируем оцененный критерий, чтобы получить оценку оптимального размера блока. В методе подстановки мы сначала выражаем оптимальный размер блока через модельные параметры, а затем подставляем оценки этих параметров для получения конкретного значения размера блока. Метод перекрестной проверки обычно требует меньше аналитической работы, но требует больше вычислений.

Идея блочного бутстрапа лучше всего иллюстрируется на примере. Предположим у нас есть временной ряд:

X_1, X_2, X_3, X_4, X_5, X_6, X_7, X_8, X_9, X_{10}

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

\overbrace{X_1, X_2}^{\text {block } 1}, \overbrace{X_3, X_4}^{\text {block } 2}, \overbrace{X_5, X_6}^{\text {block3 }}, \quad \overbrace{X_7, X_8}^{\text {block } 4}, \overbrace{X_9, X_{10}}^{\text {block5 }}.

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

block3, \quad block 1, \quad block 5 , \quad block 2, \quad block 5

или, с точки зрения исходного ряда,

X_5, X_6, \quad X_1, X_2, \quad X_9, X_{10}, \quad X_3, X_4, \quad X_9, X_{10}

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

Один из ранних вариантов блочного бутстрапа, предложенный Карльштейном (1986), ограничивает внимание на наборе непересекающихся блоков данных и формирует выборку из этого набора для генерации бутстрап-наблюдений. Это метод также называют блочным бутстрапом без перекрытия (NBB). Чтобы кратко описать его, предположим, что \ell — целое число в промежутке (1, n). Также, для упрощения, предположим, что \ell делит n и установим b=n / \ell. Бутстрап-наблюдения NBB генерируются путем случайного выбора b блоков с возвращением из коллекции \left\{\tilde{\mathbb{B}}_1, \ldots, \tilde{\mathbb{B}}_b\right\}, где

\begin{aligned}    & \tilde{\mathbb{B}}_1=\left(X_1, \ldots, X_{\ell}\right), \\    & \tilde{\mathbb{B}}_2=\qquad\qquad\qquad\left(X_{\ell+1}, \ldots, X_{2 \ell}\right), \\  & \tilde{\mathbb{B}}_b= \qquad\qquad\qquad\qquad\qquad\qquad\qquad\left(X_{(b-1) \ell+1}, \ldots, X_n\right) \text {. } \\\end{aligned}
Иллюстрация метода NBB
Иллюстрация метода NBB

Оценка методом non-overlapping bootstrap может быть лучше, чем оценка методом block bootstrap, если данные имеют слабую или среднюю зависимость и длина блока достаточно большая, чтобы уловить структуру зависимости.

NBB имеет следующие свойства:

  1. Несмещенность: NBB является несмещенным методом оценки, если длина блока b стремится к бесконечности при увеличении размера выборки n. Это означает, что среднее значение бутстрапных оценок сходится к истинному значению статистики при n \rightarrow \infty и b \rightarrow \infty.

  2. Состоятельность: NBB является состоятельным методом оценки, если отношение b/n стремится к нулю при увеличении n. Это означает, что распределение бутстрапных оценок сходится к истинному распределению статистики при n \rightarrow \infty и b/n \rightarrow \infty.

  3. Эффективность: NBB считается эффективным методом оценки, когда он обеспечивает наименьшую дисперсию оценки по сравнению с другими методами бутстрапа для конкретной статистики. Это достигается через выбор оптимальной длины блока b, которая наилучшим образом учитывает структуру зависимости в данных. Для определения оптимальной длины блока используются различные критерии, включая правило Холла-Хоровица, правило Политиса-Романо и правило Люнгбокса.

Оценка параметров временного ряда c использованием NBB
Оценка параметров временного ряда c использованием NBB

Однако, несмотря на многие преимущества метода, одна из обычных проблем при использовании NBB заключается в том, что общее количество блоков часто оказывается недостаточным. В результате бутстрап-выборки не отражают полную статистическую вариативность данных, что может привести к некорректным статистическим выводам. Для устранения этой проблемы можно применить метод moving-block bootstrap.

Moving block bootstrap

Для применения метода moving-block bootstrap (MBB) сначала нужно определить размер блока, который будет использоваться для формирования новых выборок. Затем исходная выборка разбивается на несколько блоков равного размера. Каждый блок перемещается по выборке с некоторым шагом и используется для формирования новой выборки путем случайного выбора элементов из данного блока с возвращением. Таким образом, каждая новая выборка будет иметь тот же размер и ту же структуру блоков, что и исходная выборка. Затем, на основе каждой новой выборки вычисляется интересующий нас статистический параметр. Повторяя этот процесс множество раз, мы получаем распределение данного параметра.

Иллюстрация moving-block bootstrap
Иллюстрация moving-block bootstrap

Использование метода moving-block bootstrap позволяет учесть зависимость между данными в выборке и получить более реалистичные оценки распределения параметров. Вспомним, что блоки в конструкции NBB не перекрываются, из-за этого теоретические свойства оценок NBB легче анализировать, чем свойства оценок MBB. Однако оценки NBB обычно имеют более высокую среднеквадратическую ошибку при любом размере блока \ell по сравнению с MBB.

Свойства оценки методом MBB зависят от выбора длины блока \ell и количества повторений B. В общем случае, можно сказать следующее:

  1. Несмещённость: Оценка методом MBB будет несмещенной, если длина блока \ell увеличивается с ростом объема выборки n таким образом, что \ell растет медленнее, чем квадратный корень из n.

  2. Состоятельность: Оценка методом MBB будет состоятельной, если длина блока \ell увеличивается быстрее, чем логарифм от объема выборки n.

  3. Эффективность: Оценка методом MBB считается эффективной, когда длина блока \ell выбирается таким образом, чтобы наилучшим образом учесть характеристики временного ряда. Для определения оптимальной длины блока могут применяться различные методы, включая анализ политропического спектра или выбор на основе минимизации среднеквадратичной ошибки.

Оценка параметров временного ряда c использованием MBB
Оценка параметров временного ряда c использованием MBB

Circular block bootstrap

Цель метода Circular block bootstrap (CBB) — это устранение эффекта неравномерного взвешивания наблюдений на краях, который возникает в методе MBB. В MBB, например, первое значение выборки появляется только в первом блоке, в то время как пятое значение будет встречаться в пяти блоках, если b > 5. В CBB блоки формируются аналогично MBB, с перекрытием, но дополнительно используется 'циклическое' обертывание данных, когда конец последовательности соединяется с её началом, позволяя таким образом сформировать дополнительные блоки. Этот метод особенно ценен для данных с круговой природой, как, например, годовые циклы температур.

Иллюстрация circular block bootstrap
Иллюстрация circular block bootstrap

"Обернутый круг" определяется как X_i \equiv X_j для i>n, где j=i \bmod n, и X_0=X_n.

Как упоминалось выше, пусть ряд будет иметь вид \left\{X_1, \ldots, X_n\right\}, тогда весь дуговой сегмент "обернутого круга" будет \left\{X_1, \ldots, X_n, X_1, \ldots, X_{n+l-1}\right\}; таким образом, установив длину блока равной l, получим в целом n блоков, и каждый блок определяется как Y_i=\left\{X_i, \ldots, X_{i+l-1}\right\}.

Множество блоков представлено как \left\{Y_1, \ldots, Y_n\right\}. С вероятностью 1 / n производится выборка с возвращением, случайным образом выбираются b независимых и одинаково распределенных блоков из этого множества, где l b=m \approx n. Затем у нас появляется новая последовательность X_1^*, \ldots, X_m^*. С использованием этой новой последовательности также можно рассчитать доверительный интервал.

Свойства оценки методом CBB зависят от длины блока, которая должна быть достаточной, чтобы учитывать корреляцию между наблюдениями. В общем случае, оценка методом CBB обладает следующими свойствами:

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

  2. Состоятельность: оценка методом CBB является состоятельной, если длина блока стремится к бесконечности со скоростью большей, чем логарифм от объема выборки. Это условие гарантирует, что вероятность того, что блок содержит два исходных наблюдения, стремится к нулю.

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

Преимущества метода CBB по сравнению с NBB и MBB заключаются в следующем:

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

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

  3. Метод CBB подходит для данных с периодической или почти периодической структурой, так как он сохраняет эту структуру при оборачивании данных на круг. Методы NBB и MBB могут нарушить эту структуру при разбиении данных на блоки.

Оценка параметров с использованием CBB
Оценка параметров с использованием CBB

Stationary bootstrap

Стационарный бутстрап (SB) отличается от других трёх методов, упомянутых ранее, тем, что длина блока не фиксирована, а является случайной величиной с некоторым распределением (обычно геометрическим) и ожидаемым значением l. Из-за случайной длины блока количество блоков также является случайным. Подобно CBB, метод SB был разработан с целью устранения неравномерного взвешивания наблюдений в начале и в конце ряда. Дополнительный источник случайности возникает из выбора точек начала блоков. Это аналогично случайному выбору блоков в методе MBB, поскольку каждая точка данных в последовательности является началом блока в MBB.

Предположим, что \left\{X_t\right\} - заданный слабо стационарный временной ряд, как уже упоминалось ранее. Определим каждый блок как B_{i, l}=\left\{X_i, X_{i+1}, \cdots, X_{i+l-1}\right\}, где l - длина блока. Поскольку, когда мы случайным образом выбираем начальную точку каждого блока и длина блока имеет геометрическое распределение, наблюдения в некоторых блоках могут выходить за пределы данной выборки. Чтобы гарантировать, что все наблюдения будут включены в блоки и не будет 'обрезания' данных на краях, можно использовать технику 'обертывания' ряда. 'Обертывание' означает, что конец временного ряда соединяется с его началом, образуя замкнутую циклическую структуру:

X_i \equiv X_j\;\text{для}\; i>n,\;\text{где}\;j=i(\bmod n),\;\text{и}\;X_0=X_n

Поскольку длины блоков могут быть разными, определим l_1, l_2, \cdots как длину блока, независимую от X_i, и \left\{l_i\right\} — последовательность независимых одинаково распределенных случайных величин, имеющих геометрическое распределение:

P\left\{l_i=u\right\}=\left(1-p_S\right)^{u-1} p_S \;\text{для}\;u=1,2, \ldots.

Тогда для любого заданного p_S мы можем сгенерировать последовательность \left\{l_i\right\}. Что касается индекса начальной точки каждого блока i_i, независимого от X_i и l_i, случайно выбираем из \{1, \ldots, n\}, то есть \left\{i_i\right\} — последовательность независимых одинаково распределенных случайных величин, имеющих дискретное равномерное распределение на \{1, \cdots, n\}. Таким образом, выборочная последовательность блоков с случайной длиной будет иметь следующий вид:

\left\{B_{i_1, l_1}, B_{i_2, l_2}, \ldots\right\}=\left\{X_{i_1}, X_{i_1+1}, \ldots, X_{i_1+l_1-1}, X_{i_2}, X_{i_2+1}, \ldots, X_{i_2+l_2-1}, \ldots\right\}

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

Свойства оценки методом stationary bootstrap зависят от выбора параметра p, который определяет среднюю длину блока. В общем случае, оценка методом stationary bootstrap обладает следующими свойствами:

  1. Несмещенность: оценка методом stationary bootstrap является несмещенной, если p стремится к нулю с увеличением объема выборки n. Это означает, что среднее значение оценки сходится к истинному значению параметра.

  2. Состоятельность: оценка методом stationary bootstrap является состоятельной, если p удовлетворяет условию np \rightarrow \infty \text { при } n \rightarrow \infty. Это означает, что дисперсия оценки уменьшается с увеличением объема выборки и стремится к нулю.

  3. Эффективность: оценка методом stationary bootstrap является асимптотически эффективной, если p выбирается таким образом, чтобы минимизировать среднеквадратичную ошибку оценки. Оптимальное значение p зависит от конкретной статистики и степени зависимости данных. В некоторых случаях, оптимальное значение p может быть получено из данных с помощью метода максимального правдоподобия или кросс-валидации.

Оценка параметров с использованием SB
Оценка параметров с использованием SB

Model-Based Bootstrap

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

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

Пайплайн Model-Based бутстрапа
Пайплайн Model-Based бутстрапа

Мы предполагаем, что существует линейная зависимость между вектором x и переменной отклика y.

\epsilon— последовательность независимых одинаково распределенных случайных величин (i.i.d.) нормально распределенных с нулевым средним и стандартным отклонением \sigma, т.ч. y = x\beta + \epsilon, \;\; \epsilon \sim N(0, \sigma^2) . Обучая модель получим параметр \hat\betaи остаток \epsilon:

\begin{equation*} \begin{cases} \hat\beta = \underset{\beta}{argmin} \sum^n_{i=1} (y_i - x\beta)^2 \\ \epsilon = y - x\hat\beta \end{cases} \end{equation*}

Мы получаем параметр \hat\beta путем выбора такого \beta, что квадрат разницы между наблюдаемыми и предсказанными значениями будет минимален. Далее вычитаем предсказание из значений, и таким образом получаем остатки. Исследуем данные остатки, чтобы убедиться, что они соответствуют нормальному распределению, и оценим дисперсию шума:

\sigma^2 = \frac{RSS(\hat\beta)}{n - 1}

Исходя из того, что модель адекватно описывает данные, мы создаём новые синтетические наборы данных y^{*}, используя параметры и результаты диагностики исходной модели. Далее, на этих синтетических данных обучаем модель, получая новые оценки параметров \hat\beta^{*}. Этот процесс повторяется многократно, создавая бутстрап-выборки, что позволяет нам сформировать распределение оценок для каждого параметра модели. Таким образом, на этапе генерации синтетических данных и последующего их анализа осуществляется бутстрап.

Вариация этих оценок затем используется для аппроксимации стандартной ошибки SE исходной оценки ( \hat{\beta}).

SE(\hat\beta) = \sqrt{Var(\hat\beta^{*})}

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

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

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

Оценка параметров с использованием MB из постановки задачи
Оценка параметров с использованием MB из постановки задачи

Frequency Domain Bootstrap

Для этого метода перейдём от временного домена (точки во времени) к частотному домену при помощи быстрого преобразования Фурье:

Для справки:

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

I_n(\lambda) = \frac{1}{2\pi n} \bigg | \sum_{t=1}^n X_t e^{-i\lambda t} \bigg | ^ 2, \;\; \lambda \in [-\pi, \pi]

Такой переход позволяет использовать сэмплирование периодограммы (например, предположив, что частоты независимы).

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

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

Посмотрим на реальный пример такого преобразования:

Пример преобразования Фурье временного ряда
Пример преобразования Фурье временного ряда

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

Оценка параметров с использованием FD из постановки задачи
Оценка параметров с использованием FD из постановки задачи

Сравнения методов

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

Мы сгенерировали данные при помощи авторегрессионного процесса второго порядка, взяв \phi_1 = 0.2 и \phi_2 = 0.4 при помощи библиотеки statsmodels.api, а основные методы бутстрапа были реализованы в библиотеке arch.bootstrap.

Классический бутстрап метод не позволяет сделать оценку параметров, если у нас есть созависимость данных, поскольку он не учитывает структуру ряда при построении бутстрап семпла
Классический бутстрап метод не позволяет сделать оценку параметров, если у нас есть созависимость данных, поскольку он не учитывает структуру ряда при построении бутстрап семпла
При размере блока = 1 блочные методы равносильны обычному бутстрапу
При размере блока = 1 блочные методы равносильны обычному бутстрапу
При низком значении автоковариации поведение методов и классического бустрапа совпадает
При низком значении автоковариации поведение методов и классического бустрапа совпадает
Усредненный MSE для различных методов бутстрапа
Усредненный MSE для различных методов бутстрапа

Для разных значений параметров \hat{\phi}_1 и \hat{\phi}_2 мы проведели серию из 500 симуляций для определения усредненной среднеквадратической ошибки (MSE) различных методов бутстрапа. Результаты показали, что методы NBB и CBB обеспечивают наиболее низкую MSE для обоих параметров модели, что указывает на их высокую эффективность в оценке. Однако результаты других четырех методов: SB, ModelBased, MBB, FrequencyDomain не сильно хуже. А вот метод BM показал значительно большие значения MSE, подчеркивая низкую точность в данной задаче оценки.

Распределение MSE для различных методов бутстрапа
Распределение MSE для различных методов бутстрапа

Методы ModelBased и FrequencyDomain демонстрируют меньшую дисперсию распределения MSE, что свидетельствует о их высокой точности и низкой вариативности оценок параметров модели.

NBB, MBB, SB и CBB показывают большую вариативность в оценках, с MSE, достигающими более высоких значений, особенно для параметра \phi_2. Однако наиболее частое значение MSE этих методов находится ближе к 0.

BM показал себя как метод с наибольшей вариативностью и наибольшими значениями MSE, что делает его наименее подходящим для точной оценки параметров модели AR(2).

Заключение

Вот и конец!

Мы рассмотрели основные методы для бутстрапа временных рядов. Выводы, к которым мы пришли по итогам работы:

  • Когда временной ряд не соответствует предположению о независимости и одинаковом распределении наблюдений (i.i.d), становится актуальным использование блочных методов бутстрапа. Эти методы программно реализованы и позволяют учитывать структурные особенности данных, хотя они чувствительны к выбору размера окна.

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

  • Frequency Domain Bootstrap представляет собой экспериментальный подход, который позволит выявлять сложные паттерны, которые трудно определить во временной области, хотя его применение ограничивается из-за отсутствия эффективных программных реализаций.

Проведённые симуляции для MSE различных методов бутстрапа показали, что NBB и CBB наиболее эффективны, а методы ModelBased и FrequencyDomain более стабильны и точны в оценках параметров модели. Однако мы советуем опираться на свойства ваших данных, и подбирать метод исходя из характеристик конкретно вашей задачи.

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

Над статьёй работали члены команды Ad-Hoc X5 Tech:

Комментарии (0)