孤立森林(Isolation Forest)

 著名的，人手一本的西瓜书(就是这本)的作者周志华老师,于2008年在第八届IEEE数据挖掘国际会议上提出孤立森林(Isolation Forest) 算法,

先简单解释一下什么是孤立森林： 

Forest算法是由南京大学的周志华和澳大利亚莫纳什大学的Fei Tony Liu，Kai Ming Ting等人共同移除，用于挖掘数据，它是适用于连续数据（Continuous numerical data）的异常检测，将异常定义为“容易被孤立的离群点（more  likely to be separated）”——可以理解为分布稀疏且离密度高的群体较远的点。用统计学来解释，在数据空间里面，分布稀疏的区域表示数据发生在此区域的概率很低，因此可以认为落在这些区域里的数据是异常的。

假设我们用一个随机超平面来切割（split）数据空间（data space）, 切一次可以生成两个子空间（想象拿刀切蛋糕一分为二）。之后我们再继续用一个随机超平面来切割每个子空间，循环下去，直到每子空间里面只有一个数据点为止。直观上来讲，我们可以发现那些密度很高的簇是可以被切很多次才会停止切割，但是那些密度很低的点很容易很早的就停到一个子空间里了。

和随机森林一样，孤立森林由 iTree(isolation tree) 组成，iTree树和随机森林的决策树不太一样，构建过程只是一个完全随机的过程。下面详细解释一下构建过程:

现有数据集中有n条数据，先从这n条数据中抽取一批样本(一般是无放回抽样)，假设样本个数 ψ 。随机选择一个特征作为起始节点，并在该特征的值域里随机选择一个值，对ψ个样本进行二叉划分，将样本中小于该取值的样本划到左分支，样本中大于该取值的划到右分支。然后在左右两个分支重复这样的二叉划分操作，直到数据集只有一条记录或者达到了树的限定高度。

由于异常数据较小且特征值和正常数据差别很大。因此，构建 iTree的时候，异常数据离根更近，而正常数据离根更远。一颗ITree的结果往往不可信，iForest算法通过多次抽样，构建多颗二叉树。最后整合所有树的结果，并取平均深度作为最终的输出深度，由此计算数据点的异常分支。

　　下图为iForest 构建 iTree 示例，异常数据点（17, 17）通常离根节点很近。

算法简介

孤立森林算法属非监督学习算法，不需要定义参数模型和进行历史训练样本，通过采用多次迭代的方式构建二叉搜索树(Binary Search Tree)，然后将这些二叉树组成森林，默认二叉搜索树的高度为 8，树的高度限制 l 与子样本数量ψ的关系为  l=ceiling(log2(ψ))，它近似等于树的平均高度。每 t=100 棵树组成一个森林，每次最多生成 ψ= 256 个森林。算法主要构建思想如下：

构建二叉树 iTree

1. 随机选择一个属性Attr
2. 随机选择该属性的一个值Value
3. 根据Attr 对每条记录进行分类，把Attr小于Value的记录放在左节点，把大于等于Value的记录放在右节点。
4. 然后递归的构造左节点和右节点，直到满足以下条件：

• (1)树达到了限制的高度；
• (2)节点上只有一个样本；
• (3)节点上的样本所有特征都相同。

Algorithm iTree(X,e,h):

Input: X-input data; e-current height; h-height limit; 
Output: an iTree;
if e >= h OR |X| <= 1 then 
    return exNode{Size <- |X|}
else
    //随机选择一个样本 q
    l <- filter(X, q<p)
    r <- filter(X, q>p)
    return inNode{Left <- iTree(l, e+1, h), Right <- iTree(r, e+1, h), SplitAttr q, SplitValue p}
end if

构建二叉树森林 iForest

构造一个iForest，iForest和Random Forest的方法有点类似，都是随机采样一部分数据集去构造一棵树，保证不同树之间的差异性，不过iForest与RF不同，采样的数据量Psi不需要等于n，可以远远小于n，论文提到采样大小超过256效果就提升不大了，并且越大还会造成计算时间上的浪费，为什么不像其他算法一样，数据越多效果越好呢？可以看看下面这两个图：

 左边是原始数据，右边是采样了数据，蓝色是正常样本，红色是异常样本。可以看到，在采样之前，正常样本和异常样本出现重叠，因此很难分开，但我们采样之后，异常样本和正常样本可以明显的分开。

