CADE — интересный способ поиска аномалий в многомерных данных

Введение

Один из способов искать аномалии в наборе данных — использовать соответствующую данным плотность случайной величины как оценку «аномальности» точек. Интуиция проста — чем меньше значение плотности распределения в точке, тем более эта точка «аномальна». И наоборот, области с высокими значениями плотности распределения соответствуют наборам более «нормальных» точек. То есть зная плотность распределения, мы, в одной из возможных постановок задачи, можем находить и аномалии. И именно в обладании (а реальной жизни в приближении) плотности вероятности и возникает трудность, особенно когда мы работаем с данными большой размерности, и стандартные способы приближения плотности становятся вычислительно дорогими. На помощь приходят более «хитрые» способы приблизить плотность вероятности. Существует целый класс таких методов, и один из них произвёл на меня особое впечатление оригинальностью своей мысли. Он показался мне достаточно интересным для того, чтобы более популярно изложить его на Хабре. Этот метод называется CADE (Classifier Adjusted Density Estimation), и впервые он был изложен еще в 2001 году в статье [1], хотя его реальный потенциал был исследован позже.

CADE — это метод приближения плотности распределения, который хорошо справляется с большими размерностями и неинформативными признаками. Что важно, для корректной работы он требует предположения о независимости признаков. Однако по утверждениям этой статьи [2], метод ранжирует точки по аномальности на том же уровне, а при некоторых условиях даже лучше, чем LOC (Local Outlier Factor) — один из популярных density‑based методов для поиска аномалий. Суть работы CADE заключается в использовании любого любимого вами алгоритма классификации, обученного различать известные данные от искусственно сгенерированных точек. И как это нам поможет приблизить плотность распределения? Итак, давайте рассмотрим.

Как работает CADE

Если вкратце, то CADE предлагает следующее:

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

Как это работает? Давайте разбираться.

Я буду использовать нотации из вышеупомянутой статьи [2], на которую я и опирался, чтобы вам было проще её перечитать, если у вас возникнет такое желание.

Пусть T — это имеющиеся данные (с небольшой долей аномалий). P(X|T) — их плотность распределения, которую мы хотим приблизить. Первое, что нам надо сделать — это выбрать то самое «произвольное» распределение, из которого мы можем сгенерировать данные. Обозначим его как A, а соответствующую плотность распределения — P(X|A). Авторы статьи утверждают, что в качестве такого распределения хорошо подходит равномерное или нормальное распределения, хотя ограничений на выбор распределения нет. Желательно только, чтобы оно полностью покрывало данные, то есть была ненулевая вероятность выбрать точку, близкую к изначальным данным. Что касается размера выборки, то авторы использовали размер, равный размеру выборки T, однако я поделюсь своими соображениями на этот счёт позже. Точки из A иногда называют искусственными аномалиями.

Следующий шаг — выбрать классификатор. В качестве него мы можем использовать всё, что способно рассчитать вероятность принадлежность точки к изначальному классу T, а именно P(C = T|X). Классификатором может быть Random Forest Classifier, KNN, Decision Tree Classifier, подходящая по архитектуре нейронная сеть и так далее.

Естественно, для каждой точки сумма вероятностей быть в классе T и классе A равна 1. Другими словами:

Теперь приступим к выводу формулы, которая позволит нам использовать предсказание классификатора для приближения изначальной плотности распределения P(X|T). Для этого мы применим теорему Байеса к P(C=T|X):

Далее, решив это уравнение относительно P(X|T), получим следующую формулу:

Это и есть приближение искомой плотности распределения по CADE! Оно состоит из трёх основных множителей.

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

Известная плотность распределения у сгенерированных данных:

Отношение, вычисляемое по предсказаниям классификатора:

Если P(X|A) — константа, то есть A было выбрано из равномерного распределения, то ранжирование точек по значениям плотности распределения будет совпадать с ранжированием по значениям классификатора (умножение на положительную константу не повлияет на порядок).

Если же мы угадали и выбрали P(X|A) близким к P(X|T), то классификатор не сможет различить распределения и для всех точек предскажет значение 0.5. Тогда ранжирование логичным образом будет совпадать с ранжированием по P(X|A).

