Обнаружение аномалий – это поиск непредвиденных значений или шаблонов последовательностей значений в потоках данных.

Аномалия в этом контексте – это отклонение поведения системы от ожидаемого.

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

Виды аномалий

Аномалии обычно относят к одному из трёх типов¹:

• Точечные аномалии
• Контекстные аномалии
• Коллективные аномалии

Точечные аномалии

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

Контекстные аномалии

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

Коллективные аномалии

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

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

Обнаружение точечных аномалий

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

Распознавание аномалий с учителем

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

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

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

Основная сложность таких методик – формирование данных для обучения. Кроме очевидных трудозатрат, часто аномальный класс к тому же хуже представлен, чем нормальный, что может приводить к неточностям в полученной модели.

Распознавание аномалий частично с учителем

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

Распознавание без учителя

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

Методы распознавания аномалий

Классификация

Данный метод основан на том, что нормальное поведение системы может определяться одним или несколькими классами. Тогда экземпляр, не принадлежащий ни одному из классов, является аномальным. Этот метод обычно использует подход “частично с учителем”. Основные механизмы: нейронные сети, Байесовские сети, метод опорных векторов, на основе правил.

Кластеризация

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

• Нормальные экземпляры данных относятся к какому-то кластеру
• Нормальные данные ближе к центру кластера, аномальные – дальше
• Нормальные данные образуют большие плотные кластеры, а аномальные – маленькие и разрозненные.

Один из простейших алгоритмов на основе кластеризации – алгоритм k-means.

Статистический анализ

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

• Параметрические методы. Предполагается, что нормальные данные имеют функцию плотности вероятности \(ρ(x,θ)\), где \(θ\) – вектор параметров, \(x\) – экземпляр данных (наблюдение).

• Непараметрические методы. Структура модели не определена a priori, а определяется из предоставленных данных.

Алгоритм ближайшего соседа

При использовании этой методики вводится метрика (мера похожести между объектами). Тогда возможны два подхода:

• Расстояние до k-го ближайшего соседа. Аномальные данные наиболее удалены от всех других данных.

• Использование относительной плотности. Экземпляры в областях с низкой относительной плотностью оцениваются как аномальные. (напр. метод локального уровня выброса)

Спектральные методы

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

Обнаружение аномалий во временных рядах

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

Оконные методы

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

В самом обобщённом описании, берутся обучающие временные ряды \(S_i\), из каждого из них извлекается \(p\) окон \(s_i^k\), аналогично из тестовых временных рядов \(T_j\) извлекается \(r\) окон \(t_j^l\).

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

К недостаткам оконных методов можно отнести чувствительность к размеру окна (аномалии должны укладываться в окно), и относительную вычислительную “дороговизну”. Действительно, каждое тестовое окно должно сравниваться с каждым обучающим окном, поэтому оценка сложности не лучше, чем \(O(nprw)\) для каждого тестового ряда, где \(w\) – размер окна, \(n\) – число обучающих рядов.

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

Прямые метрические методы

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

Два основных подхода, используемых в данном случае это

• Мера аномальности – расстояние от тестового ряда до k-го ближайшего соседа в обучающем наборе.
• Кластеризация. Мера аномальности – расстояние от тестового ряда до центра ближайшего кластера.

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

В таком случае применяется алгоритм динамической трансформации временной шкалы (Dynamic time wrapping, DTW). Идея алгоритма в том, чтобы минимизировать отличия по оси значений через трансформацию по оси времени. Следует однако учитывать, что DTW может давать некорректные результаты: например, выровнять ряды так, что одной точке одного ряда ставится в соответствие большая группа точек второго ряда. DTW в общем имеет квадратичную сложность, но есть линейные алгоритмы, такие как fastDTW, вычисляющие “достаточно хорошее приближение”.

DTW хуже работает для несинхронизированных периодических временных рядов. Для таких случаев анализ автокорреляционной функции позволяет найти задержку и выровнять ряды.

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

Методы на основе прогнозирования

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

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

В качестве прогностической модели могут использоваться спектральные модели, авторегрессионные модели типа ARIMA и т.п.

Прогностические методы могут захватывать все типы аномалий.

Алгоритмы на основе скрытых марковских моделей

Основное предположение СММ – наблюдаемый ряд \(O=\{o_1,o_2,\ldots,o_n\}\) является косвенным наблюдением “скрытого” ряда \(Q=\{q_1,q_2,\ldots,q_n\}\), причём процесс, описывающий ряд \(Q\) является марковским, т.е. описывается авторегрессионной моделью первого порядка AR(1): \(X_{t}=c+\alpha X_{t-1}+\varepsilon _{t}\)