其构造 iForest 的步骤如下：

• 1，从训练数据中随机选择 n 个点样本作为subsample。
• 2，根据样本数据容量迭代重复步骤(1)过程创建二叉搜索树 iTree，并将生成的 iTree 组成二叉树森林。

计算森林中二叉树的路径长度，当二叉树森林 iForest 构建完成后，就可以对样本进行预测了，预测的过程就是把测试记录在ITree上走一下（递归中序遍历），看看测试记录在哪个叶子节点。iTree能有效检测异常的假设是：异常点一般都是非常稀有的，在iTree中会很快被划分到叶子节点，并记录从根结点到叶子结点的路径长度 h(x)（也可以叫深度）。

Algorithm pathLength(x,T,e): 

Input: x-an instance;
       T-an iTree,
       e-current path length;
Output: path length of x;
if T is an external node then
    return e+c(T.size) //c(n) 为二叉树森林平均路径长度 
end if
//递归中序遍历二叉搜索树 iTree 
a <- T.splitAttr
if x.a < T.splitValue then
    return pathLength(x, T.left, e+1) 
else {x.a >= T.splitValue}
    return pathLength(x, T.right, e+1) 
end if

因此可以把叶子节点到根节点的路径h(x)长度来判断一条记录x是否是异常点（也就是根据h(x)判断x是否是异常点）。

假设 iTree 的训练样本中同样落在 x 所在叶子节点的样本数为 T.size，则数据 x 在这棵 iTree 上的路径长度 h(x)，可以用下面这个公式计算：

    

公式中，e 表示数据 x 从 iTree 的根节点到叶节点过程中经过的边的数目，C(T.size) 可以认为是一个修正值，它表示在一棵用 T.size 条样本数据构建的二叉树的平均路径长度。 

计算离群点偏离值，当森林中所有样本路径长度 h(x) 计算完毕后，通过运用统计学的方法计算得出所有数据样本期望值 E(h(x)) 和方差 S(h(x))，进而得到偏离期望和方差的异常数据点。

常见机器学习聚类算法通常根据空间距离或者密度来寻找异常数据，孤立森林算法独辟蹊径，采用构建二叉树森林再进行中序遍历计算叶子结点平均高度的方式来寻找异常数据，算法实现了对于海量数据的异常检测仅需 O(n) 的线性时间复杂度，能够在短暂的批处理时间间隔内有效检测出离群数据点。

对于一个包含n条记录的数据集，其构造的树的高度最小值为log(n)，最大值为n-1，论文中提到说用log(n)和 n-1 归一化不能保证有界和不方便比较，所以下面用一个稍微复杂一点的归一化公式：

其中H(n-1)为调和数，该值可以被估计为ln(n-1)+0.5772156649，这里的常数是欧拉常数。 E(h(x)) 表示数据 x 在多棵 iTree 的路径长度的均值，c(n)表示用来归一化E(h(x))。

S(x,n) 就是记录x在n个样本的训练数据构成的iTree的异常指数，S(x,n)取值范围为[0,1]。

下图给出了s和E(h(x))的关系：

由上图可以得到一些结论： 

1. 当E(h(x))→c(n)时，s→0.5，即样本x的路径平均长度与树的平均路径长度相近时，则不能区分是不是异常。
2. 当E(h(x))→0时，s→1，即x的异常分数接近1时，被判定为异常。
3. 当E(h(x))→n−1时，s→0，被判定为正常。

从异常分值的公式看，如果数据 x 在多棵 iTree 中的平均路径长度越短，得分越接近 1，表明数据 x 越异常；如果数据 x 在多棵 iTree 中的平均路径长度越长，得分越接近 0，表示数据 x 越正常；如果数据 x 在多棵 iTree 中的平均路径长度接近整体均值，则打分会在 0.5 附近。

注意：