Покончим с формулами! Давайте посмотрим на небольшую анимацию о процессе работы CADE без вычислений и финальной формулы.

Как я упомянул выше, исследователи из [2] показали, что метод действительно хорошо работает с данными высокой размерности и при наличии «шумных» атрибутов. А также использование множителя P(X|A), что является отличительной особенностью CADE, улучшает работу метода при наличии коррелирующих признаков. Если у вас будет желание, я рекомендую ознакомиться с результатами из этого исследования. В нем вы также увидите значения метрики ROC AUC для разных наборов данных и классификаторов.

Работаем с CADE на практике

На момент написания статьи, я не обнаружил реализации CADE на Python (не путать с другим CADE — Compass‑aligned Distributional Embeddings, реализация для которого есть в PyPI). Поэтому в этой статье, а также на Github (cade‑outliers), я оставлю пример кода, который реализует поиск аномалий с помощью этого метода.

Предложения по улучшению кода приветствуются!

Класс CADEOutliers отвечает за то, чтобы приближать плотность распределения исходных данных и возвращать точки, ранжированные по значениям плотности.

import numpy as np


class CADEOutliers():
    """Класс для обнаружения выбросов с использованием CADE (Classifier Adjusted Density Estimation)"""
    
    _supported_A_dist = ["uniform"]
    
    def __init__(self, classifier, A_dist, A_size):
        """Создать объект обнаружения выбросов CADE.
        
        Параметры
        ----------
        classifier : model
            Модель классификатора sklearn, поддерживающая методы fit и predict_proba.
        
        A_dist : string
            Распределение для генерации искусственных аномалий.
        
        A_size : float | int
            Если float, то количество искусственных аномалий равно |X| * A_size, где X - анализируемые данные. 
            Если int, то количество искусственных аномалий равно A_size.
        """
        # Установка классификатора
        self.classifier = classifier
        
        # Установка типа распределения
        if A_dist not in self._supported_A_dist:
            raise ValueError(f"{A_dist} распределение пока не поддерживается для генерации искусственных аномалий. Выберите из следующего списка: {self._supported_A_dist}")
        self.A_dist = A_dist
        
        # Установка количества искусственных аномалий для генерации
        if isinstance(A_size, int) and A_size >= 1:
            self.sample_size = A_size
        elif isinstance(A_size, float) and A_size >= 0:
            self.sample_size = None
        else:
            raise ValueError("A_size должен быть либо int[1, inf) либо float[0, inf)")
        self.A_size = A_size
        
    def _generate_A(self, dist, attrs):
        """Генерация искусственных аномалий.
        
        Параметры
        ----------
        dist : string
            Распределение для генерации искусственных аномалий.
        attrs : dict
            Параметры для генерации искусственных аномалий. 
            Если dist равен "uniform", attrs = {'low': array, 'high': array, 'size': (sample_size, n_dim)} с минимальными и максимальными порогами для каждого измерения.
        """
        if dist == 'uniform':
            generator = np.random.uniform
        else:
            raise ValueError(f"{dist} распределение пока не поддерживается для генерации искусственных аномалий. Выберите из следующего списка: {self._supported_A_dist}")
        
        A = generator(**attrs)
        
        return A
    
    def fit(self, X):
        """Оценка распределения X.
        
        Параметры
        ----------
        X : np.ndarray
            Массив с образцами
        """
        if self.sample_size is None:
            sample_size = int(len(X) * self.A_size) + 1
        else:
            sample_size = int(self.sample_size)
        
        if self.A_dist == 'uniform':
            n_dim = X.shape[1]
            attrs = {
                'low': X.min(axis=0),
                'high': X.max(axis=0),
                'size': (sample_size, n_dim)
            }
            # вычисление плотности вероятности для равномерного распределения
            P_A = 1
            for s in X.max(axis=0) - X.min(axis=0):
                if s == 0:
                    s = 1
                P_A /= s
            self.P_A = lambda x: np.array([P_A] * x.shape[0])
        else:
            attrs = None
            
        A = self._generate_A(self.A_dist, attrs)
        
        combined_data = np.vstack([X, A])
        target = np.hstack([np.zeros(X.shape[0]), np.ones(A.shape[0])])
        
        self.classifier.fit(combined_data, target)
        self.A_X_ratio = A.shape[0] / X.shape[0]
        
        
    def get_ranking(self, X):
        """Возвращает ранжирование для каждого образца в X.
        Более высокие значения указывают на большую уверенность в том, что точка является выбросом.
        
        Параметры
        ----------
        X : np.ndarray
            Массив с образцами
        """
        predictions = self.classifier.predict_proba(X)[:, 0]
        
        anomaly_score = self.P_A(X) * self.A_X_ratio * (predictions / (1 - predictions + 1e-5))
        
        return anomaly_score

