异常检测系列：主成分分析PCA异常值分数检测

数智笔记

于 2024-01-25 13:39:28 发布

阅读量1.9k

点赞数 19

分类专栏：数据挖掘数据挖掘案例文章标签：机器学习数据挖掘算法人工智能

本文链接：https://blog.csdn.net/wjjc1017/article/details/135842327

版权

数据挖掘同时被 2 个专栏收录

166 篇文章 33 订阅

订阅专栏

数据挖掘案例

20 篇文章 2 订阅

订阅专栏

主成分分析（PCA）是数据分析中的一种流行技术之一。它执行线性转换，将可能相关的变量的一组数据点转换为线性不相关的变量的一组数据点。

在本章中，我将解释为什么降维技术可以找到异常值。我将解释PCA如何应用于异常检测并生成异常值分数。然后，我将演示如何在PyOD中使用PCA异常值分数来检测异常值。

(A) PCA是如何工作的？

高维数据集是具有大量变量的数据集。这样的数据集带来了“维度灾难”，并经常出现计算挑战。尽管高性能计算可以在某种程度上处理高维数据，但在许多应用中，仍然需要降低原始数据的维度。PCA降低了由大量相关变量组成的数据集的维度，同时尽可能保留方差。PCA找到了原始变量的新变量，这些新变量只是它们的线性组合。这些新变量被称为主成分（PCs）。这些主成分是正交的，意味着主成分彼此垂直。在三维情况下，X不能由Y表示，Y也不能由Z表示。正交一词来自希腊语orthogōnios，意思是“直角的”。

线性转换可以在图（A）中显示。原始的X和Y轴被“旋转”以更好地与数据体（以红色显示）对齐。第一个主成分（PC1）将捕获数据中的最大方差。第二个主成分捕获了PC1未捕获的最大方差。第三、第四、…和其他成分将继续捕获之前主成分未捕获的内容，直到解释了所有方差。主成分的数量应等于原始变量的数量。

图（A）：PCA（作者提供的图像）

在线性变换中，数据的协方差矩阵被分解为与特征值相关联的正交向量，称为特征向量。特征值是特征向量缩放的因子。具有较高特征值的特征向量捕捉了数据中大部分的方差。“Eigen”是德语单词，意思是“自己”或“典型”。特征向量也可以称为“特征向量”。

(B) 降维可以找出异常值

当维度降低到几个主要维度时，可以识别出模式，然后揭示出异常值。这种直觉来自于异常值本身。D.M.Hawkins [1]曾经对异常值给出了一个很好的定义：“异常值是与其他观察值相差很大，以至于引起怀疑它是由不同机制生成的观察值。”因为异常值往往遵循不同的规律，它们通常不在前几个主成分中。当它们投影到低维超平面时，它们将落在一些具有较小特征值的唯一特征向量上。我们可以说异常值检测是降维的副产品。具有这个属性的数据点在PCA中的异常值得分如下：

异常值得分 = 每个观察值到所选特征向量构成的超平面的加权欧氏距离之和

© 在运行PCA之前记得对数据进行标准化

在执行PCA之前，对数据进行标准化是必要的。如果对数据进行标准化，所有变量都具有相同的标准差，因此所有变量具有相同的权重。如果忘记对数据进行标准化，具有较高标准差的变量在轴的计算中将具有较高的权重，而具有较低标准差的变量则具有较低的权重。另一种思考标准化需求的方式是，数据集中的不同变量可能具有不同的测量单位。例如，一个变量可能是美元金额，另一个变量可能是单位等。需要对所有变量进行数据标准化。PyOD中的PCA类具有内置的标准化过程。它在执行PCA之前对数据进行标准化。

(D) 建模过程

我遵循第1章提出的建模过程。这个过程让您专注于模型的开发、评估和结果的解释。

模型开发
阈值确定
正常组和异常组的描述统计信息

(D.1) 第一步 — 构建模型

我生成了一个包含500个观测值和六个变量的模拟数据集。我将异常值的百分比设置为5%，即“contamination=0.05”。无监督模型只使用X变量。模拟数据集中的目标变量Y仅用于验证。


import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from pyod.utils.data import generate_data
contamination = 0.05 # percentage of outliers
n_train = 500       # number of training points
n_test = 500        # number of testing points
n_features = 6      # number of features
X_train, X_test, y_train, y_test = generate_data(
    n_train=n_train, 
    n_test=n_test, 
    n_features= n_features, 
    contamination=contamination, 
    random_state=123)