• 训练样本中异常样本的比例比较高，可能最终的效果会受影响
• 异常检测跟具体的应用场景紧密相关，算法检测出的“异常”不一定是我们实际想要的，所以，在特征选择时，要注意过滤不太相关的特征。个人觉得较好的实践是，在预处理时把不关注的特征normalize

优缺点

　　具有线性时间复杂度。因为是随机森林的方法，所以可以用在含有海量数据的数据集上，通常树的数量越多，算法越稳定。由于每棵树都是互相独立生成的，因此可以部署在大规模分布式系统上来加速运算。

　　不适用于特别高维的数据。由于每次切数据空间都是随机选取一个维度和该维度的随机一个特征，建完树后仍然有大量的维度没有被使用，导致算法可靠性降低。此推荐降维后使用，或者考虑使用One Class SVM 。

　　仅对即全局稀疏点敏感。不擅长处理局部的相对稀疏点，这样在某些局部的异常点较多的时候检测可能不是很准。

　　而One Class SVM对于中小型数据分析，尤其是训练样本不是特别海量的时候用起来经常会比IForest顺手，因此比较适合做原型分析。

scikit-learn   Isolation Forest算法库概述

调用模块：

class sklearn.ensemble.IsolationForest（n_estimators = 100，max_samples ='auto'，contamination ='legacy'，max_features = 1.0，bootstrap = False，n_jobs = None，behavior ='old'，random_state = None，verbose = 0 ） 

  iForest常用参数解释：

n_estimators ： int，optional（默认值= 100）
整体中基本估算器的数量。
 
max_samples ： int或float，optional（default =“auto”）
从X中抽取的样本数量，用于训练每个基本估算器。
如果是int，则绘制max_samples样本。
如果是float，则绘制max_samples * X.shape [0]样本。
如果是“auto”，则max_samples = min（256，n_samples）。
如果max_samples大于提供的样本数，则所有样本将用于所有树（无采样）。
 
contamination ： float（0.，0.5），可选（默认值= 0.1）
数据集的污染量，即数据集中异常值的比例。在拟合时用于定义决策函数的阈值。如果是“自动”，则确定决策函数阈值，如原始论文中所示。
 
在版本0.20中更改：默认值contamination将从0.20更改为'auto'0.22。
 
max_features ： int或float，可选（默认值= 1.0）
从X绘制以训练每个基本估计器的特征数。
 
如果是int，则绘制max_features特征。
如果是float，则绘制max_features * X.shape [1]特征。
bootstrap ： boolean，optional（default = False）
如果为True，则单个树适合于通过替换采样的训练数据的随机子集。如果为假，则执行未更换的采样。
 
n_jobs ： int或None，可选（默认=无）
适合和预测并行运行的作业数。 None除非在joblib.parallel_backend上下文中，否则表示1 。 -1表示使用所有处理器。
 
random_state ： int，RandomState实例或None，可选（默认=无）
如果是int，则random_state是随机数生成器使用的种子; 如果是RandomState实例，则random_state是随机数生成器; 如果没有，随机数生成器所使用的RandomState实例np.random。
 
verbose ： int，optional（默认值= 0）
控制树构建过程的详细程度。

 iForest常用属性说明：

• estimators_：构造好的子树的集合
• estimators_samples：每个子树抽取的样本的子集
• max_samples_：样本的真正数量
• offset_：float

       offset用来从原始分数开始定义决策函数，其关系是decision_function = score_samples - offset_。假设behaviour == ‘new’，则offset_如下定义：

　　当contamination参数设置为'auto'，当inliers的得分接近0且outliers的得分接近-1时，偏移量等于-0.5；
       当提供与“auto”不同的contamination参数时，则以在训练中获取期望的异常个数的方式来定义偏移量(决策函数< 0的样本)。
       假设behaviour ==“old”，我们总是有offset_ = -0.5，使得决策函数独立于contamination参数。

注意：其实现是基于一系列的ExtraTreeRegressor。每个树的最大深度设置为ceil(log_2(n))，其中n是用于构建树的样本数量

  iForest常用方法介绍：

• score_samples（X）：与原文定义的异常值相反
• set_params（**params）：设置该森林的参数

实例一（iForest算法检验数据）

