目录
1. 引言
真实数据集中不同维度的数据通常具有高度的相关性，这是因为不同的属性往往是由相同的基础过程以密切相关的方式产生的。在古典统计学中，这被称为——回归建模，一种参数化的相关性分析
- 一类相关性分析试图通过其他变量预测单独的属性值，如线性回归
- 另一类方法用一些潜在变量来代表整个数据，例如主成分分析
需要明确的是，这里有两个重要的假设：
- 假设一：近似线性相关假设。线性相关假设是使用两种模型进行一场检测的重要理论基础。
- 假设2：子空间假设。子空间假设认为数据是镶嵌在低维子空间中的，线性方法的目的是找到合适的低维子空间使得异常点（o）在其中区别于正常点（n）
基于这两点假设，在异常检测的第一阶段，为了确定模型是否适合特定的数据集，对数据进行探索性和可视化分析是非常关键的
2. 数据可视化
数据集下载https://download.csdn.net/download/weixin_39940512/14121801本例是用的其中的breast-cancer-unsupervised-ad.csv数据集
##导入相关包
import warnings
warnings.filterwarnings('ignore')
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
import missingno as msno
##载入训练集与测试集
path='./data/'
Train_data=pd.read_csv(path+'breast-cancer-unsupervised-ad.csv')
##简单观察数据
print(Train_data.shape)
Train_data.head()
有30个特征和一个标签，数据集大小为367*31
##使用describe()查看相关统计量
Train_data.describe()
##通过info（）来查看数据类型
Train_data.info()
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 367 entries, 0 to 366
'''
Data columns (total 31 columns):
# Column Non-Null Count Dtype
--- ------ -------------- -----
0 f0 367 non-null float64
1 f1 367 non-null float64
2 f2 367 non-null float64
3 f3 367 non-null float64
4 f4 367 non-null float64
5 f5 367 non-null float64
6 f6 367 non-null float64
7 f7 367 non-null float64
8 f8 367 non-null float64
9 f9 367 non-null float64
10 f10 367 non-null float64
11 f11 367 non-null float64
12 f12 367 non-null float64
13 f13 367 non-null float64
14 f14 367 non-null float64
15 f15 367 non-null float64
16 f16 367 non-null float64
17 f17 367 non-null float64
18 f18 367 non-null float64
19 f19 367 non-null float64
20 f20 367 non-null float64
21 f21 367 non-null float64
22 f22 367 non-null float64
23 f23 367 non-null float64
24 f24 367 non-null float64
25 f25 367 non-null float64
26 f26 367 non-null float64
27 f27 367 non-null float64
28 f28 367 non-null float64
29 f29 367 non-null float64
30 label 367 non-null object
dtypes: float64(30), object(1)
memory usage: 89.0+ KB
'''
可以看出除标签外都为float64类型
numeric_features=['f'+str(i) for i in range(30)]#特征名
##相关性分析
#查看特征之间的相关性
numeric=Train_data[numeric_features]
correlation=numeric.corr()
#绘制热力图
f , ax = plt.subplots(figsize = (14, 14))
sns.heatmap(correlation)
plt.title('Correlation of Numeric Features with Price',y=1,size=16)
plt.show()
##每个数字特征分布可视化
'''
绘图 https://blog.csdn.net/weixin_42398658/article/details/82960379
绘图 https://www.jianshu.com/p/6a210c2ad3ad
'''
'''
数据分析的时候经常要把宽数据--->>长数据，有点像你们用excel 做透视跟逆透视的过程，直接看下面例子，希望有助于理解.
pandas.melt 使用参数：
pandas.melt(frame, id_vars=None, value_vars=None, var_name=None, value_name='value', col_level=None)
参数解释：
frame:要处理的数据集。
id_vars:不需要被转换的列名。
value_vars:需要转换的列名，如果剩下的列全部都要转换，就不用写了。
var_name和value_name是自定义设置对应的列名。
col_level :如果列是MultiIndex，则使用此级别。
'''
print(Train_data['f0'].max(),Train_data['f0'].min())
f=pd.melt(Train_data,value_vars=numeric_features)
g = sns.FacetGrid(f, col="variable", col_wrap=6, sharex=False, sharey=False)
g = g.map(sns.distplot, "value", hist=False, rug=True)
##变量两两之间的相关性
sns.set()
sns.pairplot(Train_data[numeric_features],size=2,
kind='scatter',diag_kind='kde') #kde为密度图
##数据降维可视化
from sklearn.manifold import TSNE
tsne = TSNE(n_components=2, init='pca', random_state=0)
result = tsne.fit_transform(numeric)
x_min, x_max = np.min(result, 0), np.max(result, 0)
result = (result - x_min) / (x_max - x_min)
label = Train_data['label']
fig=plt.figure(figsize=(7,7))
color = {'o':0, 'n':7}
for i in range(result.shape[0]):
plt.text(result[i, 0], result[i, 1], str(label[i]),
color=plt.cm.Set1(color[label[i]] / 10.),
fontdict={'weight': 'bold', 'size': 9})
plt.xticks([]) #不显示刻度线
plt.yticks([])
plt.title('Visualization of data dimension reduction')
3. 线性回归
在线性回归中，我们假设不同维度的变量具有一样的相关性，并可以通过一个相关系数矩阵进行衡量。因此对于特定的观测值，可以通过线性方程组来建模。在实际应用中，观测值的数量往往大于数据的维度，导致线性方程组一定是超定方程，不能直接求解。因此需要通过优化的方法，最小化模型预测值与真实数据点的误差。
线性回归是统计学中的一个重要应用，它能够通过一系列的自变量去预测一个特殊的因变量的值。在这种情况下，异常值是根据其他自变量对因变量的影响来定义的，而自变量之间相互关系中的异常则不那么重要。这里异常点检测主要用于数据降噪，避免异常点的出现对模型性能的影响，因而这里关注的兴趣点主要是正常值(n)。
而异常检测中并不会对任何变量给与特殊对待，异常值的定义是基于基础数据点的整体分布，这里我们关注的兴趣点主要是异常值(o)。
广义的回归建模只是一种工具，这种工具既可以用来进行数据降噪也可以进行异常点检测。
3.1 基于自变量和因变量的线性回归
3.1.1 最小二乘法
此处以一元线性回归为例：
Y
=
∑
i
=
1
d
a
i
⋅
X
i
+
a
d
+
1
Y=\sum_{i=1}^{d} a_{i} \cdot X_{i}+a_{d+1}
Y=i=1∑dai⋅Xi+ad+1
变量Y为因变量，也就是我们要预测的值； X 1 . . . X d X_{1}...X_{d} X1...Xd为一系列因变量，也就是输入值。系数 a 1 . . . a d + 1 a_{1}...a_{d+1} a1...ad+1为要学习的参数。假设数据共包含 N N N个样本，第 j j j个样本包含的数据为 x j 1 . . . x j d x_{j1}...x_{jd} xj1...xjd和 y j y_{j} yj，带入式(1)如下式所示：
y
j
=
∑
i
=
1
d
a
i
⋅
x
j
i
+
a
d
+
1
+
ϵ
j
y_{j}=\sum_{i=1}^{d} a_{i} \cdot x_{j i}+a_{d+1}+\epsilon_{j}
yj=i=1∑dai⋅xji+ad+1+ϵj
这里
ϵ
j
\epsilon_{j}
ϵj为第
j
j
j个样本的误差。以
Y
Y
Y 代表
N
×
1
N \times 1
N×1 的因变量矩阵
(
y
1
.
.
.
y
N
)
T
{(y_{1}...y_{N})}^{T}
(y1...yN)T，即样本中的真实值；以
U
U
U代表
N
×
(
d
+
1
)
N \times (d+1)
N×(d+1)的自变量矩阵，其中第
j
j
j行为
(
x
j
1
.
.
.
x
j
d
,
1
)
(x_{j1}...x_{jd}, 1)
(xj1...xjd,1)；以
A
A
A 代表
(
d
+
1
)
×
1
(d+1) \times 1
(d+1)×1 的系数矩阵
(
a
1
.
.
.
a
d
+
1
)
T
(a_{1}...a_{d+1})^{T}
(a1...ad+1)T。则模型可表示为：
f
(
U
,
A
)
=
U
⋅
A
f(U, A) = U \cdot A
f(U,A)=U⋅A
定义目标函数为：
L ( A ) = 1 2 ∥ Y − U ⋅ A ∥ 2 L(A) = \frac{1}{2}{\left\| {Y - U \cdot A} \right\|^2} L(A)=21∥Y−U⋅A∥2
目标函数是关于 A A A的凸函数，其对 A A A求偏导为：
∂ L ( A ) ∂ A = 1 2 ∂ ∥ Y − U ⋅ A ∥ 2 ∂ A = − U T ( Y − U ⋅ A ) \frac{{\partial L(A)}}{{\partial A}} = \frac{1}{2}\frac{{\partial {{\left\| {Y - U \cdot A} \right\|}^2}}}{{\partial A}} = - {U^T}(Y - U \cdot A) ∂A∂L(A)=21∂A∂∥Y−U⋅A∥2=−UT(Y−U⋅A)
令 ∂ L ( A ) ∂ A = 0 \frac{{\partial L(A)}}{{\partial A}}=0 ∂A∂L(A)=0，得到最优参数为：
A = ( U T ⋅ U ) − 1 ⋅ ( U T ⋅ Y ) A=\left(U^{T} \cdot U\right)^{-1} \cdot\left(U^{T} \cdot Y\right) A=(UT⋅U)−1⋅(UT⋅Y)
这种求解线性回归参数的方法也叫最小二乘法。
最小二乘法要求矩阵 U T ⋅ U U^{T} \cdot U UT⋅U 可逆，即 U T ⋅ U U^{T} \cdot U UT⋅U是满秩的。当 U T ⋅ U U^{T} \cdot U UT⋅U不可逆时可以通过两种方法进行参数估计
- 一种先使用主成分分析等方法来预处理数据，消除不同特征之间的相关性，然后再使用最小二乘法
- 第二种方法是使用梯度下降法。
3.1.2 梯度下降法
数据集
- 训练集
- 测试集
- 验证集
损失函数
如果把线性回归看作是一个优化问题，那么我们要优化的目标就是损失函数。损失函数是用来衡量样本误差的函数，我们的优化目标是要求得在误差最小的情况下模型参数的值。
【注意】损失函数和代价函数的区别
Loss Function(损失函数)：the error for single training example;单个训练样本的误差
Cost Function(代价函数)：the average of the loss functions of the entire training set;整个训练集损失函数的平均值；
线性回归中常用的是均方误差
l
(
i
)
(
w
,
b
)
=
1
2
(
y
^
(
i
)
−
y
(
i
)
)
2
l^{(i)}(\mathbf{w}, b)=\frac{1}{2}\left(\hat{y}^{(i)}-y^{(i)}\right)^{2}
l(i)(w,b)=21(y^(i)−y(i))2
L
(
w
,
b
)
=
1
n
∑
i
=
1
n
l
(
i
)
(
w
,
b
)
=
1
n
∑
i
=
1
n
1
2
(
w
⊤
x
(
i
)
+
b
−
y
(
i
)
)
2
L(\mathbf{w}, b)=\frac{1}{n} \sum_{i=1}^{n} l^{(i)}(\mathbf{w}, b)=\frac{1}{n} \sum_{i=1}^{n} \frac{1}{2}\left(\mathbf{w}^{\top} \mathbf{x}^{(i)}+b-y^{(i)}\right)^{2}
L(w,b)=n1i=1∑nl(i)(w,b)=n1i=1∑n21(w⊤x(i)+b−y(i))2
其中
y
^
\hat{y}
y^ 为预测值，
y
y
y 为真实值。
优化算法 - 随机梯度下降
当模型和损失函数形式较为简单时，上面的误差最小化问题的解可以直接用公式表达出来。这类解叫作解析解（analytical solution）。本节使用的线性回归和平方误差刚好属于这个范畴。
然而，大多数深度学习模型并没有解析解，只能通过优化算法有限次迭代模型参数来尽可能降低损失函数的值。这类解叫作数值解（numerical solution）。
在求数值解的优化算法中，==小批量随机梯度下降（mini-batch stochastic gradient descent）==被广泛使用。它的算法很简单
- 先选取一组模型参数的初始值，如随机选取；
- 接下来对参数进行多次迭代，使每次迭代都可能降低损失函数的值。
- 在每次迭代中，先随机均匀采样一个由固定数目训练数据样本所组成的小批量（mini-batch），然后求小批量中数据样本的平均损失和有关模型参数的导数（梯度）
- 最后用此结果与预先设定的学习率的乘积作为模型参数在本次迭代的减小量。如下式所示：
( w , b ) ← ( w , b ) − η ∣ B ∣ ∑ i ∈ B ∂ ( w , b ) l ( i ) ( w , b ) (\mathbf{w}, b) \leftarrow(\mathbf{w}, b)-\frac{\eta}{|\mathcal{B}|} \sum_{i \in \mathcal{B}} \partial_{(\mathbf{w}, b)} l^{(i)}(\mathbf{w}, b) (w,b)←(w,b)−∣B∣η∑i∈B∂(w,b)l(i)(w,b)
学习率( η \eta η): 代表在每次优化中，能够学习的步长的大小
批量大小( B B B): 是小批量计算中的批量大小batch size
3.2 基于异常检测的线性回归
一个特定的变量被认为是特殊的，最优平面是通过最小化该特殊变量的均方误差而确定的。
而我们通常所说的异常检测中并不会对任何变量给与特殊对待，异常值的定义是基于基础数据点的整体分布，因此需要采用一种更一般的回归建模：即以相似的方式对待所有变量，通过最小化数据对该平面的投影误差确定最佳回归平面。在这种情况下，假设我们有一组变量
X
1
…
X
d
X_{1}… X_{d}
X1…Xd， 对应的回归平面如下：
a 1 ⋅ X 1 + … + a d ⋅ X d + a d + 1 = 0 a_{1} \cdot X_{1}+\ldots+a_{d} \cdot X_{d}+a_{d+1}=0 a1⋅X1+…+ad⋅Xd+ad+1=0
为了后续计算的方便，对参数进行如下约束：
∑
i
=
1
d
a
i
2
=
1
\sum\limits_{i = 1}^d {a_i^2 = 1}
i=1∑dai2=1
以
L
2
L_{2}
L2范数作为目标函数：
L
=
∥
U
⋅
A
∥
2
L = {\left\| {U \cdot A} \right\|_2}
L=∥U⋅A∥2
&emsp这样的一个问题可以通过主成分分析方法得到有效解决
4. 主成分分析
上一节的最小二乘法试图找到一个与数据具有最佳匹配 ( d − 1 ) (d−1) (d−1) 维超平面。主成分分析方法可用于解决这一问题的广义版本。具体来说，它可以找到任意 k ( k < d ) k( k<d ) k(k<d) 维的最优表示超平面，从而使平方投影误差最小化。
4.1 原理推导
对于
d
d
d 维，包含
N
N
N 个样本的数据，用
R
i
R_{i}
Ri 表示其中第
i
i
i 行为：
[
x
i
1
.
.
.
x
i
d
]
[x_{i1}... x_{id}]
[xi1...xid]。由此可以得到
d
×
d
d \times d
d×d 的协方差矩阵（标准的PCA应当计算相关系数矩阵，即对数据进行均值为0方差为1的标准化处理，而协方差矩阵只需要减去均值即可）:
Σ
=
(
R
−
R
ˉ
)
T
⋅
(
R
−
R
ˉ
)
Σ = (R - \bar{R})^{T} \cdot (R - \bar{R})
Σ=(R−Rˉ)T⋅(R−Rˉ)
易知协方差矩阵
Σ
Σ
Σ 是对称并且半正定的，因此可以进行相似对角化：
Σ
=
P
⋅
D
⋅
P
T
Σ = P \cdot D \cdot P^{T}
Σ=P⋅D⋅PT
这里的
D
D
D 为对角矩阵，对角元素为特征值；
P
P
P 为标准正交矩阵，每一行为对应的特征向量；这些标准正交向量提供了数据应该投影的轴线方向。与异常检测相关的主成分分析的主要性质如下：
- 如果前 k k k 的特征向量选定之后（根据最大的 k k k个特征值），由这些特征向量定义的 k k k 维超平面是在所有维度为 k k k 的超平面中，所有数据点到它的均方距离尽可能小的平面。
- 如果将数据转换为与正交特征向量对应的轴系，则转换后的数据沿每个特征向量维的方差等于相应的特征值。在这种新表示中，转换后的数据的协方差为0。
- 由于沿特征值小的特征向量的转换数据的方差很低，因此沿这些方向的变换数据与平均值的显著偏差可能表示离群值。
在得到这些特征值和特征向量之后，可以将数据转换到新的坐标系中。以
Y
1
.
.
.
Y
N
Y_{1}...Y_{N}
Y1...YN 表示新坐标系中的数据，这些数据可以通过原始向量
R
i
R_{i}
Ri 与包含新轴系的标准正交特征向量矩阵
P
P
P 的乘积来实现。
Y
i
=
R
i
⋅
P
{Y_i} = {R_i} \cdot P
Yi=Ri⋅P
在许多涉及高维数据集的真实场景中，很大一部分特征值往往非常接近于零。这意味着大多数数据都沿着一个低维的子空间排列。从异常检测的角度来看，这是非常方便的，因为离这些投影方向非常远的观测值可以被假定为离群值。例如，对于特征值较小（方差较小）的特征向量
j
j
j，第
i
i
i 条记录的
y
i
j
y_{ij}
yij 与
y
k
j
y_{kj}
ykj 的其他值的偏差较大，说明有离群行为。这是因为当
j
j
j固定而
k
k
k变化时，
y
k
j
y_{kj}
ykj 的值应当变化不大。因此，
y
i
j
y_{ij}
yij 值是不常见的。
在不选取任何特定的
k
k
k 维集合的情况下，一种更精确的异常检测建模方法是使用特征值来计算数据点沿每个主分量方向到质心的归一化距离。设
e
j
e_{j}
ej为第
j
j
j 个特征向量，
λ
j
λ_{j}
λj 为沿该方向的方差(特征值)。数据点
X
ˉ
\bar{X}
Xˉ相对于对数据质心
μ
ˉ
\bar{\mu}
μˉ的总体归一化异常得分可以由下式给出:
S
core
(
X
ˉ
)
=
∑
j
=
1
d
∣
(
X
ˉ
−
μ
ˉ
)
⋅
e
ˉ
j
∣
2
λ
j
S \operatorname{core}(\bar{X})=\sum_{j=1}^{d} \frac{\left|(\bar{X}-\bar{\mu}) \cdot \bar{e}_{j}\right|^{2}}{\lambda_{j}}
Score(Xˉ)=j=1∑dλj∣∣(Xˉ−μˉ)⋅eˉj∣∣2
值得注意的是，对异常得分的大部分贡献是由
λ
j
λ_{j}
λj 值较小的主成分的偏差提供的，这一点上文中有提及过。主成分分析比因变量回归能更稳定地处理少数异常值的存在。这是因为主成分分析是根据最优超平面来计算误差的，而不是一个特定的变量。当数据中加入更多的离群点时，最优超平面的变化通常不会大到影响离群点的选择。因此，这种方法更有可能选择正确的异常值，因为回归模型一开始就更准确。
4.2 归一化问题
当不同维度的尺度差别较大时，使用
P
C
A
PCA
PCA 有时并不能得到直观有效的结果。
例如，考虑一个包含年龄和工资等属性的人口统计数据集。工资属性的范围可能是几万，而年龄属性几乎总是小于100，使用主成分分析会导致主成分被高方差属性所控制。对于一个只包含年龄和工资的二维数据集，最大的特征向量几乎与工资轴平行，这会降低异常点检测过程的有效性。因此，一个自然的解决方案是对数据进行均值为0方差为1的标准化处理。这隐含地导致在主成分分析中使用相关矩阵而不是协方差矩阵。当然，这个问题并不是线性建模所独有的，对于大多数异常检测算法，都需要使用这样的预处理。
5.回归分析的局限性
回归分析作为检测离群值的工具有一些局限性：
- 了使回归分析技术有效，数据需要高度相关，并沿着低维子空间对齐
- 当数据不相关，但在某些区域高度聚集时，这种方法可能不会有效。
另一个相关的问题是，数据中的相关性在本质上可能不是全局性的。最近的一些分析观察**[1]**表明，子空间相关性是特定于数据的特定位置的。在这种情况下，由主成分分析发现的全局子空间对于异常检测是次优的。因此，为了创建更一般的局部子空间模型，有时将线性模型与邻近模型结合起来是有用的。
6. 总结
真实数据中，数据不同属性之间往往具有显著的相关性。在这种情况下，线性建模可以提供一种有效的工具来从底层数据中移除异常值或者进行异常检测。对于其他基于因变量回归的应用，线性建模是一种工具，去除异常值对于提高此类应用的性能是非常重要的。在大多数情况下，主成分分析提供了去除异常值和进行异常检测最有效的方法，因为它对存在少数异常值的数据更有鲁棒性。
鲁棒性：https://blog.csdn.net/zhouxinxin0202/article/details/79915873
Huber从稳健统计的角度系统地给出了鲁棒性3个层面的概念：
- 一是模型具有较高的精度或有效性，这也是对于机器学习中所有学习模型的基本要求；
- 二是对于模型假设出现的较小偏差，只能对算法性能产生较小的影响；
主要是：噪声- 三是对于模型假设出现的较大偏差，不可对算法性能产生“灾难性”的影响。
主要是：离群点（outlier）
7. 练习
使用pyod库生成example并使用该库的pca模块进行检测
https://mybinder.org/v2/gh/yzhao062/pyod/master
##pca example
from pyod.models.pca import PCA
from pyod.utils.data import generate_data
from pyod.utils.data import evaluate_print
from pyod.utils.example import visualize
if __name__=='__main__':
contamination=0.1 #异常值百分比
n_train=200
n_test=100
X_train,y_train,X_test,y_test=generate_data(n_train=n_train,
n_test=n_test,
n_features=20,
contamination=contamination,
random_state=42)
clf_name='PCA'
clf=PCA(n_components=3)
clf.fit(X_train)
#得到训练数据的预测标签和异常值得分
y_train_pred=clf.predict(X_train)
y_train_scores=clf.decision_scores_
#对测试数据进行预测
y_test_pred = clf.predict(X_test) # outlier labels (0 or 1)
y_test_scores = clf.decision_function(X_test) # outlier scores
# evaluate and print the results
print("\nOn Training Data:")
evaluate_print(clf_name, y_train, y_train_scores)
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)
参考文献
[1]:Outlier Analysis Charu C. Aggarwal
[2]:Anomaly Detection: A Survey VARUN CHANDOLA, ARINDAM BANERJEE, and VIPIN KUMAR University of Minnesota
[3]:Outlier Analysis Charu C. Aggarwal
[4]:Anomaly Detection: A Tutorial
[5]:Data Mining Concepts and Techniques Third Edition
133
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