X_train_pd = pd.DataFrame(X_train)
X_train_pd.head()

以下是一些观察结果：

图表 (D.1) 绘制了散点图中的前两个变量。黄色点是异常值，紫色点是正常数据点。

# Plot
plt.scatter(X_train_pd[0], X_train_pd[1], c=y_train, alpha=0.8)
plt.title('Scatter plot')
plt.xlabel('x0')
plt.ylabel('x1')
plt.show()

图(D.1)

以下代码构建模型并对训练和测试数据进行评分。让我们逐行查看：

label_：这是训练数据的标签向量。如果您使用.predict()训练数据，则相同。
decision_scores_：这是训练数据的分数向量。如果您使用.decision_functions()训练数据，则相同。
decision_score()：此评分函数为每个观察分配异常值分数。
predict()：这是分配1或0的预测函数，基于分配的阈值。
contamination：这是异常值的百分比。在大多数情况下，我们不知道异常值的百分比，因此可以根据任何先前知识分配一个值。 PyOD默认将污染率设置为10％。此参数不影响异常值分数的计算。


from pyod.models.pca import PCA
pca = PCA(contamination=0.05) 
pca.fit(X_train)

# get the prediction labels and outlier scores of the training data
y_train_pred = pca.labels_  # binary labels (0: inliers, 1: outliers)
y_train_scores = pca.decision_scores_  # .decision_scores_ yields the raw outlier scores for the training data
y_train_scores = pca.decision_function(X_train) # You also can use .decision_function()
y_train_pred = pca.predict(X_train) # You also can use .decision_function()

y_test_scores = pca.decision_function(X_test) # You also can use .decision_function()
y_test_pred = pca.predict(X_test) # You also can use .decision_function()

让我们看看PCA的默认参数。标准化的默认参数为true，这意味着在执行PCA之前执行归一化。污染率设置为5％。


pca.get_params()

在PCA中一个非常重要的指标是解释方差比。它解释了每个主成分解释了多少总方差。PyOD通过其函数.explained_variance_ratio_传递解释方差。下面的输出显示第一个PC解释了73.4%的方差，第二个PC解释了额外的7.4%，第三个解释了5.6%，以此类推。六个主成分的解释方差之和应该为100%。


[pca.explained_variance_,
pca.explained_variance_ratio_]

(D.2) 第二步 — 确定一个合理的阈值

PyOD有一个内置函数threshold_，它根据污染率计算训练数据的阈值。如果我们不指定污染率，默认值为10%。因为我们在模型中指定了5%，所以阈值是在5%处的值。


print("The threshold for the defined comtanimation rate:" , pca.threshold_)

在大多数情况下，我们不知道异常值的百分比。我们可以使用异常值得分的直方图来选择一个合理的阈值。如果有任何先前的知识表明异常值为1％，则应选择一个导致大约1％异常值的阈值。图（D.2）中PCA异常值得分的直方图表明阈值为200.0，因为直方图中存在一个自然的切割点。


# get the prediction on the test data
y_test_pred = pca.predict(X_test)  # outlier labels (0 or 1)
y_test_scores = pca.decision_function(X_test)  # outlier scores

import matplotlib.pyplot as plt
plt.hist(y_train_scores, bins='auto')  # arguments are passed to np.histogram
plt.title("Histogram with 'auto' bins")
plt.xlabel('PCA outlier score')
plt.show()

图(D.2): PCA异常值得分的直方图

(D.3) 第三步 — 展示正常组和异常组的描述统计

对正常组和异常组进行分析是展示模型的可靠性的关键步骤。正常组和异常组的特征统计应与领域知识一致。如果异常组中某个特征的均值应该很高，但结果却相反，建议您检查、修改或丢弃该特征。您应该迭代建模过程，直到所有特征与先前的知识一致。另一方面，建议您验证数据的先前知识是否提供了新的见解。


threshold = pca.threshold_
threshold

