Введение

Один из способов искать аномалии в наборе данных — использовать соответствующую данным плотность случайной величины как оценку «аномальности» точек. Интуиция проста — чем меньше значение плотности распределения в точке, тем более эта точка «аномальна». И наоборот, области с высокими значениями плотности распределения соответствуют наборам более «нормальных» точек. То есть зная плотность распределения, мы, в одной из возможных постановок задачи, можем находить и аномалии. И именно в обладании (а реальной жизни в приближении) плотности вероятности и возникает трудность, особенно когда мы работаем с данными большой размерности, и стандартные способы приближения плотности становятся вычислительно дорогими. На помощь приходят более «хитрые» способы приблизить плотность вероятности. Существует целый класс таких методов, и один из них произвёл на меня особое впечатление оригинальностью своей мысли. Он показался мне достаточно интересным для того, чтобы более популярно изложить его на Хабре. Этот метод называется 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(C=T|X)=\frac{P(X|T)P(C=T)}{P(X)}=\frac{P(X|T)P(C=T)}{P(X|T)P(C=T)+P(X|A)P(C=A)}$

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

$P(X|T)=\frac{P(C=A)}{P(C=T)}P(X|A)\frac{P(C=T|X)}{1-P(C=T|X)}$

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

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

$\frac{P(C=A)}{P(C=T)}$

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

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

$\frac{P(C=T|X)}{1-P(C=T|X)}$

Если 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():
    """Class for outliers detection with CADE (Classifier Adjusted Density Estimation)"""
    
    _supported_A_dist = ["uniform"]
    
    def __init__(self, classifier, A_dist, A_size):
        """Create CADE outlier detection object.
        
        Parameters
        ----------
        classifier : model
            sklearn classifier model supporting fit and predict_proba methods.
        
        A_dist : string
            Distribution for generating artificial anomalies.
        
        A_size : float | int
            If float, then the number of artifical anomalies is |X| * A_size, where X - analyzed data. 
            If int, then the number of artifical anomalies is just A_size.
        """
        # Set classifier
        self.classifier = classifier
        
        # Set the distribution kind
        if A_dist not in self._supported_A_dist:
            raise ValueError(f"{A_dist} distribution is not supported yet for generating artifical anomalies. Choose from the following list: {self._supported_A_dist}")
        self.A_dist = A_dist
        
        # Set the number of artifical anomalies to generate
        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 must be either int[1, inf) of float[0, inf)")
        self.A_size = A_size
        
    def _generate_A(self, dist, attrs):
        """Generate artifical anomalies.
        
        Parameters
        ----------
        dist : string
            Distribution for generating artificial anomalies.
        attrs : dict
            Key-value arguments for generating artifical anomalies. 
            If dist is "uniform", attrs = {'low': array, 'high': array, 'size': (sample_size, n_dim)} with min and high thresholds for each dimension.
        """
        if dist == 'uniform':
            generator = np.random.uniform
        else:
            raise ValueError(f"{dist} distribution is not supported yet for generating artifical anomalies. Choose from the following list: {self._supported_A_dist}")
        
        A = generator(**attrs)
        
        return A
    
    def fit(self, X):
        """Fit the distribution of X.
        
        Parameters
        ----------
        X : np.ndarray
            Array with samples
        """
        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)
            }
            # calculate uniform probability density
            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):
        """Return rankings for each sample in X.
        Higher values indicate higher certainty that the point is an outlier.
        
        Parameters
        ----------
        X : np.ndarray
            Array with samples
        """
        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

import CADEOutliers


# Generate data
# X consists of a mix of 3 normal distributions
# The 3rd is a small portion of anomalies
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('The shape of X: {}'.format(X.shape))

# Split the data into train and test part
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.33)

# Create CADE outliers object with uniform distribution of artifical anomalies and with size of 100% of the given dataset
cade = CADEOutliers(
    classifier=RandomForestClassifier(max_depth=3), 
    A_dist='uniform', 
    A_size=1.
)

# Fit probability density of X using X_train
cade.fit(X_train)

# Get rankings for X_test. Lower values indicate higher similarity to outliers
ranking_by_dens = cade.get_ranking(X_test)

print('ROC AUC of ranking anomalies by CADE on test data: {}'.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('Points with rankings derived with CADE (low values indicate anomalies)')
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

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