Генерирование данных и применение к ним CADE:

import numpy as np
import seaborn as sns
import matplotlib.pyplot as plt

from sklearn.metrics import roc_auc_score
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import train_test_split

from my_module import CADEOutliers  # Убедитесь, что CADEOutliers импортируется из правильного модуля


# Генерация данных
# X состоит из смеси 3 нормальных распределений
# 3-е распределение - это небольшая часть аномалий
X = np.vstack([2 * np.random.randn(5000, 2), 7 + np.random.randn(3000, 2), 25 + np.random.randn(100, 2)])
y = np.hstack([np.zeros(8000), np.ones(100)])

print('Размерность X: {}'.format(X.shape))

# Разделение данных на обучающую и тестовую части
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.33)

# Создание объекта CADE для обнаружения выбросов с использованием равномерного распределения искусственных аномалий размером 100% от заданного набора данных
cade = CADEOutliers(
    classifier=RandomForestClassifier(max_depth=3), 
    A_dist='uniform', 
    A_size=1.
)

# Обучение оценки плотности вероятности X с использованием X_train
cade.fit(X_train)

# Получение ранжирования для X_test. Меньшие значения указывают на более высокое сходство с аномалиями
ranking_by_dens = cade.get_ranking(X_test)

print('ROC AUC ранжирования аномалий с использованием CADE на тестовых данных: {}'.format(roc_auc_score(y_test, -ranking_by_dens)))

fig, ax = plt.subplots(1, 1, figsize=(7, 7))
sns.scatterplot(x=X_test[:, 0], y=X_test[:, 1], hue=ranking_by_dens, ax=ax, palette='magma')
plt.title('Точки с ранжированием, полученным с помощью CADE (низкие значения указывают на аномалии)')
plt.show()

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

Переобученный классификатор более склонен считать любые увиденные в процессе обучения точки «нормальными» (даже те, что кажутся аномалиями), а любые области без исходных данных при обучении — аномальными (даже те, что лежат рядом с увиденными точками). Другими словами, переобучение классификатора ведет к более детальному описанию плотности вероятности и большей чувствительностью к обнаружению аномалий.
Бóльшая плотность заполнения пространства искусственными аномалиями также ведет к более детальному описанию плотности распределения, т. е. к тому же эффекту, что и переобученный классификатор.
Из статьи [2] — в качестве классификатора неплохо подходит Random Forest Classifier и KNN, а в качестве распределения для искусственных аномалий — равномерное.

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

Источники

[1] T. Hastie, R. Tibshirani, and J. H. Friedman. The Elements of Statistical Learning: Data Mining, Inference, and Prediction. New York: Springer-Verlag, 2001.

[2] https://seppe.net/aa/papers/CADE.pdf

PreviousФункции потерь в Python — простая реализация NextКак использовать Python для проверки нормальности

Last updated 7 months ago