　　isolation Forest 通过随机选择一个特征然后随机选择所选特征的最大值和最小值之间的分割值来“隔离”观察。

　　由于递归分区可以表示为一个树结构，分裂需要隔离一个样本的数量相当于从根节点到终止节点路径长度。

　　随机划分为异常生产明显更短的路径。因此，当一个随机森林产生短路径长度，他们极有可能是异常点。

　　代码：

#_*_coding:utf-8_*_
import numpy as np
import matplotlib.pyplot as plt
from sklearn.ensemble import IsolationForest
 
rng = np.random.RandomState(42)
 
# Generate  train data
X = 0.3 * rng.randn(100, 2)
X_train = np.r_[X + 2, X - 2]
X = 0.3 * rng.randn(20, 2)
X_test = np.r_[X + 2, X - 2]
X_outliers = rng.uniform(low=-4, high=4, size=(20, 2))
 
# fit the model
clf = IsolationForest(behaviour='new', max_samples=100,
                      random_state=rng, contamination='auto')
clf.fit(X_train)
y_pred_train = clf.predict(X_train)
y_pred_test = clf.predict(X_outliers)
 
xx, yy = np.meshgrid(np.linspace(-5, 5, 50), np.linspace(-5, 5, 50))
Z = clf.decision_function(np.c_[xx.ravel(), yy.ravel()])
Z = Z.reshape(xx.shape)
 
plt.title("IsolationForest")
plt.contourf(xx, yy, Z, camp=plt.cm.Blues_r)
b1 = plt.scatter(X_train[:, 0], X_train[:, 1], c='white',
                 s=20, edgecolor='k')
b2 = plt.scatter(X_test[:, 0], X_test[:, 1], c='green',
                 s=20, edgecolor='k')
c = plt.scatter(X_outliers[:, 0], X_outliers[:, 1], c='red',
                s=20, edgecolor='k')
plt.axis('tight')
plt.xlim((-5, 5))
plt.ylim((-5, 5))
plt.legend([b1, b2, c],
           ["training observations",
            "new regular observations", "new abnormal observations"],
           loc="upper left")
plt.show()

　　结果展示：

 实例二（多种异常检测算法比较）

　　此示例显示了2D数据集上不同异常检测算法的特征。数据集包含一种或两种模式（高密度区域），以说明算法处理多模态数据的能力。

 　　对于每个数据集，15%的样本被生成为随机均匀噪声。该比例是给予OneClassSVM的nu参数的值和其他异常值检测算法的污染参数。除了局部异常因子（LOF）之外，内部和异常值之间的决策边界以黑色显示，因为当用于异常值检测时，它没有预测方法应用于新数据。

　　svm.OneClassSVM被称为是对异常值敏感，并因此对异常值检测不执行的非常好。当训练集未被异常值污染时，该估计器最适合于新颖性检测。也就是说，高维中的离群检测，或者对上层数据的分布没有任何假设是非常具有挑战性的，并且单类SVM可能在这些情况下根据其超参数的值给出有用的结果。

　　covariance.EllipticEnvelope假设数据是高斯数据并学习椭圆。因此，当数据不是单峰时，它会降级。但请注意，此估算器对异常值很稳健。

　　ensemble.IsolationForest并且neighbors.LocalOutlierFactor 似乎对多模态数据集表现得相当好。neighbors.LocalOutlierFactor对于第三数据集示出了优于其他估计器的优点 ，其中两种模式具有不同的密度。这一优势可以通过LOF的局部方面来解释，这意味着它只将一个样本的异常得分与其邻居的得分进行比较。

　　最后，对于最后一个数据集，很难说一个样本比另一个样本更异常，因为它们均匀分布在超立方体中。除了svm.OneClassSVM稍微过度拟合之外，所有估算者都为这种情况提供了不错的解决方案。在这种情况下，更仔细地观察样本的异常分数是明智的，因为良好的估计器应该为所有样本分配相似的分数。

　　虽然这些例子给出了一些关于算法的直觉，但这种直觉可能不适用于非常高维的数据。

　　最后，请注意模型的参数已经在这里精心挑选，但实际上它们需要进行调整。在没有标记数据的情况下，问题完全没有监督，因此模型选择可能是一个挑战。