def descriptive_stat_threshold(df,pred_score, threshold):
    # Let's see how many '0's and '1's.
    df = pd.DataFrame(df)
    df['Anomaly_Score'] = pred_score
    df['Group'] = np.where(df['Anomaly_Score']< threshold, 'Normal', 'Outlier')

    # Now let's show the summary statistics:
    cnt = df.groupby('Group')['Anomaly_Score'].count().reset_index().rename(columns={'Anomaly_Score':'Count'})
    cnt['Count %'] = (cnt['Count'] / cnt['Count'].sum()) * 100 # The count and count %
    stat = df.groupby('Group').mean().round(2).reset_index() # The avg.
    stat = cnt.merge(stat, left_on='Group',right_on='Group') # Put the count and the avg. together
    return (stat)

descriptive_stat_threshold(X_train,y_train_scores, threshold)

表格（D.3）：正常组和异常组的描述统计

上述表格展示了正常组和异常组的特征。它显示了正常组和异常组的计数和计数百分比。请记住，为了有效展示，您需要使用特征名称标记特征。该表格告诉我们几个重要的结果：

异常组的大小： 一旦确定了阈值，大小就确定了。如果阈值是从图（D.2）中得出的，并且没有先验知识，大小统计量成为一个很好的起点参考。
每个组中的特征统计： 所有均值必须与领域知识一致。在我们的案例中，异常组的均值小于正常组的均值。
平均异常分数： 异常组的平均分数应该高于正常组的平均分数（844.33 > 124.59）。您不需要对分数进行过多解释。

因为我们在数据生成中有真实值，所以我们可以生成一个混淆矩阵来了解模型的性能。该模型表现良好，并且识别出了所有的25个异常值。


Actual_pred = pd.DataFrame({'Actual': y_test, 'Anomaly_Score': y_test_scores})
Actual_pred['Pred'] = np.where(Actual_pred['Anomaly_Score']< threshold,0,1)
pd.crosstab(Actual_pred['Actual'],Actual_pred['Pred'])

(E) 模型比较

从第2章到现在，我们已经建立了几个模型。看到这些模型的重叠部分是很有趣的。让我们检查一下是否有多个模型都识别出了异常值。

########
# HBOS #
########
from pyod.models.hbos import HBOS
n_bins = 50
hbos = HBOS(n_bins=n_bins,contamination=0.05) 
hbos.fit(X_train)
y_train_hbos_scores = hbos.decision_function(X_train)
y_test_hbos_scores = hbos.decision_function(X_test)

########
# ECOD #
########
from pyod.models.ecod import ECOD
clf_name = 'ECOD'
ecod = ECOD(contamination=0.05) 
ecod.fit(X_train)
y_test_ecod_scores = ecod.decision_function(X_test)

########
# PCA  #
########
from pyod.models.pca import PCA
pca = PCA(contamination=0.05) 
pca.fit(X_train)

# get the prediction labels and outlier scores of the training data
y_test_pca_scores = pca.decision_function(X_test) # You also can use .decision_function()

# Thresholds
[ecod.threshold_, hbos.threshold_, pca.threshold_]

我将实际的Y值、HBO分数、ECOD分数和PCA分数放在一个数据框中。


# Put the actual, the HBO score and the ECOD score together
Actual_preds = pd.DataFrame({'Actual': y_test, 
                             'HBO_Score': y_test_hbos_scores, 
                             'ECOD_Score': y_test_ecod_scores, 
                             'PCA_Score': y_test_pca_scores})
Actual_preds['HBOS_pred'] = np.where(Actual_preds['HBO_Score']>hbos.threshold_,1,0)
Actual_preds['ECOD_pred'] = np.where(Actual_preds['ECOD_Score']>ecod.threshold_,1,0)
Actual_preds['PCA_pred'] = np.where(Actual_preds['PCA_Score']>pca.threshold_,1,0)
Actual_preds.head()

当我对HBOS和PCA的预测进行交叉表格分析时，有25个观测值被两个模型都识别为异常值。您可以交叉分析其他分数。如果没有基准数据，这种交叉检验尤其有帮助。

pd.crosstab(Actual_preds['HBOS_pred'],Actual_preds['PCA_pred'])

(F) PCA算法概述

异常值与正常数据点不同。当它们投影到低维超平面时，它们将落在具有较小特征值的特征向量上。
数据点的异常值分数是数据点到低维超平面的加权欧氏距离的总和。

参考文献

Hawkins D. M. (1980), Identification of Outliers, Chapman & Hall, London.
pyod.models.pca — pyod 1.0.5 documentation: https://pyod.readthedocs.io/en/latest/_modules/pyod/models/pca.html#PCA