import numpy as np class CADEOutliers(): """Класс для обнаружения выбросов с использованием CADE (Classifier Adjusted Density Estimation)""" _supported_A_dist = ["uniform"] def __init__(self, classifier, A_dist, A_size): """Создать объект обнаружения выбросов CADE. Параметры ---------- classifier : model Модель классификатора sklearn, поддерживающая методы fit и predict_proba. A_dist : string Распределение для генерации искусственных аномалий. A_size : float | int Если float, то количество искусственных аномалий равно |X| * A_size, где X - анализируемые данные. Если int, то количество искусственных аномалий равно A_size. """ # Установка классификатора self.classifier = classifier # Установка типа распределения if A_dist not in self._supported_A_dist: raise ValueError(f"{A_dist} распределение пока не поддерживается для генерации искусственных аномалий. Выберите из следующего списка: {self._supported_A_dist}") self.A_dist = A_dist # Установка количества искусственных аномалий для генерации if isinstance(A_size, int) and A_size >= 1: self.sample_size = A_size elif isinstance(A_size, float) and A_size >= 0: self.sample_size = None else: raise ValueError("A_size должен быть либо int[1, inf) либо float[0, inf)") self.A_size = A_size def _generate_A(self, dist, attrs): """Генерация искусственных аномалий. Параметры ---------- dist : string Распределение для генерации искусственных аномалий. attrs : dict Параметры для генерации искусственных аномалий. Если dist равен "uniform", attrs = {'low': array, 'high': array, 'size': (sample_size, n_dim)} с минимальными и максимальными порогами для каждого измерения. """ if dist == 'uniform': generator = np.random.uniform else: raise ValueError(f"{dist} распределение пока не поддерживается для генерации искусственных аномалий. Выберите из следующего списка: {self._supported_A_dist}") A = generator(**attrs) return A def fit(self, X): """Оценка распределения X. Параметры ---------- X : np.ndarray Массив с образцами """ if self.sample_size is None: sample_size = int(len(X) * self.A_size) + 1 else: sample_size = int(self.sample_size) if self.A_dist == 'uniform': n_dim = X.shape[1] attrs = { 'low': X.min(axis=0), 'high': X.max(axis=0), 'size': (sample_size, n_dim) } # вычисление плотности вероятности для равномерного распределения P_A = 1 for s in X.max(axis=0) - X.min(axis=0): if s == 0: s = 1 P_A /= s self.P_A = lambda x: np.array([P_A] * x.shape[0]) else: attrs = None A = self._generate_A(self.A_dist, attrs) combined_data = np.vstack([X, A]) target = np.hstack([np.zeros(X.shape[0]), np.ones(A.shape[0])]) self.classifier.fit(combined_data, target) self.A_X_ratio = A.shape[0] / X.shape[0] def get_ranking(self, X): """Возвращает ранжирование для каждого образца в X. Более высокие значения указывают на большую уверенность в том, что точка является выбросом. Параметры ---------- X : np.ndarray Массив с образцами """ predictions = self.classifier.predict_proba(X)[:, 0] anomaly_score = self.P_A(X) * self.A_X_ratio * (predictions / (1 - predictions + 1e-5)) return anomaly_score

import numpy as np import seaborn as sns import matplotlib.pyplot as plt from sklearn.metrics import roc_auc_score from sklearn.ensemble import RandomForestClassifier from sklearn.model_selection import train_test_split from my_module import CADEOutliers # Убедитесь, что CADEOutliers импортируется из правильного модуля # Генерация данных # X состоит из смеси 3 нормальных распределений # 3-е распределение - это небольшая часть аномалий X = np.vstack([2 * np.random.randn(5000, 2), 7 + np.random.randn(3000, 2), 25 + np.random.randn(100, 2)]) y = np.hstack([np.zeros(8000), np.ones(100)]) print('Размерность X: {}'.format(X.shape)) # Разделение данных на обучающую и тестовую части X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.33) # Создание объекта CADE для обнаружения выбросов с использованием равномерного распределения искусственных аномалий размером 100% от заданного набора данных cade = CADEOutliers( classifier=RandomForestClassifier(max_depth=3), A_dist='uniform', A_size=1. ) # Обучение оценки плотности вероятности X с использованием X_train cade.fit(X_train) # Получение ранжирования для X_test. Меньшие значения указывают на более высокое сходство с аномалиями ranking_by_dens = cade.get_ranking(X_test) print('ROC AUC ранжирования аномалий с использованием CADE на тестовых данных: {}'.format(roc_auc_score(y_test, -ranking_by_dens))) fig, ax = plt.subplots(1, 1, figsize=(7, 7)) sns.scatterplot(x=X_test[:, 0], y=X_test[:, 1], hue=ranking_by_dens, ax=ax, palette='magma') plt.title('Точки с ранжированием, полученным с помощью CADE (низкие значения указывают на аномалии)') plt.show()