Для построения модели СММ на основе обучающих данных используется алгоритм Баума-Велша. Затем для найденных параметров скрытой модели, вычисляется вероятность порождения данного тестового ряда моделью с такими параметрами.

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

Алгоритм k-means

Алгоритм k-means на основе k средних значений выделяет k кластеров.

Наивный алгоритм работает следующим образом. Если даны начальные значения средних \(m_1, \ldots, m_k\), то:

1. Шаг присваивания. Каждое наблюдение \(x_i\) присваивается кластеру с ближайшим целым: \[S_i^{(t)} = \{x_p : \lVert x_p-m_i^{(t)}\rVert^2 \le \lVert x_p-m_j^{(t)}\rVert^2,\, \forall p, \forall j \in [1,k]\}\]

   При этом каждый \(x_p\) присваивается только одному кластеру.

2. Шаг обновления. Средние значения (центры) кластеров пересчитываются с учётом нового распределения: \[m_i^{(t+1)} = \frac{1}{|S_i^{(t)}|}\sum_{x\in S_i^{(t)}} x\]

3. Повторять 1, 2 пока кластеры не перестанут изменяться.

Для инициализации обычно используются подходы Forgy или случайного разделения.

Forgy случайным образом выбирает k наблюдений и использует их в качестве начальных.

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

Более сложный алгоритм выбора начальных средних, дающий лучшие результаты – алгоритм k-means++:

1. Случайно (с равномерным распределением) выбрать одно наблюдение в качестве первого среднего
2. Для каждого из оставшихся наблюдений \(x\), вычислить \(D(x)\) – расстояние между \(x\) и ближайшим средним, из уже выбранных.
3. Выбрать новое наблюдение в качестве нового среднего с вероятностью, пропорциональной \(D(x)^2\).
4. Повторять шаги 2, 3 пока не выбрано \(k\) средних.

Хотя k-means++ требует дополнительных затрат времени на генерацию начального приближения, он значительно уменьшает число итераций, необходимых для основного алгоритма.

Кроме того, k-means++ имеет лучшие оценки оптимальности распределения (k-means не обязательно сходится к глобальному оптимуму, k-means++ чаще сходится к глобальному оптимуму, чем стандартный)

Метод опорных векторов

Каждый объект данных (наблюдение) представляется как точка в \(p\)-мерном пространстве. Каждая точка принадлежит одному из двух классов. Вопрос: можно ли разделить эти точки гиперплоскостью размерности \(p-1\). Это задача линейной разделимости.

Таких гиперплоскостей, если есть хотя бы одна, может быть бесконечно много, поэтому полагают, что оптимальная гиперплоскость – та, которая максимизирует расстояние до ближайших точек с обеих сторон (“зазор”).

Формально, мы полагаем, что имеется обучающий набор вида \[\{(\vec x_1, c_1), \ldots, (\vec x_n, c_n)\},\] где \(c_i \in \{-1, 1\}\) показывает, к какому классу принадлежит точка \(\vec x_i\). Каждый \(\vec x_i\) представляет собой \(p\)-мерный вектор, обычно нормализованный.

Чтобы классифицировать новые точки аналогично прочим, строится гиперплоскость, описываемая уравнением \[\vec w \cdot \vec x - b = 0,\] где \(\vec w\) – перпендикуляр к разделяющей гиперплоскости, а \(\frac{b}{|\vec w|}\) – расстояние от начала координат до гиперплоскости.

Размер зазора – это расстояние между гиперплоскостями, параллельными оптимальной, описываемыми уравнениями \[\vec w \cdot \vec x - b = 1\] \[\vec w \cdot \vec x - b = -1\] (с точностью до нормировки).

Тогда зазор \(d = \frac{2}{|\vec w|}\), и чтобы его максимизировать необходимо минимизировать \(|\vec w|\).

Также, чтобы в зазоре не было ни одной точки, должно выполняться условие \[c_i(\vec w \cdot \vec x_i - b) \ge 1\]

Это задача ограниченной оптимизации.

На практике, обычно решается “смягчённая” задача, которая смягчает требования неравенства и вводит в функцию минимизации “штраф”:

\[|\vec w|^2 + C\sum_{i=1}^n\xi_i \to \min_{\vec w, b, \vec \xi}\] \[c_i (\vec w\cdot \vec x_i - b) \ge 1 -\xi_i\] \[\xi_i \ge 0\]

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

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