　　代码：

import time
 
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
 
from sklearn import svm
from sklearn.datasets import make_blobs, make_moons
from sklearn.covariance import EllipticEnvelope
from sklearn.ensemble import IsolationForest
from sklearn.neighbors import LocalOutlierFactor
 
matplotlib.rcParams['contour.negative_linestyle'] = 'solid'
 
# Example settings
n_samples = 300
outliers_fraction = 0.15
n_outliers = int(outliers_fraction * n_samples)
n_inliers = n_samples - n_outliers
 
# define outlier/ anomaly detection methods to be compared
anomaly_algorithms = [
    ("Robust covariance", EllipticEnvelope(contamination=outliers_fraction)),
    ("One-Class SVM", svm.OneClassSVM(nu=outliers_fraction, kernel='rbf',gamma=0.1)),
    ("Isolation Forest", IsolationForest(behaviour='new', contamination=outliers_fraction, random_state=42)),
    ("Local Outlier Factor", LocalOutlierFactor(n_neighbors=35, contamination=outliers_fraction))
]
 
# define datasets
blobs_params = dict(random_state=0, n_samples=n_inliers, n_features=2)
datasets = [
    make_blobs(centers=[[0, 0], [0, 0]], cluster_std=0.5, **blobs_params)[0],
    make_blobs(centers=[[2, 2], [-2, -2]], cluster_std=[0.5, 0.5], **blobs_params)[0],
    make_blobs(centers=[[2, 2], [-2, -2]], cluster_std=[1.5, 0.3], **blobs_params)[0],
    4. * (make_moons(n_samples=n_samples, noise=0.05, random_state=0)[0] - np.array([0.5, 0.25])),
    14. * (np.random.RandomState(42).rand(n_samples, 2) - 0.5)
]
 
# Compare given classifiers under given settings
xx, yy = np.meshgrid(np.linspace(-7, 7, 150), np.linspace(-7, 7, 150))
 
plt.figure(figsize=(len(anomaly_algorithms) * 2 + 3, 12.5))
plt.subplots_adjust(left=0.02, right=0.98, bottom=0.001, top=0.96, wspace=0.05, hspace=0.01)
 
plot_num = 1
rng = np.random.RandomState(42)
 
 
for i_dataset, X in enumerate(datasets):
    # add outliers
    X = np.concatenate([X, rng.uniform(low=-6, high=6, size=(n_outliers, 2))], axis=0)
    for name, algorithm in anomaly_algorithms:
        print(name , algorithm)
        t0 = time.time()
        algorithm.fit(X)
        t1 = time.time()
        plt.subplot(len(datasets), len(anomaly_algorithms), plot_num)
        if i_dataset == 0:
            plt.title(name, size=18)
 
        # fit the data and tag outliers
        if name == 'Local Outlier Factor':
            y_pred = algorithm.fit_predict(X)
        else:
            y_pred = algorithm.fit(X).predict(X)
 
        # plot the levels lines and the points
        if name != "Local Outlier Factor":
            Z = algorithm.predict(np.c_[xx.ravel(), yy.ravel()])
            Z = Z.reshape(xx.shape)
            plt.contour(xx, yy, Z, levels=[0], linewidths=2, colors='black')
 
        colors = np.array(["#377eb8", '#ff7f00'])
        plt.scatter(X[:, 0], X[:, 1], s=10, color=colors[(y_pred + 1) // 2])
 
        plt.xlim(-7, 7)
        plt.ylim(-7, 7)
        plt.xticks(())
        plt.yticks(())
        plt.text(0.99, 0.01, ('%.2fs' % (t1 - t0)).lstrip('0'),
                 transform=plt.gca().transAxes, size=15,
                 horizontalalignment='right')
        plot_num += 1
 
plt.show()

　　结果

参考资料：

• https://www.kesci.com/home/project/5ab37cefba5a9d0f1b5f9078
• https://github.com/lishiyun19/AnomalyFilter
• https://blog.csdn.net/extremebingo/article/details/80108247
• https://www.cnblogs.com/wj-1314/p/10461816.html