六、异常检测
在本章中，我们将讨论无监督学习的实际应用。 我们的目标是训练模型，这些模型要么能够重现特定数据生成过程的概率密度函数，要么能够识别给定的新样本是内部数据还是外部数据。 一般而言，我们可以说，我们要追求的特定目标是发现异常，这些异常通常是在模型下不太可能出现的样本（也就是说，给定概率分布p(x) << λ，其中λ是预定义的阈值），或者离主分布的质心很远。
特别是，本章将包含以下主题：
- 概率密度函数及其基本性质简介
- 直方图及其局限性
- 核密度估计（KDE）
- 带宽选择标准
- 异常检测的单变量示例
- 使用 KDD Cup 99 数据集的 HTTP 攻击异常检测示例
- 单类支持向量机
- 隔离森林的异常检测
技术要求
本章中提供的代码要求：
- Python3.5+（强烈建议使用 Anaconda 发行版）
- 库：
- SciPy 0.19+
- NumPy 1.10+
- Scikit-Learn 0.20+
- Pandas 0.22+
- Matplotlib 2.0+
- Seaborn 0.9+
概率密度函数
在所有先前的章节中，我们一直认为我们的数据集是从隐式数据生成过程p_data以及所有算法假设x[i] ∈ X为独立同分布的（IID）并进行均匀采样。 我们假设X足够准确地表示p_data，以便算法可以学习使用有限的初始知识进行概括。 相反，在本章中，我们感兴趣的是直接建模p_data，而没有任何具体限制（例如，高斯混合模型通过对数据结构施加约束来实现此目标分布）。 在讨论一些非常有效的方法之前，简要回顾一下在可测量子集X包含于ℜ^n上定义的通用连续概率密度函数p(x)的性质很有帮助（为了避免混淆，我们将用p(x)表示密度函数，用P(x)表示实际概率）：
例如，单变量高斯分布完全由均值μ和方差σ^2来表征：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-GbGW0bw5-1681652675127)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/b3f2250c-4a29-4b19-9652-1af23703d1d2.png)]
因此，x ∈ (a, b)的概率如下：
即使连续空间（例如，高斯）中某个事件的绝对概率为零（因为积分具有相同的极值），概率密度函数还是一种非常有用的度量，可以用来将一个样本与另一个对比来了解它。 例如：考虑高斯分布N(0, 1)，密度p(1)= 0.4，而对于x = 2密度降低到大约0.05。 这意味着1的可能性比2高0.4 / 0.05 = 8倍。 同样，我们可以设置可接受阈值α并定义所有x[i]样本，p(x[i]) < α的这些样本为异常（例如，在我们的情况下，α = 0.01）。 这种选择是异常检测过程中的关键步骤，正如我们将要讨论的那样，它还必须包括潜在的异常值，但是这些异常值仍然是常规样本。
在许多情况下，特征向量是使用多维随机变量建模的。 例如：数据集X包含于R^3可以用联合概率密度函数p(x, y, z)表示。 在一般情况下，实际概率需要三重积分：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-KE1EJsKK-1681652675128)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/04050c08-87c3-4985-8fa1-45242b98a498.png)]
不难理解，任何使用这种联合概率的算法都会受到复杂性的负面影响。 通过假设单个组件的统计独立性可以大大简化：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-5bsdj9x2-1681652675128)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/35903294-84ca-486e-856c-930c74b81e27.png)]
不熟悉此概念的读者可以想象考试前会有大量学生。 用随机变量建模的特征是学习时间（x）和完成课程的数量（y），鉴于这些因素，我们希望找出成功的可能性p(Success | x, y)（此类示例基于条件概率，但主要概念始终相同）。 我们可以假设一个完成所有课程的学生需要在家少学习； 但是，这样的选择意味着两个因素之间的依赖性（和相关性），不能再单独评估了。 相反，我们可以通过假设不存在任何相关性来简化程序，并根据给定的上课次数和作业时间与成功的边际概率进行比较。 重要的是要记住，特征之间的独立性不同于随后从分布中抽取的样本的独立性。 当我们说数据集由 IID 样本组成时，是指概率p(x[i] | x[i-1], x[i-2], ..., p[1]) = p(x[i])。 换句话说，我们假设样本之间没有相关性。 这样的条件更容易实现，因为通常足以洗净数据集以删除任何残余相关性。 取而代之的是，特征之间的相关性是数据生成过程的特殊属性，无法删除。 因此，在某些情况下，我们假定独立性是因为我们知道其影响可以忽略不计，并且最终结果不会受到严重影响，而在其他情况下，我们将基于整个多维特征向量训练模型。 现在，我们可以定义将在其余部分中使用的异常的概念。
作为异常值或新颖性的异常
本章的主题是在没有任何监督的情况下自动检测异常。 由于模型不是基于标记样本提供的反馈，因此我们只能依靠整个数据集的属性来找出相似之处并突出显示不同之处。 特别是，我们从一个非常简单但有效的假设开始：常见事件为正常，而不太可能发生的事件通常被视为异常。 当然，此定义意味着我们正在监视的过程运行正常，并且大多数结果都被认为是有效的。 例如：一家硅加工厂必须将晶圆切成相等的块。 我们知道它们每个都是0.2×0.2英寸（约0.5×0.5厘米），每侧的标准差为 0.001 英寸。 此措施是在 1,000,000 个处理步骤后确定的。 我们是否被授权将0.25×0.25英寸芯片视为异常？ 当然可以。 实际上，我们假设每边的长度都建模为高斯分布（一个非常合理的选择），其中μ = 0.2和σ = 0.001；在经过三个标准差后，概率下降到几乎为零。 因此，例如：P(edge > 0.23)≈ 0，具有这种大小的芯片必须清楚地视为异常。
显然，这是一个非常简单的示例，不需要任何模型。 但是，在现实生活中，密度的结构可能非常复杂，几个高概率区域被低概率区域包围。 这就是为什么必须采用更通用的方法来对整个样本空间进行建模的原因。
当然，异常的语义无法标准化，并且始终取决于所分析的特定问题。 因此，定义异常概念的常见方法是在异常值和新奇之间进行区分。 前者是数据集中包含的样本，即使它们与其他样本之间的距离大于平均值。 因此，离群值检测过程旨在找出此类新奇的样本（例如：考虑之前的示例，如果将0.25×0.25英寸的芯片包含在数据集中，则显然是一个离群值）。 相反，新奇检测的目标略有不同，因为在这种情况下，我们假定使用仅包含正常样本的数据集； 因此，给定一个新的芯片，我们有兴趣了解我们是否可以将其视为来自原始数据生成过程还是离群值（例如：新手技术人员向我们提出以下问题：0.25×0.25英寸的芯片是否是正常芯片）。如果我们已经收集了数据集，则可以使用我们的模型来回答问题。
描述这种情况的另一种方法是将样本视为一系列可能受可变噪声影响的值：y(t) = x(t) + n(t)。 当||n(t)|| << ||x(t)||可以分类为干净：y(t) ≈ x(t)。 相反，当||n(t)|| ≈ ||x(t)||时（或更高），它们是离群值，不能代表真实的基础过程p_data。 由于噪声的平均大小通常比信号小得多，因此P(||n(t)|| ≈ ||x(t)||)的概率接近于零。 因此，我们可以将异常想象为受异常外部噪声影响的正常样本。 异常和噪声样本管理之间真正的主要区别通常在于检测真实异常并相应地标记样本的能力。 实际上，虽然嘈杂的信号肯定已损坏，然后目标是最大程度地减少噪声的影响，但是异常现象可以很容易地被人类识别并正确标记。 但是，正如已经讨论过的，在本章中，我们有兴趣找出不依赖现有标签的发现方法。 此外，为避免混淆，我们总是引用异常，每次定义数据集的内容（仅内部数据或内部数据及异常值）以及我们的分析目标。 在下一部分中，我们将简要讨论数据集的预期结构。
数据集的结构
在标准监督（通常也是非监督）任务中，数据集有望达到平衡。 换句话说，属于每个类别的样本数量应该几乎相同。 相反，在本章要讨论的任务中，我们假设数据集X非常不平衡（包含N个样本）：
N[outlier] << N，如果存在离群检测（即数据集部分为污垢； 因此找出一种方法将所有异常值过滤掉）N[outlier] = 0（或更实际地，P(N[outlier] > 0) → 0，如果存在新颖性检测（也就是说，我们通常可以信任现有样本，而将注意力集中在新样本上）
这些标准的原因很明显：让我们考虑前面讨论的示例。 如果在 1,000,000 个处理步骤后观察到的异常率等于 0.2%，则表示存在 2,000 个异常，这对于一个工作过程而言可能是一个合理的值。 如果这个数字大得多，则意味着系统中应该存在一个更严重的问题，这超出了数据科学家的职责范围。 因此，在这种情况下，我们期望一个数据集包含大量正确的样本和非常少的异常（甚至为零）。 在许多情况下，经验法则是反映潜在的数据生成过程，因此，如果专家可以确认例如发生 0.2% 的异常，则该比率应为1000÷2来找出现实的概率密度函数。 实际上，在这种情况下，更重要的是找出确定异常值可区分性的因素。 另一方面，如果要求我们仅执行新颖性检测（例如：区分有效和恶意网络请求），则必须对数据集进行验证，以便不包含异常，但同时要进行反映负责所有可能有效样本的真实数据生成过程。
实际上，如果正确样本的数量是详尽无遗的，则与高概率区域的任何较大偏差都足以触发警报。 相反，真实数据生成过程的有限区域可能会导致假阳性结果（也就是说，尚未包含在训练集中并被错误标识为异常值的有效样本）。 在最坏的情况下，如果特征发生更改（即，错误地识别为有效样本的离群值），则噪声很大的子集也可能确定假阴性。 但是，在大多数现实生活中，最重要的因素是样本的数量和收集样本的环境。 毋庸置疑，任何模型都必须使用将要测试的相同类型的元素进行训练。 例如：如果使用低精度的仪器在化工厂内部进行测量，则高精度采集的测试可能无法代表总体（当然，它们比数据集可靠得多）。 因此，在进行分析之前，我强烈建议您仔细检查数据的性质，并询问是否所有测试样本均来自同一数据生成过程。
现在，我们可以介绍直方图的概念，这是估计包含观测值的数据集分布的最简单方法。
直方图
找出概率密度函数近似值的最简单方法是基于频率计数。 如果我们有一个包含m样本的数据集x[i] ∈ X（为简单起见，我们仅考虑单变量分布，但是过程对于多维样本完全相同），我们可以如下定义m和M：
间隔(m, M)可以分为固定数量的b个桶（它们可以具有相同或不同的宽度，表示为w(b[j])，因此n[p](b[j])对应于箱b[j]中包含的样本数。此时，给定测试样本x[t]，很容易理解，通过检测包含x[t]的桶可以很容易地获得概率的近似值 ，并使用以下公式：
在分析这种方法的利弊之前，让我们考虑一个简单的示例，该示例基于细分为 10 个不同类别的人群的年龄分布：
import numpy as np nb_samples = [1000, 800, 500, 380, 280, 150, 120, 100, 50, 30] ages = [] for n in nb_samples: i = np.random.uniform(10, 80, size=2) a = np.random.uniform(i[0], i[1], size=n).astype(np.int32) ages.append(a) ages = np.concatenate(ages)
只能使用随机种子1000（即，设置np.random.seed(1000)）来复制数据集。
ages 数组包含所有样本，我们想创建一个直方图以初步了解分布。 我们将使用 NumPy np.histrogram()函数，该函数提供所有必需的工具。 要解决的第一个问题是找出最佳箱数。 对于标准分布，这可能很容易，但是如果没有关于概率密度的先验知识，则变得非常困难。 原因很简单：因为我们需要用一个逐步函数来近似一个连续函数，所以 bin 的宽度决定了最终精度。 例如：如果密度是平坦的（例如：均匀分布），那么几个箱就足以达到良好的效果。 相反，当存在峰时，在函数的一阶导数较大时将更多（较短）的 bin 放在区域中，在导数接近零（表示平坦区域）时将较小的数目放置在区域中会很有帮助。 正如我们将要讨论的，使用更复杂的技术可以使此过程变得更容易，而直方图通常基于对最佳仓数的更粗略计算。 特别是，NumPy 允许设置bins='auto'参数，该参数将强制算法根据明确定义的统计方法（基于 Freedman Diaconis Estimator 和 Sturges 公式）自动选择数字：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-nwdpIFPC-1681652675129)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/ad9c7ed6-68a9-47f3-9993-9ed18def14f8.png)]
在上式中，四分位数范围（IQR）对应于第 75^个和第 25^个百分位数。 由于我们对分布没有一个清晰的概念，因此我们希望依靠自动选择，如以下代码片段所示：
import numpy as np h, e = np.histogram(ages, bins='auto') print('Histograms counts: {}'.format(h)) print('Bin edges: {}'.format(e))
上一个代码段的输出如下：
Histograms counts: [177 86 122 165 236 266 262 173 269 258 241 116 458 257 311 1 1 5 6] Bin edges: [16\. 18.73684211 21.47368421 24.21052632 26.94736842 29.68421053 32.42105263 35.15789474 37.89473684 40.63157895 43.36842105 46.10526316 48.84210526 51.57894737 54.31578947 57.05263158 59.78947368 62.52631579 65.26315789 68\. ]
因此，该算法定义了 19 个 bin，并已输出频率计数和边缘（即，最小值为16，最大值为68）。 现在，我们可以显示直方图的图：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-rirVoFrE-1681652675129)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/ab4c2cc0-deca-42f2-a783-a8ed3259ff51.png)]
测试分布的直方图
该图证实该分布是非常不规则的，并且一些区域的峰被平坦区域包围。 如前所述，当查询基于样本属于特定区域的概率时，直方图会很有帮助。 例如，在这种情况下，我们可能有兴趣确定某个人的年龄在 48.84 和 51.58 之间（对应于第 12 个桶）的概率从 0 开始）。 由于所有箱的宽度相同，因此我们可以简单地用n[p](b[12])（h[12]）和m（ages.shape[0]）：
d = e[1] - e[0] p50 = float(h[12]) / float(ages.shape[0]) print('P(48.84 < x < 51.58) = {:.2f} ({:.2f}%)'.format(p50, p50 * 100.0))
输出如下：
P(48.84 < x < 51.58) = 0.13 (13.43%)
因此，概率的近似值约为 13.5%，这也由直方图的结构证实。 但是，读者应该清楚地了解到这种方法有明显的局限性。 首先，也是最明显的是关于箱的数量和宽度。 实际上，一小部分产生的粗略结果无法考虑快速振荡。 另一方面，非常大的数量会产生带孔的直方图，因为大多数桶都没有样本。 因此，考虑到现实生活中可能遇到的所有动态因素，需要一种更可靠的方法。 这是我们将在下一节中讨论的内容。
核密度估计（KDE）
直方图不连续性问题的解决方案可以通过一种简单的方法有效地解决。 给定样本x[i] ∈ X，假设我们使用的是中心为x[i]多元分布，则可以考虑超体积（通常是超立方体或超球体）。 通过一个称为带宽的常数h定义了这样一个区域的扩展（已选择名称以支持该值为正的有限区域的含义）。 但是，我们现在不只是简单地计算属于超体积的样本数量，而是使用具有一些重要特征的平滑核函数K(x[i]; h)来近似估计该值：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-zAqSsgwV-1681652675130)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/623f8390-126c-4574-a155-21296102eb0d.png)]
此外，出于统计和实际原因，还必须强制执行以下整数约束（为简单起见，仅在单变量情况下显示，但扩展很简单）：
在讨论称为核密度估计（KDE）的技术之前，显示K(·)的一些常见选择将很有帮助。
高斯核
这是最常用的内核之一，其结构如下：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-LNMNB3YN-1681652675130)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/96057fde-dc28-4cce-945a-c791662917e2.png)]
以下屏幕截图显示了图形表示：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-zz6t4lY5-1681652675130)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/8f5b2bb9-40c8-42f2-a728-cdf8638c4e33.png)]
高斯核
鉴于其规律性，高斯核是许多密度估计任务的常见选择。 但是，由于该方法不允许混合不同的内核，因此选择时必须考虑所有属性。 从统计数据中，我们知道高斯分布可被视为峰度的平均参考值（峰度与峰高和尾巴的重量成正比）。 为了最大化内核的选择性，我们需要减少带宽。 这意味着即使最小的振荡也会改变密度，并且结果是非常不规则的估计。 另一方面，当h大时（即高斯的方差），近似变得非常平滑，并且可能失去捕获所有峰的能力。 因此，结合选择最合适的带宽，考虑其他可以自然简化流程的内核也会很有帮助。
Epanechnikov 核
已经提出该内核以最小化均方误差，并且它还具有非常规则的性质（实际上，可以想象为倒抛物线）。 计算公式如下：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-A4qQhiCJ-1681652675131)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/802044f2-601f-4896-b783-f68372d52094.png)]
引入常数ε可以使内核规范化并满足所有要求（以类似的方式，可以在范围内扩展内核(-h, h)，以便与其他函数更加一致）。 以下屏幕截图显示了图形表示：
Epanechnikov 核
当h → 0时，内核会变得非常尖峰。但是，由于其数学结构，它将始终保持非常规则； 因此，在大多数情况下，无需用它代替高斯核（即使后者的均方误差稍大）。 此外，由于函数在x = ±h（对于|x| > h，K(x; h) = 0）不连续，因此可能会导密集度估计值迅速下降，特别是在边界处，例如高斯函数非常缓慢地下降。
指数核
指数核是一个非常高峰的内核，其通用表达式如下：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-dEJcl60l-1681652675131)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/00116f24-e19c-42be-b5cd-868d44bcd7c3.png)]
与高斯核相反，该核的尾巴很重，峰尖尖。 以下屏幕截图显示了一个图：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-FuJ50puX-1681652675136)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/ec0ecff9-d103-47f7-ae02-f4c43a4192cf.png)]
指数核
可以看到，这样的函数适合于建模非常不规则的分布，其密度高度集中在某些特定点周围。 另一方面，当数据生成过程非常规则且表面光滑时，误差可能会非常高。 平均积分平方误差（MISE）可以用来评估内核（和带宽）的表现，是一种很好的理论方法，其定义如下：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-f7oB2vWE-1681652675136)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/2f177691-5433-440d-b548-5620efa7b7af.png)]
在上一个公式中，p[K](x)是估计的密度，而p(x)是实际的密度。 不幸的是，p(x)是未知的（否则，我们不需要任何估计）。 因此，这种方法只能用于理论评估（例如：Epanechnikov 核的最优性）。 但是，很容易理解，只要内核无法保持接近实际表面，MISE 就会更大。 由于指数突然跃升至峰值，因此仅在特定情况下才适用。 在所有其他情况下，它的行为会导致更大的 MISE，因此最好使用其他内核。
均匀（或 Tophat）内核
这是最简单且不太平滑的内核函数，其用法类似于构建直方图的标准过程。 它等于以下内容：
显然，这是一个在带宽界定的范围内恒定的步骤，仅在估计不需要平滑时才有用。
估计密度
一旦选择了核函数，就可以使用 k 最近邻方法建立概率密度函数的完全近似值。 实际上，给定数据集X（为简单起见，X ∈ R^m，所以这些值是实数），例如，通过创建球形树，很容易（如第 2 章，“聚类基础知识”中所述）以有效的方式对数据进行分区。 当数据结构准备就绪时，可以在带宽定义的半径范围内获得查询点x[j]的所有邻居。 假设这样的集合是X[j] = {x[1], ..., x[t]}，点数是N[j]。 概率密度的估计如下：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-Bae1MZgy-1681652675137)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/93a535f7-f555-4e22-b5da-bed6fa8503ba.png)]
不难证明，如果适当地选择了带宽（根据邻域中包含的样本数量而定），p[K]的概率就会收敛到实际的p(x)。 换句话说，如果粒度足够大，则近似值与真实密度之间的绝对误差将收敛为零。 下图显示了p[K](x[j])的构建过程：
x[j]的密度估计。 在属于x[j]邻域的每个点中评估内核函数
在这一点上，自然会问为什么不为每个查询使用整个数据集而不是 KNN 方法？ 答案很简单，它基于这样的假设：可以使用局部行为轻松地插值以x[j]计算的密度函数的值（即，对于多变量分布，以x[j]为中心的球和远点对估计没有影响。 因此，我们可以将计算限制为X的较小子集，避免包含接近零的贡献。
在讨论如何确定最佳带宽之前，让我们展示一下先前定义的数据集的密度估计（使用 scikit-learn）。 由于我们没有任何特定的先验知识，因此我们将使用具有不同带宽（0.1、0.5 和 1.5）的高斯核。 所有其他参数均保留为其默认值。 但是，KernelDensity类允许设置度量（默认为metric='euclidean'），数据结构（默认为algorithm='auto'，它根据维度在球树和 kd 树之间执行自动选择），以及绝对公差和相对公差（分别为 0 和10^(-8)）。 在许多情况下，无需更改默认值。 但是，对于具有特定特征的超大型数据集，例如，更改leaf_size参数以提高性能可能会有所帮助（如第 2 章，“聚类基础知识”中讨论的 ）。 此外，默认度量标准不能满足所有任务的要求（例如：标准文档显示了一个基于 Haversine 距离的示例，在使用纬度和经度时可以使用该示例）。 在其他情况下，最好使用超立方体而不是球（曼哈顿距离的情况就是这样）。
让我们首先实例化类并拟合模型：
from sklearn.neighbors import KernelDensity kd_01 = KernelDensity(kernel='gaussian', bandwidth=0.1) kd_05 = KernelDensity(kernel='gaussian', bandwidth=0.5) kd_15 = KernelDensity(kernel='gaussian', bandwidth=1.5) kd_01.fit(ages.reshape(-1, 1)) kd_05.fit(ages.reshape(-1, 1)) kd_15.fit(ages.reshape(-1, 1))
此时，可以调用score_samples()方法来获取一组数据点的对数密度估计值（在我们的示例中，我们正在考虑以 0.05 为增量的范围(10, 70)）。 由于值是log(p)，因此有必要计算exp(log(p))以获得实际概率。
生成的图显示在以下屏幕截图中：
带宽的高斯密度估计：0.1（顶部），0.5（中间）和 1.5（底部）
可能会注意到，当带宽很小（0.1）时，由于缺少特定子范围的样本，因此密度具有强烈的振荡。 当h = 0.5时，轮廓（由于数据集是单变量的）变得更加稳定，但是仍然存在一些由邻居的内部方差引起的残留快速变化。 当h变大（在我们的情况下为 1.5）时，几乎完全消除了这种行为。 一个明显的问题是：如何确定最合适的带宽？ 当然，最自然的选择是使 MISE 最小的h值，但是，正如所讨论的，只有在知道真实的概率密度时才可以使用此方法。 但是，有一些经验标准已经被证实是非常可靠的。 给定完整的数据集X ∈ R^m，第一个数据集基于以下公式：
在我们的案例中，我们获得以下信息：
import numpy as np N = float(ages.shape[0]) h = 1.06 * np.std(ages) * np.power(N, -0.2) print('h = {:.3f}'.format(h))
输出如下：
h = 2.415
因此，建议是增加带宽，甚至超过我们上一个实验中的带宽。 因此，第二种方法基于四分位数间距（IQR = Q3-Q1或等效地，第 75 个百分位减去第 25 个百分位），并且对于非常强大的内部变化，它更加健壮：
计算如下：
import numpy as np IQR = np.percentile(ages, 75) - np.percentile(ages, 25) h = 0.9 * np.min([np.std(ages), IQR / 1.34]) * np.power(N, -0.2) print('h = {:.3f}'.format(h))
现在的输出是这样的：
h = 2.051
该值比上一个值小，表明p[K](x)可以使用较小的超体积来更精确。 根据经验，我建议选择带宽最小的方法，即使第二种方法通常在不同情况下也能提供最佳结果。 现在让我们使用h = 2.0以及高斯，Epanechnikov 和指数核（我们将统一的数排除在外，因为最终结果与直方图等效）来重新执行估计：
from sklearn.neighbors import KernelDensity kd_gaussian = KernelDensity(kernel='gaussian', bandwidth=2.0) kd_epanechnikov = KernelDensity(kernel='epanechnikov', bandwidth=2.0) kd_exponential = KernelDensity(kernel='exponential', bandwidth=2.0) kd_gaussian.fit(ages.reshape(-1, 1)) kd_epanechnikov.fit(ages.reshape(-1, 1)) kd_exponential.fit(ages.reshape(-1, 1))
图形输出显示在以下屏幕截图中：
带宽等于 2.0 的密度估计，高斯核（上），Epanechnikov 核（中）和指数核（下）
不出所料，Epanechnikov 和指数核都比高斯核振荡（因为当h较小时，它们倾向于更趋于峰值）； 但是，很明显，中心图肯定是最准确的（就 MISE 而言）。 以前使用高斯核和h = 0.5时已经获得了相似的结果，但是在那种情况下，振荡极为不规则。 如所解释的， Epanechnikov 内核在值达到带宽边界时具有非常强的不连续趋势。 通过查看估计的极端现象可以立即理解该现象，该估计值几乎垂直下降到零。 相反，h = 2的高斯估计似乎非常平滑，并且无法捕获 50 到 60 年之间的变化。 指数核也发生了同样的情况，它也显示出其独特的行为：极端尖刺的极端。 在下面的示例中，我们将使用 Epanechnikov 内核； 但是，我邀请读者也检查带宽不同的高斯过滤器的结果。 这种选择有一个精确的理由（没有充分的理由就不能丢弃）：我们认为数据集是详尽无遗的，并且我们希望对克服自然极端的所有样本进行惩罚。 在所有其他情况下，可以选择非常小的残差概率。 但是，必须考虑每个特定目标做出这样的选择。
异常检测
现在，我们使用 Epanechnikov 密度估计来执行异常检测的示例。 根据概率密度的结构，我们决定在p(x) < 0.005处设置一个截止点。 以下屏幕快照中显示了这种情况：
具有异常截止的 Epanechnikov 密度估计
红点表示将样本归类为异常的年龄限制。 让我们计算一些测试点的概率密度：
import numpy as np test_data = np.array([12, 15, 18, 20, 25, 30, 40, 50, 55, 60, 65, 70, 75, 80, 85, 90]).reshape(-1, 1) test_densities_epanechnikov = np.exp(kd_epanechnikov.score_samples(test_data)) test_densities_gaussian = np.exp(kd_gaussian.score_samples(test_data)) for age, density in zip(np.squeeze(test_data), test_densities_epanechnikov): print('p(Age = {:d}) = {:.7f} ({})'.format(age, density, 'Anomaly' if density < 0.005 else 'Normal'))
上一个代码片段的输出是这样的：
p(Age = 12) = 0.0000000 (Anomaly) p(Age = 15) = 0.0049487 (Anomaly) p(Age = 18) = 0.0131965 (Normal) p(Age = 20) = 0.0078079 (Normal) p(Age = 25) = 0.0202346 (Normal) p(Age = 30) = 0.0238636 (Normal) p(Age = 40) = 0.0262830 (Normal) p(Age = 50) = 0.0396169 (Normal) p(Age = 55) = 0.0249084 (Normal) p(Age = 60) = 0.0000825 (Anomaly) p(Age = 65) = 0.0006598 (Anomaly) p(Age = 70) = 0.0000000 (Anomaly) p(Age = 75) = 0.0000000 (Anomaly) p(Age = 80) = 0.0000000 (Anomaly) p(Age = 85) = 0.0000000 (Anomaly) p(Age = 90) = 0.0000000 (Anomaly)
可以看到，函数的突然下降造成了某种垂直分离。 年龄15的人几乎处于边界（p(15) ≈ 0.0049），而行为的上限更加剧烈。 截止日期约为 58 年，但年龄60的样本比年龄 57 岁的样本低约 10 倍（这也由初始直方图证实）。 由于这只是一个教学示例，因此很容易检测到异常。 但是，如果没有标准化的算法，即使是稍微更复杂的分布也会产生一些问题。 特别地，在这种简单的单变量分布的特定情况下，异常通常位于尾部。
因此，我们假设给定整体密度估计p[K](x)：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-KoLmz3cu-1681652675138)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/749c0092-8f96-4900-8bb3-062de85e7b6a.png)]
当考虑包含所有样本（正常样本和异常样本）的数据集时，这种行为通常是不正确的，并且数据科学家在确定阈值时必须小心。 即使很明显，也最好通过从数据集中删除所有异常来学习正态分布，以展开异常所在的区域（p[K](x) → 0）。 这样，先前的标准仍然有效，并且可以轻松比较不同的密度以进行区分。
在继续下一个示例之前，我建议通过创建人工漏洞并设置不同的检测阈值来修改初始分布。 此外，我邀请读者根据年龄和身高生成双变量分布（例如：基于一些高斯的总和），并创建一个简单的模型，该模型能够检测所有参数不太可能出现的人。
将 KDD Cup 99 数据集用于异常检测
本示例基于 KDD Cup 99 数据集，该数据集收集了一系列正常和恶意的互联网活动。 特别是，我们将重点放在 HTTP 请求的子集上，该子集具有四个属性：持续时间，源字节，目标字节和行为（这是一个分类元素，但是对我们而言，可以立即访问某些特定的属性很有帮助。 攻击）。 由于原始值是非常小的零附近的数字，因此所有版本（包括 scikit-learn 在内）都使用公式log(x + 0.1)（因此，在用新样本模拟异常检测时必须使用它）。 当然，逆变换如下：
让我们首先使用 scikit-learn 内置函数fetch_kddcup99()加载并准备数据集，然后选择percent10=True将数据限制为原始集合的 10% （非常大）。 当然，我邀请读者也使用整个数据集和完整的参数列表（包含 34 个数值）进行测试。
在这种情况下，我们还选择subset='http'，它已经准备好包含大量的正常连接和一些特定的攻击（如在标准期刊日志中）：
from sklearn.datasets import fetch_kddcup99 kddcup99 = fetch_kddcup99(subset='http', percent10=True, random_state=1000) X = kddcup99['data'].astype(np.float64) Y = kddcup99['target'] print('Statuses: {}'.format(np.unique(Y))) print('Normal samples: {}'.format(X[Y == b'normal.'].shape[0])) print('Anomalies: {}'.format(X[Y != b'normal.'].shape[0]))
输出如下：
Statuses: [b'back.' b'ipsweep.' b'normal.' b'phf.' b'satan.'] Normal samples: 56516 Anomalies: 2209
因此，使用2209恶意样本和56516正常连接有四种类型的攻击（在此情况下，其详细信息不重要）。 为了进行密度估计，为了进行一些初步考虑，我们将把这三个分量视为独立的随机变量（虽然不完全正确，但是可以作为一个合理的起点），但是最终估计是基于完整的联合分布 。 当我们要确定最佳带宽时，让我们执行基本的统计分析：
import numpy as np means = np.mean(X, axis=0) stds = np.std(X, axis=0) IQRs = np.percentile(X, 75, axis=0) - np.percentile(X, 25, axis=0)
上一个代码段的输出如下：
Means: [-2.26381954 5.73573107 7.53879208] Standard devations: [0.49261436 1.06024947 1.32979463] IQRs: [0\. 0.34871118 1.99673381]
持续时间的 IQR（第一个部分）为空； 因此，大多数值是相等的。 让我们绘制一个直方图来确认这一点：
第一部分的直方图（持续时间）
不出所料，这种成分不是很重要，因为只有一小部分样本具有不同的值。 因此，在此示例中，我们将跳过它，仅使用源字节和目标字节。 现在，如前所述，计算带宽：
import numpy as np N = float(X.shape[0]) h0 = 0.9 * np.min([stds[0], IQRs[0] / 1.34]) * np.power(N, -0.2) h1 = 0.9 * np.min([stds[1], IQRs[1] / 1.34]) * np.power(N, -0.2) h2 = 0.9 * np.min([stds[2], IQRs[2] / 1.34]) * np.power(N, -0.2) print('h0 = {:.3f}, h1 = {:.3f}, h2 = {:.3f}'.format(h0, h1, h2))
输出如下：
h0 = 0.000, h1 = 0.026, h2 = 0.133
除了第一个值，我们需要在h1和h2之间进行选择。 由于值的大小不大并且我们希望具有较高的选择性，因此我们将设置h = 0.025，并使用高斯核，该核提供了良好的平滑度。 下面的屏幕快照显示了分割输出（使用包含一个内部 KDE 模块的 seaborn 可视化库获得），其中还包含第一个组件。
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-129AhkA6-1681652675139)(https://gitcode.net/apachecn/apachecn-dl-zh/-/raw/master/docs/handson-unsup-learn-py/img/99fdb6a4-0ddc-48d4-8396-53b95f4be661.png)]
正常连接（上面一行）和恶意攻击（下面一行）的密度估计
第一行显示了正常连接的密度，而下一行是恶意攻击。 正如预期的那样，两种情况下的第一部分（持续时间）几乎相同，可以将其丢弃。 相反，源字节和目标字节都表现出非常不同的行为。 在不考虑对数变换的情况下，普通连接平均发送 5 个字节，其方差很小，从而将电位范围扩展到间隔（ 4 ， 6 ） 。 响应具有较大的方差，其值在 4 和 10 之间，并且从 10 开始具有非常低的密度。 相反，恶意攻击的源字节和目标字节都具有两个峰值：一个较短的峰值对应于 -2 ，一个较高的峰值分别对应于大约 11 和。 9 （与正常区域的重叠最小）。 即使不考虑全部联合概率密度，也不难理解大多数攻击会发送更多的输入数据并获得更长的响应（而连接持续时间并没有受到很大影响）。
现在，我们可以通过仅选择正常样本（即，对应于Y == b'normal.'）来训练估计器：
from sklearn.neighbors import KernelDensity X = X[:, 1:] kd = KernelDensity(kernel='gaussian', bandwidth=0.025) kd.fit(X[Y == b'normal.'])
让我们计算正常样本和异常样本的密度：
Yn = np.exp(kd.score_samples(X[Y == b'normal.'])) Ya = np.exp(kd.score_samples(X[Y != b'normal.'])) print('Mean normal: {:.5f} - Std: {:.5f}'.format(np.mean(Yn), np.std(Yn))) print('Mean anomalies: {:.5f} - Std: {:.5f}'.format(np.mean(Ya), np.std(Ya)))
输出如下：
Mean normal: 0.39588 - Std: 0.25755 Mean anomalies: 0.00008 - Std: 0.00374
显然，当例如p[K](x) < 0.05（考虑三个标准差），我们得到p时，我们可以预期到异常。 [K] (x) ∈ (0, 0.01)，而Yn的中位数约为 0.35。 这意味着至少一半的样本具有p[K](x) > 0.35。 但是，通过简单的计数检查，我们得到以下信息：
print(np.sum(Yn < 0.05)) print(np.sum(Yn < 0.03)) print(np.sum(Yn < 0.02)) print(np.sum(Yn < 0.015))
输出如下：
3147 1778 1037 702
Python 无监督学习实用指南：6~10（2）https://developer.aliyun.com/article/1426887
