异常检测——基于统计学的方法
主要内容:高斯分布、箱线图、HBOS
1、概述
1.1、基于统计学方法的一般思想:学习一个拟合给定数据集的生成模型,然后识别该模型低概率区域中的对象,把他们作为异常点
1.2、基于统计学方法的主要类型:参数方法和非参数方法
-
参数方法:
-
非参数方法:
2、参数方法
2.1、基于正态分布的一元异常点检测
- 原始数据->正态分布参数->概率密度函数->阈值->异常点
# 构造数据
import numpy as np
data1 = np.random.randn(5000)*3+2.5
data2 = np.random.randn(5000)*4+3
# data2 = np.random.randint(1,size=(5000))
data = np.r_[data1,data2]
# 基于箱线图可视化数据
'''中位数、上四分之一位数、下四分之一位数'''
import matplotlib.pyplot as plt
import seaborn as sns
sns.boxplot(data=data)
# 求得所构造数据的均值和标准差
data_u = data.mean()
data_sigma = data.std()
# 基于ks检验查看数据是否符合正态分布
from scipy import stats
"""
kstest方法:KS检验,参数分别是:待检验的数据,检验方法(这里设置成norm正态分布),均值与标准差
结果返回两个值:statistic → D值,pvalue → P值
p值大于0.05,为正态分布
H0:样本符合
H1:样本不符合
如何p>0.05接受H0 ,反之
"""
ks_result = stats.kstest(data,'norm',(data_u,data_sigma))
if ks_result[1]>0.05: # pvalue
print('构造的数据data符合正态分布')
else:
print('构造的数据data不符合正态分布')
# 基于阈值通过概率密度函数找到异常点
import math
def norm_density(x,u,std):
norm_den_res = 1/(std*pow(2*math.pi,0.5))*np.exp(-(x-u)**2/(2*std**2))
return norm_den_res
data_abnormal=[]
for data_value in data:
if norm_density(data_value,data_u,data_sigma)<10**(-4):
data_abnormal.append(data_value)
print('数据集中的异常值为:'+ str(data_abnormal))
构造的数据data符合正态分布
数据集中的异常值为:[-11.46365186764682, 16.416823416708013, 16.327736773836584, 16.308623544015518, 17.745661790920828, -10.685020690566537, 16.26914075409187]
2.2、多元异常点检测
1、各个维度的特征之间相互独立:一元高斯分布
- 核心思想:将多元异常点检测任务转化为一元异常点检测任务
2、各个维度的特征之家有相关性:多元高斯分布
2.2.1、独立高斯分布的异常点检测
# 构造数据——二维随机变量
import numpy as np
data1 = np.random.randn(5000)*3+2.5
data2 = np.random.randn(5000)*4+3
# data2 = np.random.randint(1,size=(5000))
data = np.c_[data1,data2]
# 求得所构造数据的均值和标准差
data_u1 = data[:,0].mean()
data_u2 = data[:,1].mean()
data_sigma1 = data[:,0].std()
data_sigma2 = data[:,0].std()
# 基于阈值通过概率密度函数找到异常点
import math
def norm_density(x,u,std):
norm_den_res = 1/(std*pow(2*math.pi,0.5))*np.exp(-(x-u)**2/(2*std**2))
return norm_den_res
data_abnormal=[]
for i in range(data.shape[0]):
if norm_density(data[i,0],data_u1,data_sigma1)*norm_density(data[i,1],data_u2,data_sigma1)<10**(-6):
data_abnormal.append(data[i,:])
print('数据集中的异常值为:'+ str(data_abnormal))
数据集中的异常值为:[array([ 5.22234993, -10.72594557]), array([ 3.29865232, 16.81809963]), array([ 1.78894018, -10.44736964]), array([ 6.19805727, 15.6589662 ]), array([ 5.44369811, -10.3625621 ]), array([ -2.68441286, -10.14012666])]
对于多维随机变量,将所有维度上的概率密度函数相乘,然后在判断其值是否小于阈值epsilon。若大于或等于于阈值则该点是正常点;若小于阈值则该点是异常点。
2.2.2、多元高斯分布的异常点检测
# 构造数据——二维随机变量
# import numpy as np
# data1 = np.random.randn(5000)*3+2.5
# data2 = np.random.randn(5000)*4+3
# # data2 = np.random.randint(1,size=(5000))
# data = np.c_[data1,data2]
# 求得所构造数据的均值和标准差
data_u1 = data[:,0].mean()
data_u2 = data[:,1].mean()
data_u = [data_u1,data_u2]
for i in range(data.shape[0]):
data_sigma = data_sigma + np.mat(data[i,:]-data_u).T*np.mat(data[i,:]-data_u)
data_sigma = data_sigma/data.shape[0]
# 基于阈值通过概率密度函数找到异常点
import math
def norm_density(x,dimensional,u,std):
norm_den_res = 1/(pow(2*math.pi,dimensional/2)*pow(np.linalg.det(std),0.5))*np.exp(-0.5*np.mat(x-u)*np.linalg.inv(std)*np.mat(x-u).T)
return norm_den_res
data_abnormal=[]
for i in range(data.shape[0]):
if norm_density(data[i,:],data.shape[1],data_u,data_sigma)<10**(-6):
data_abnormal.append(data[i,:])
print('数据集中的异常值为:'+ str(data_abnormal))
数据集中的异常值为:[]
可以看到,当将独立高斯分布修改为多元高斯分布后,相同的数据集中不再产生异常值,这是因为多元高斯分布加入了变量之间的相关性信息
2.3、使用混合参数分布
在许多情况下假定数据是由正态分布产生的。当实际数据很复杂时,这种假定过于简单,可以假定数据是被混合参数分布产生的。
3、非参数方法
- 输入数据->学习“正常”模型(而非先验假设)->适用情况更多
使用直方图检测异常点
1、构造直方图:使用输入数据构造直方图——一元or多元
2、检测异常点:简单方法:(1)对照直方图检查是否落入直方图箱体(2)使用直方图赋予每个对象一个异常点得分
缺点:很难选择合适的箱的尺寸(太小:误识为异常点;太大:假扮为正常点)
# 构造数据——二维随机变量
# import numpy as np
# data1 = np.random.randn(5000)*3+2.5
# data2 = np.random.randn(5000)*4+3
# # data2 = np.random.randint(1,size=(5000))
# data = np.c_[data1,data2]
import matplotlib.pyplot as plt
plt.hist(data)
plt.show()
4、基于角度的方法
- 核心思想:数据边界上的数据很可能将整个数据包围在一个较小的角度内,而内部的数据点则可以不同的角度围绕着他们
# 导入库函数
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
%matplotlib inline
import matplotlib.font_manager
from pyod.models.abod import ABOD
from pyod.models.knn import KNN
# 创建随机数据集
from pyod.utils.data import generate_data, get_outliers_inliers
## 生成2features的随机数据
X_train, Y_train = generate_data(n_train=200,train_only=True, n_features=2)
## 默认情况下,Generate数据功能中的异常分数为0.1
outlier_fraction = 0.1
## 在不同的numpy阵列中储存outliers和inliers
x_outliers, x_inliers = get_outliers_inliers(X_train,Y_train)
n_inliers = len(x_inliers)
n_outliers = len(x_outliers)
## 分隔2features并使用它来绘制数据
F1 = X_train[:,[0]].reshape(-1,1)
F2 = X_train[:,[1]].reshape(-1,1)
## 创建meshgrid
xx , yy = np.meshgrid(np.linspace(-10, 10, 200), np.linspace(-10, 10, 200))
## 绘制散点图
plt.scatter(F1,F2)
plt.xlabel('F1')
plt.ylabel('F2')
# 创建一个dictionary并添加要用于检测异常值的所有模型:
classifiers = {
'Angle-based Outlier Detector (ABOD)' : ABOD(contamination=outlier_fraction),
'K Nearest Neighbors (KNN)' : KNN(contamination=outlier_fraction)}
# 将数据拟合到我们在dictionary中添加的每个模型,然后,查看每个模型如何检测异常值:
## 设置features的大小
plt.figure(figsize=(10, 10))
for i,(clf_name,clf) in enumerate(classifiers.items()) :
## 将数据集适合模型
clf.fit(X_train)
## 预测原始异常分数
scores_pred = clf.decision_function(X_train)*-1
## 预测DataPoint类别outlier或Inlier
y_pred = clf.predict(X_train)
## 预测中没有错误
n_errors = (y_pred != Y_train).sum()
print('No of Errors : ',clf_name,n_errors)
# 创建可视化
## 数据点为正常值还是异常值的阈值
threshold = stats.scoreatpercentile(scores_pred,100 *outlier_fraction)
## 决策功能计算每个点的原始异常分数
Z = clf.decision_function(np.c_[xx.ravel(), yy.ravel()]) * -1
Z = Z.reshape(xx.shape)
subplot = plt.subplot(1,2,i+1)
## 将蓝色ColorMAP从最小异常分数填充到阈值
subplot.contourf(xx, yy, Z, levels = np.linspace(Z.min(), threshold, 10),cmap=plt.cm.Blues_r)
## 绘制红轮廓线,其中异常分数等于阈值
a = subplot.contour(xx, yy, Z, levels=[threshold],linewidths=2, colors='red')
## 填充橙色轮廓线,其中异常分数范围从阈值到最大的异常分数
subplot.contourf(xx, yy, Z, levels=[threshold, Z.max()],colors='orange')
## 用白色点散布inliers
b = subplot.scatter(X_train[:-n_outliers, 0], X_train[:-n_outliers, 1], c='white',s=20, edgecolor='k')
## 用黑色点散布outliers
c = subplot.scatter(X_train[-n_outliers:, 0], X_train[-n_outliers:, 1], c='black',s=20, edgecolor='k')
subplot.axis('tight')
subplot.legend([a.collections[0], b, c],['learned decision function', 'true inliers', 'true outliers'],
prop=matplotlib.font_manager.FontProperties(size=10),loc='lower right')
subplot.set_title(clf_name)
subplot.set_xlim((-10, 10))
subplot.set_ylim((-10, 10))
plt.show()
No of Errors : Angle-based Outlier Detector (ABOD) 7
No of Errors : K Nearest Neighbors (KNN) 1
5、HBOS
单变量方法的组合,不能对特征之间的依赖关系进行建模,但计算速度快,对大数据集友好
-
基本假设:数据集的每个维度相互独立,对每个维度进行bin(区间)划分,区间密度越高,异常评分越低
-
算法流程:
1、为每个数据维度做出数据直方图,对分类数据统计每个值的频数并计算相对频率。对数值数据根据分布的不同采用以下两种方法:
i)静态宽度直方图:标准的直方图构建方法,在值范围内使用k个等宽箱。样本落入每个桶的频率(相对数量)作为密度(箱子高度)的估计。时间复杂度:O(n)
ii)动态宽度直方图:首先对所有值进行排序,然后固定数量的N/k个连续值装进一个箱里,其中N是总实例数,k是箱个数;直方图中的箱面积表示实例数。因为箱的宽度是由箱中第一个值和最后一个值决定的,所有箱的面积都一样,因此每一个箱的高度都是可计算的。这意味着跨度大的箱的高度低,即密度小,只有一种情况例外,超过k个数相等,此时允许在同一个箱里超过N/k值。时间复杂度:O(n×log(n))
2、对每个维度计算独立直方图,其中每个箱子的高度表示密度的估计。然后为了使得最大高度为1(确保了每个特征与异常值得分的权重相等),对直方图进行归一化处理。最后,每一个实例的HBOS值由以下公式计算:
# 生成toy example
from pyod.utils.data import generate_data,evaluate_print
contamination = 0.1 # percentage of outliers
n_train = 200 # number of training points
n_test = 100 # number of testing points
X_train, y_train, X_test, y_test = generate_data(
n_train=n_train, n_test=n_test, contamination=contamination)
# 拟合模型
from pyod.models import hbos
from pyod.utils.example import visualize
clf = hbos.HBOS()
clf.fit(X_train)
y_train_pred = clf.labels_
y_train_socres = clf.decision_scores_
y_test_pred = clf.predict(X_test) # 返回未知数据上的分类标签 (0: 正常值, 1: 异常值)
y_test_scores = clf.decision_function(X_test) # 返回未知数据上的异常值 (分值越大越异常)
# 可视化
clf_name = 'HBOS'
print("\nOn Test Data:")
evaluate_print(clf_name, y_test, y_test_scores)
visualize(clf_name, X_train, y_train, X_test, y_test, y_train_pred,y_test_pred, show_figure=True, save_figure=False)
On Test Data:
HBOS ROC:0.9667, precision @ rank n:0.7
5、总结
1.异常检测的统计学方法由数据学习模型,以区别正常的数据对象和异常点。使用统计学方法的一个优点是,异常检测可以是统计上无可非议的。当然,仅当对数据所做的统计假定满足实际约束时才为真。
2.HBOS在全局异常检测问题上表现良好,但不能检测局部异常值。但是HBOS比标准算法快得多,尤其是在大数据集上。
参考:
多元高斯分布异常点检测
OCSVW和混合高斯分布
异常检测算法演变
异常检测的N种方法
异常检测——ABOD
异常检测——HBOS
异常检测代码详解