\[x(t) = \{x(t-E+1), x(t-E+2), \ldots, x(t)\}^{\mathrm T},\]

где \(E\) – размерность пространства (число элементов вектора). Этот подход близок по сути к оконным методам.

Алгоритм Баума-Велша

Скрытая модель Маркова – вертоятностная модель множества случайных переменных \(\{Y_1,\ldots,Y_t, Q_1, \ldots, Q_t\}\). Переменные \(Y_i\) – наблюдения, переменные \(Q_i\) – скрытые величины. Два основных предположения:

• \(Q_t\) зависит только от \(Q_{t-1}\): \[P(Q_t | Q_{t-1}, Y_{t-1}, \ldots, Q_{1}, Y_{1}) =P(Q_t | Q_{t-1})\]
• \(Y_t\) зависит только от \(Q_t\): \[P(Y_t | Q_{t-1}, Y_{t-1}, \ldots, Q_{1}, Y_{1}) =P(Y_t | Q_t)\]

Пусть \(Q_t\) – дискретная переменная, принимающая одно из \(N\) значений \(q_t \in \{1,\ldots,N\}\). В предположении, что модель Маркова однородна по времени, т.е. \(P(Q_t | Q_{t-1})\) не зависит от времени, мы можем задать матрицу переходов \[A = \{a_{ij}\} = p(Q_t=j|Q_{t-1}=i).\]

Вероятности начальных состояний определяются начальным распределением \(π_i = p(Q_1 = i).\)

Наблюдения \(Y_t\) может иметь одно из \(L\) возможных значений \(y_t \in \{o_1,\ldots,o_L\}\). Вероятность данного наблюдения тогда определяется как \[b_{ij} = b_j(o_i) = p(Y_t = o_i | Q_t = j),\] и тогда можно определить матрицу \(B = \{b_{ij}\}\) размера \(L\times N\).

Тогда модель описывается тройкой \((A, B, π)\), которую мы обозначим \(λ\).

Алгоритм Баума-Велша находит \(λ^* = \arg \max_λ P(y|λ)\), параметры, максимизирующие вероятность наблюдений \(y\).

Алгоритм

В качестве исходных выбираются случайные значения \(λ=(A,B,π)\). Алгоритм итеративно обновляет \(λ\) до схождения к неподвижной точке.

Прямая процедура

Обозначим вероятность появления заданной последовательности измерений \(\{y_1, \ldots, y_t\}\) как \(α_i(t)\). \[α_i(t) = p(Y_1=y_1, \ldots, Y_t=y_t, Q_t = i | λ)\]

Эту величину можно вычислить рекурсивно:

\[α_i(1) = π_i \cdot b_i(y_1)\] \[α_i(t+1) = b_i(y_{t+1})\sum_{j=1}^N α_j(t)\cdot a_{ji}\]

Обратная процедура

Обозначим вероятность “будущей” последовательности измерений, при условии что в “текущий” момент \(t\) скрытое состояние – \(i\): \[β_i(t) = p(Y_{t+1}=y_{t+1}, \ldots, Y_T=y_T | Q_t=i, λ)\]

Можно вычислить \(β_i(t)\) как

\[β_i(T) = p(Y_T=y_T | Q_t=i, λ) = 1\] \[β_i(t) = \sum_{j=1}^N β_j(t+1)a_{ij}b_j(y_{t+1})\]

Теперь, зная \(α\) и \(β\), можно вычислить

\[γ_i(t) = p(Q_t=i|y, λ) = \frac{α_i(t)β_i(t)}{\sum_{j=1}^N α_j(t)β_j(t)}\] \[ξ_{ij}(t) = p(Q_t=i, Q_{t+1}=j|y, λ) = \frac{α_i(t)a_{ij}β_j(t+1)b_j(y_{t+1})}{\sum_{k=1}^N\sum_{w=1}^N α_k(t)a_{kw}β_w(t+1)b_w(y_{t+1})}\]

Из этого вычисляются новые параметры модели:

\[π_i^* = γ_i(1)\]

\[a_{ij}^* = \frac{\sum_{t=1}^{T-1}ξ_{ij}(t)} {\sum_{t=1}^{T-1} γ_i(t)}\]

\[b_i^*(o_k) = \frac{\sum_{t=1}^{T}δ_{y_t=o_k}γ_i(t)} {\sum_{t=1}^{T} γ_i(t)}\]