网站导航
> 技术文档 > Python 时间序列的机器学习（二）

Python 时间序列的机器学习（二）

网友: 技术文档 2025-09-03 00:59:19

原文：annas-archive.org/md5/6fe7868ac7b4b689aa7ceab9954a4809

译者：飞龙

协议：CC BY-NC-SA 4.0

第五章：使用移动平均和自回归模型进行预测

本章介绍基于移动平均和自回归的时间序列建模。这个主题包含了一大类模型，在不同学科中都非常流行，包括计量经济学和统计学。我们将讨论自回归和移动平均模型，以及一些将这两者结合起来的模型，如 ARMA、ARIMA、VAR、GARCH 等。

这些模型仍然受到高度评价并且有广泛的应用。然而，自那时以来，许多新模型涌现出来，证明它们与这些简单模型具有竞争力，甚至能够超越它们。然而，在其主要应用领域，即单变量预测中，简单模型往往能够提供准确或足够准确的预测，因此这些模型在时间序列建模中仍然占据主导地位。

我们将涵盖以下主题：

  • 经典模型是什么？

    • 移动平均和自回归

    • 模型选择和阶数

    • 指数平滑法

    • ARCH 和 GARCH

    • 向量自回归

  • Python 库

    • statsmodels

  • Python 实践

    • 在 Python 中建模

我们将从经典模型的介绍开始。

经典模型是什么？

在本章中，我们将讨论那些有着更长传统的模型，这些模型根植于统计学和数学领域。它们在计量经济学和统计学中有着广泛应用。

虽然统计学和机器学习方法之间有相当多的重叠，而且两个领域的学者都在吸收对方的研究成果，但仍然存在一些关键差异。统计学论文仍然主要是正式和推理性的，而机器学习研究者则更加务实，依赖于模型的预测准确性。

我们在 第一章《Python 时间序列介绍》中讨论了时间序列模型的早期历史。在本章中，我们将讨论用于预测的移动平均和自回归方法。这些方法在 20 世纪初期被提出，并在 1970 年由 George Box 和 Gwilym Jenkins 在他们的著作《时间序列分析：预测与控制》中推广。关键的是，在这本书中，Box 和 Jenkins 使 ARIMA 模型得到了形式化，并描述了如何将其应用于时间序列预测。

许多时间序列展现出趋势和季节性，而本章中的许多模型假设时间序列是平稳的。如果一个时间序列是平稳的，它的均值和标准差会随着时间保持不变。这意味着该时间序列没有趋势，也没有周期性波动。

因此，去除不规则成分、趋势和季节波动是应用这些模型的一个内在方面。模型然后预测去除季节性和趋势后的残余部分：商业周期。

因此，应用经典模型时，时间序列通常需要分解为不同的组成部分。因此，经典模型通常按以下方式应用：

  1. 平稳性检验

  2. 差分处理【如果检测到平稳性】

  3. 拟合方法和预测

  4. 加回趋势和季节性

本章中的大多数方法仅适用于单变量时间序列。尽管已经提出了扩展到多变量时间序列的方法，但它们不像单变量版本那样流行。单变量时间序列由单个向量组成，换句话说，就是一个随时间变化的值。尽管如此，我们将在本章末尾看到向量自回归（VAR），这是一种多变量时间序列的扩展。

另一个重要的考虑因素是，大多数经典模型都是线性的，这意味着它们假设时间点之间以及不同时间步长之间的值之间存在线性依赖关系。实际上，本章中的模型在处理一系列平稳时间序列时效果良好。平稳性意味着分布随时间保持不变。一个例子就是温度随时间的变化。这些模型在数据量较小的情况下尤其有价值，因为在这种情况下，非线性模型中的额外估计误差会超过精度上的潜在增益。

然而，平稳性假设意味着本章中模型的应用仅限于具有这一属性的时间序列。否则，我们就需要对时间序列进行预处理，以强制其平稳性。相比之下，非线性时间序列分析与预测的统计方法发展较少受到关注；然而，仍然存在一些模型，例如阈值自回归模型（我们在这里不讨论）。

最后，需要指出的是，虽然这是一个合理的初步方法，但许多时间序列（如温度）的预测通过基于物理的大气高维模型，比通过统计模型更为准确。这说明了复杂性的问题：本质上，建模是将一组假设凝练并与参数一起形式化的过程。

现实世界的时间序列来自复杂的过程，这些过程可能是非线性和非平稳的，描述它们的方式有很多种，每种方法都有其优缺点。因此，我们可以从很多参数的角度看待建模问题，或者仅将其视为单一或几个参数。在下面的专门部分，我们将讨论如何根据参数数量和精度，从一组备选模型中选择一个模型的问题。

目前，非线性模型来源于不同的研究方向，主要是神经网络或更广泛的机器学习领域。我们将在第十章中讨论神经网络，《时间序列的深度学习》，并将在第七章中讨论并应用最先进的机器学习方法，《时间序列的机器学习模型》。

移动平均和自回归

经典模型可以分为几类模型——移动平均（MA）、自回归（AR）模型、ARMA 和 ARIMA。这些模型经过时间的推敲，由许多数学家和统计学家在书籍和论文中正式化并普及开来，包括彼得·惠特尔（1951 年）和乔治·博克斯与吉威尔姆·詹金斯（1970 年）。不过，让我们从更早的时候开始。

移动平均标志着现代时间序列预测的开始。在移动平均中，通常会对过去一段时间内（时间框架）一定数量的时间点的值取平均（通常是算术平均）。

更正式地说，简单移动平均是一个无权重的平均值，计算范围为 k 个点，公式为：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_001.png

其中 x[i] 表示观察到的时间序列。

移动平均可以用来平滑时间序列，从而去除短期内发生的噪声和周期波动，实际上起到低通滤波器的作用。因此，正如数学家雷金纳德·胡克在 1902 年的一篇出版物中指出的那样，移动平均可以用来分离趋势和振荡成分。他将趋势概念化为忽略振荡后，序列前进的方向。

移动平均可以平滑时间序列中的趋势和周期；然而，作为一种模型，移动平均也可以用于预测未来。时间序列是当前序列值与观察值（误差项）之间的线性回归。移动平均模型的阶数为 q，即 MA(q)，可以表示为：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_002.png

其中 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_003.png 是 x[t] 的平均值（期望值）（通常假设为 0），https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_004.png 是参数，https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_005.png 是随机噪声。

胡克在剑桥大学接受教育，曾在英国农业、渔业和食品部的统计部门工作。他是一个业余统计学家，撰写关于气象学和社会经济话题的文章，如工资、婚姻率、贸易以及作物预测等。

AR 技术的发明可以追溯到 1927 年，英国统计学家乌德尼·尤尔（Udny Yule）的论文（《On a Method of Investigating Periodicities in Disturbed Time-Series with special reference to Wolfer’s Sunspot Numbers》），尤尔是胡克的朋友。自回归模型（autoregressive model）将变量与其自身的滞后值回归。换句话说，当前值由紧接其后的值通过线性组合驱动。

太阳黑子的变化具有高度周期性，正如这张显示太阳黑子随时间变化的图表所示（通过 statsmodels 数据工具加载）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_01.png

图 5.1：按年划分的太阳黑子观测数据

尤尔提出了一个由噪声驱动的线性模型，用来应用于太阳黑子数目，即太阳外壳上的黑斑数量。这些黑斑源自巨大的爆炸，标志着太阳的磁活动，并且与日冕物质抛射等现象相关。

下面是根据太阳黑子数目展示的低太阳活动和高太阳活动的两幅图像（来自 NASA）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_02.png

图 5.2：太阳活动

今天，我们知道太阳周期是太阳磁活动的一个几乎周期性的 11 年变化，表现为高磁活动（太阳极大期）和低磁活动（太阳极小期）之间的变化。在高点时，爆炸（太阳耀斑）会将带电粒子释放到太空中，可能危及地球上的生命。

尤尔（Yule）在伦敦大学学院（UCL）学习工程学，后来与海因里希·赫兹（Heinrich Hertz）一起在波恩工作，随后回到 UCL 与卡尔·皮尔逊（Karl Pearson）合作，并最终晋升为 UCL 的助理教授。在 UCL 担任统计学职位后，他移居剑桥。他因其关于统计学的著作《统计学理论导论》而被人们铭记，这本书首次出版于 1911 年，经过了多个版本的修订。此外，他还因描述了现在被称为优先附着过程的现象而闻名，该过程描述了网络中新节点的分配方式是根据节点已经拥有的数量来决定的；这一过程有时被称为“富者愈富”。

安德烈·柯尔莫哥洛夫（Andrey Kolmogorov）在 1931 年定义了平稳过程这一术语，尽管路易·巴谢列（Louis Bachelier）在早些时候（1900 年）就用不同的术语提出了类似的定义。平稳性由三个特征定义：

  1. 有限变差

  2. 恒定均值

  3. 恒定变差

恒定变差意味着时间序列在两个点之间的窗口内的变化随着时间的推移保持恒定：https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_006.png，尽管它可能随窗口的大小变化。

这是弱平稳性。在文献中，除非另有说明，通常平稳性指的是弱平稳性。严格平稳性意味着时间序列具有随时间不变的概率密度函数。换句话说，在严格平稳性下，联合分布在 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_007.png 上与在 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_008.png 上是相同的。

1938 年，挪威数学家赫尔曼·奥勒·安德烈亚斯·沃尔德（Herman Ole Andreas Wold）描述了平稳时间序列的分解。他观察到，平稳时间序列可以表示为一个确定性成分（自回归）与一个随机成分（噪声）的和。如今，这一分解方法以他的名字命名，称为沃尔德分解。

这导致了自回归模型的形式化，阶数为 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_009.png，AR§，表示为：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_010.png

其中 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_011.png 是模型参数，c 是一个常数，而 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_012.png 则代表噪声。在这个方程中，p 是时间序列连续值之间自相关性的度量。

此工作后来在 1951 年被新西兰人彼得·惠特尔的博士论文《时间序列假设检验》中推广到多变量时间序列中，彼得·惠特尔的导师是沃尔德。彼得·惠特尔也因将 AR 和 MA 模型整合为一个而受到赞誉，这就是自回归移动平均（ARMA）的另一个里程碑，将尤尔和胡克的工作结合在一起。

ARMA 模型包括两种类型的滞后值，一种用于自回归分量，另一种用于移动平均分量。因此，我们写成ARMA(p, q)，其中第一个参数 p 表示自回归的顺序，第二个q表示移动平均的顺序，如下所示：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_013.png

ARMA 假设该系列是稳定的。在实践中，为了确保稳定性，必须应用预处理。

模型参数是通过最小二乘法估计的，直到乔治·Box 和 Gwilym Jenkins 推广了他们的最大似然估计参数方法。

乔治·Box 是不仅在古典时间序列预测领域中最有影响力的人物之一，也在更广泛的统计领域中卓有成效。在二战期间，他未完成化学学习就被征召入伍，为军队进行毒气实验，在此过程中自学了统计学以进行分析。

战争结束后，他在伦敦大学学院学习数学和统计学，并在埃冈·皮尔逊（Karl Pearson 之子）的指导下完成了他的博士学位。后来，他在普林斯顿大学领导了一个研究团队，然后在威斯康星大学麦迪逊分校创立了统计学系。

Box 和 Jenkin 的 1970 年著作 “Time-Series Analysis: Forecasting and Control”，详细介绍了时间序列预测和季节调整的许多应用示例。所谓的 Box-Jenkins 方法是最流行的预测方法之一。他们的书籍还描述了自回归积分移动平均模型（ARIMA）。

ARIMA(p, d, q)包括一个数据预处理步骤，称为积分，使时间序列保持稳定，通过替换值来减去即时过去的值，这个转换称为差分。

模型的积分由参数 d 参数化，它是当前值和先前值之间进行差分的次数。正如提到的那样，这三个参数代表模型的三个部分。

有一些特殊情况；ARIMA(p,0,0)代表 AR§，ARIMA(0,d,0)代表 I(d)，而 ARIMA(0,0,q)则是 MA(q)。I(0)有时被用作一个约定，指的是不需要任何差分即可保持稳定的时间序列。

虽然 ARIMA 类型模型有效地考虑了平稳过程，但作为 ARMA 模型扩展的 季节性自回归积分滑动平均 模型（SARIMA）可以描述在季节内外展现非平稳行为的过程。

季节性 ARIMA 模型通常表示为 ARIMA(p,d,q)(P,D,Q)m。各个参数需要进一步解释：

  • m 表示一个季节中的周期数

  • P、D、Q 参数化季节部分的自回归、差分和移动平均组件

  • p、d、q 是我们之前讨论过的 ARIMA 项。

P 是度量时间序列中连续季节组件之间自相关性的指标。

我们可以列出季节部分来使其更清楚。季节性自回归，SAR，可以表示为：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_014.png

其中 s 是季节性周期的长度。

类似地，季节性移动平均，SMA，可以写作如下：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_015.png

请注意，这些组件中的每一个将使用一组不同的参数。

例如，模型 SARIMA(0,1,1)(0,1,1)12 过程将包含一个非季节性 MA(1) 项（对应参数 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_016.png），以及一个季节性 MA(1) 项（对应参数 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_017.png）。

模型选择与顺序

ARMA 中的参数 q 通常为 3 或更少，但这更多地反映了计算资源的限制，而非统计问题。如今，为了设置 p 和 q 参数，我们通常会查看自相关和偏自相关图，其中我们可以看到每个滞后的相关性峰值。

当我们有不同的模型，比如不同的 p 和 q 模型，每个模型都在相同的数据集上进行训练时，我们怎么知道应该使用哪个模型呢？这就是模型选择的作用。

模型选择是决定竞争模型之间选择的方法论。模型选择中的一个主要思想是奥卡姆剃刀原则，得名于英国方济各会修士、经院哲学家威廉·奥卡姆（William of Ockham），他生活在公元 1287 年至 1347 年间。

根据奥卡姆剃刀原则，在选择竞争性解决方案时，应优先选择假设最少的解释。奥卡姆基于这一思想主张，神圣干预的原则是如此简单，奇迹便是一个简洁的解释。这个规则，也叫做拉丁语中的 “lex parsimoniae”，表明一个模型应当简洁，即使它简单，却应具有较高的解释力。

在科学中，出于可证伪性原则，通常偏好更简单的解释。科学解释越简单，越容易进行测试，甚至可能被推翻——这为模型提供了科学的严谨性。

ARMA 和其他模型通常通过 最大似然估计（MLE）进行估计。在 MLE 中，这意味着最大化似然函数，以便在给定模型参数的情况下，观察到的数据是最有可能的。

最大似然法中最常用的模型选择标准之一是 赤池信息准则 (AIC)，得名于 Hirotugu Akaike，他在 1973 年首次用英语发布了这一准则。

AIC 采用最大似然法中的对数似然值 l 和模型中的参数个数 k。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_018.png

这表示 AIC 等于参数个数的两倍减去对数似然的两倍。在模型选择中，我们会选择 AIC 最小的模型，这意味着它的参数少，但对数似然高。

对于 ARIMA 模型，我们可以更具体地写出：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_019.png

我已经省略了参数 d，因为它不会引入额外的估计。

贝叶斯信息准则 (BIC) 是由 Gideon Schwarz 在几年前（1978 年）提出的用于模型选择的方法，和 AIC 非常相似。它额外考虑了 N，数据集中的样本数量：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_020.png

根据 BIC，我们希望选择一个具有较少参数和较高对数似然的模型，但同时样本量也要小。

指数平滑

指数平滑，追溯至 Siméon Poisson 的工作，是一种通过指数窗口函数平滑时间序列数据的技术，可以用于预测具有季节性和趋势的时间序列。

最简单的方法，简单指数平滑，SES，s[t] 的时间序列 x[t] 可以表示为：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_021.pnghttps://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_022.png

其中 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_023.png 是指数平滑因子（一个介于 0 和 1 之间的值）。

本质上，这是一个加权移动平均，权重为 α 和 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_024.png。你可以把第二项 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_025.png 看作递归的，在展开时， https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_026.png 会被重复相乘——这就是指数项。

参数 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_027.png 控制平滑值在当前值和前值之间的决定权重。与移动平均法相似，这个公式的效果是结果变得更加平滑。

有趣的是，John Muth 在 1960 年展示了 SES 在时间序列中提供了最佳预测，其中在每个时间步，值会随机偏离前一个值，而且步长是独立同分布的，加上噪声。这种时间序列称为随机游走，有时，波动的股票价格假定遵循这种行为。

另一种指数平滑方法是Theta 方法，它对从业者尤其具有吸引力，因为它在 2000 年的 M3 竞赛中表现出色。M3 竞赛得名于其组织者斯皮罗斯·马克里达基斯（Spyros Makridakis），他是尼科西亚大学的教授兼未来研究所的主任。该竞赛涉及 3003 个来自微观经济学、工业、金融、人口学等地方的时间序列。其主要结论之一是，非常简单的方法也可以在单变量时间序列预测中表现良好。M3 竞赛被证明是预测领域的一个转折点，提供了基准和最先进的技术（SOTA），尽管自那时以来，SOTA 已经发生了显著变化，正如我们在第七章《时间序列的机器学习模型》中将看到的那样。

Theta 方法由瓦西里斯·阿西马科普洛斯（Vassilis Assimakopoulos）和科斯塔斯·尼科洛普洛斯（Konstantinos Nikolopoulos）于 2000 年提出，并由罗布·海德曼（Rob Hyndman）和巴基·比拉（Baki Billah）在 2001 年重新阐述。Theta 模型可以理解为简单指数平滑（SES）带有漂移项。

该方法基于将去季节化数据分解为两条线。第一条所谓的“theta”线估计长期成分——趋势，然后将此趋势与 SES 的加权平均值相结合。

让我们更正式地阐述这一点！趋势成分的预测公式如下：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_028.png

在这个方程中，c是截距，https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_029.png是乘以时间步长的系数，https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_030.png是残差。https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_031.png可以通过普通最小二乘法拟合。

Theta 的公式是将这一趋势与 SES 加权求和：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_032.png

这里，https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_033.png是时间步长t时刻的X预测值。

最流行的指数平滑方法之一是霍尔茨-温特斯方法。芝加哥大学的查尔斯·霍尔茨教授于 1957 年首次发布了一种双重指数平滑方法，该方法允许基于趋势和水平进行预测。他的学生彼得·温特斯在 1960 年扩展了该方法，以捕捉季节性（“通过指数加权移动平均法预测销售额”）。即便在后期，霍尔茨-温特斯平滑法也被进一步扩展，以考虑多重季节性（n 阶平滑）。

为了应用霍尔茨-温特斯方法，我们首先去除趋势和季节性。然后我们预测时间序列并将季节性和趋势加回来。

我们可以区分方法的加性和乘法变体。趋势和季节性都可以是加性或乘法的。

加性季节性是独立于系列值添加的季节性。乘法季节性成分是按比例添加的，当季节效应随着时间序列中的值（或趋势）的增加或减少而变化时。通过视觉检查可以帮助决定使用哪种变体。

Holtz-Winters 方法也称为三重指数平滑法，因为它应用了三次指数平滑，正如我们所见。Holtz-Winters 方法捕捉了三个组成部分：

  • 对每个时间点的水平估计，L[t] —— 这可能是一个平均值

  • 趋势分量 T

  • 季节性 S[t]，具有 m 个季节（即一年中的季节数）

对于加性趋势和季节性，在数学上，Holtz-Winters 预测值 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_034.png 定义为：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_035.png

对于 乘法季节性，我们将其与季节性相乘：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_036.png

水平更新如下：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_037.png

我们基于两个项的加权平均来更新当前水平，其中 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_038.png 为两者之间的权重。这两个项分别是前一水平和去季节化后的系列值。

在这个方程中，我们通过除以季节性来进行去季节化：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_039.png

之前的趋势分量会像这样被加到前一个水平：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_040.png

趋势更新如下（对于加性趋势）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_041.png

最终，（乘法）季节性更新如下：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_042.png

我们可以根据需要将这些方程切换为加性变体。本书的范围之外会有更详细的讨论——我们就此止步。

ARCH 和 GARCH

麻省理工学院经济学教授 Robert F. Engle 提出了一个时间序列预测模型（1982），他将其命名为 ARCH（自回归条件异方差）。

对于金融机构而言，风险价值（value at risk），即在特定时间段内的金融风险水平，是风险管理中的一个重要概念。因此，考虑资产收益的协方差结构至关重要。这就是 ARCH 所做的，也解释了其重要性。

事实上，鉴于他在时间序列计量经济学领域的贡献，Engle 与之前提到的 Clive Granger 一起，于 2003 年获得了诺贝尔经济学奖（诺贝尔经济科学纪念奖）。表彰中特别提到了他在 ARCH 方面的开创性工作。

在 ARMA 型模型中，收益被建模为独立且在时间上相同分布，而 ARCH 通过对收益在不同频率下的高阶依赖性进行参数化，从而允许时间变化的（异方差）误差项。

在 ARCH 中，残差表示为由一个随机项 z[t] 和一个标准差 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_043.png 组成，这两者都是时间相关的： https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_044.png。

时间 t 处残差的标准差模型依赖于前期点的残差：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_045.png

其中 q 是方差依赖的前期时间点的数量。

模型 ARCH(q) 可以通过最小二乘法来确定。

最小二乘法算法是求解线性方程 y=X.β 以获得 β。它的原理是找到最小化误差平方和的参数，即 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_046.png。

GARCH（广义自回归条件异方差模型）是在 Tim Bollerslev（1986）和 Stephen Taylor（1986）分别扩展 Engle 模型使其更具普遍性时诞生的。GARCH 和 ARCH 之间的主要区别在于，残差来自 ARCH 模型，而不是自回归模型 AR。

通常，在应用 GARCH 或 ARCH 模型之前，会进行同方差性的统计检验，换句话说，就是检验方差是否随时间保持恒定。常用的检验方法是 ARCH-LM 检验，零假设为时间序列没有 ARCH 效应。

向量自回归

本章介绍的所有预测方法都针对单变量时间序列，即由单一时间依赖变量组成的时间序列，一个单一向量。在实践中，我们通常知道比单一测量序列更多的信息。

例如，如果我们的时间序列是关于冰淇淋销量的，我们可能还知道温度或泳衣的销量。我们可以预期冰淇淋销量与温度高度相关，实际上，当气温较高时，冰淇淋的消费可能会增加。同样，我们也可以推测泳衣的销量与冰淇淋的销量要么同时发生，要么先于冰淇淋销量，或晚于冰淇淋销量。

向量自回归模型可以追踪多个变量随时间变化的关系。它们可以捕捉时间序列与当前时间戳之前的数值向量之间的线性依赖，将 AR 模型推广到多变量时间序列。

VAR 模型的特点是其阶数，即模型中考虑的前几个时间点的数量。最简单的情况是 VAR(1)，其中模型只考虑时间序列的一个滞后项，公式如下：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_047.png

其中 c 是常数，即直线的截距，https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_048.png 是模型的系数，而 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_049.png 是 t 时刻的误差项。x[t] 和 c 是长度为 k 的向量，而 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_050.png 是一个 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_051.png 矩阵。

p 阶模型，VAR§，表示为：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_052.png

VAR 假设误差项的均值为 0，且误差项之间没有序列相关性。

就像向量自回归是自回归的多变量推广，向量 ARIMA（VARIMA）是对单变量 ARIMA 模型的扩展，用于多变量时间序列。尽管早在 1957 年就已经正式提出，但可用的软件实现直到后来才出现。

在下一节中，我们将介绍一些可以在 Python 中用于经典模型预测的库。

Python 库

Python 中有一些流行的经典时间序列建模库，但其中最流行的无疑是 statsmodels。以下图表比较了各库在 GitHub 上的星标数，显示它们的受欢迎程度：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_03.png

图 5.3：Python 经典时间序列预测库的受欢迎程度

Statsmodels 显然是这些库中最受欢迎的。我只选择了那些积极维护并直接实现算法的库，而不是从其他库导入的库。例如，SkTime 或 Darts 库提供传统的预测模型，但这些模型并未在它们那里实现，而是在 statsmodels 中实现的。

pmdarima（最初是 pyramid-arima）包含一个参数搜索功能，帮助拟合最佳 ARIMA 模型到单变量时间序列。Anticipy 包含一些模型，如指数衰减模型和阶梯模型。Arch 实现了金融计量经济学工具和 自回归条件异方差性（ARCH）的功能。

尽管 statsmodels 的活跃度不如 Scikit-Learn，且仅由少数人维护，但它仍是时间序列传统统计和计量经济学方法的首选库，尤其在参数估计和统计检验上，强调程度远高于机器学习。

Statsmodels

statsmodels 库可以帮助估计统计模型并进行统计检验。它基于 SciPy 和 NumPy，包含许多统计函数和模型。

以下表格展示了与本章相关的一些建模类：

类别 描述 ar_model.AutoReg 单变量自回归模型 arima.model.ARIMA 自回归积分滑动平均（ARIMA）模型 ExponentialSmoothing Holt Winter’s 指数平滑 SimpleExpSmoothing 简单指数平滑

图 5.4：statsmodels 中实现的几个模型

ARIMA 类还通过 seasonal_order 参数支持 SARIMA，即带有季节性组件的 ARIMA。根据定义，ARIMA 还支持 MA、AR 和差分（集成）。

还有一些其他模型，如 Markov 自回归模型，但我们不会逐一讲解，我们将选择性地讲解一些。

这里列出了其他一些有用的函数：

函数 描述 stattools.kpss Kwiatkowski-Phillips-Schmidt-Shin 平稳性检验 stattools.adfuller 扩展的 Dickey-Fuller 单位根检验 stattools.ccf 互相关函数 stattools.pacf 偏自相关估计 stats.diagnostic.het_arch Engle 的自回归条件异方差性（ARCH）检验，也称为 ARCH-LM 检验 stattools.q_stat Ljung-Box Q 统计量 tsa.seasonal.seasonal_decompose 使用移动平均法进行季节性分解 tsa.tsatools.detrend 去趋势一个向量

图 5.5：statsmodels 中的有用函数

按照惯例，我们像这样导入 statsmodels：

import statsmodels.api as sm

这些 statsmodels 算法也可以通过 SkTime 使用，它通过类似 Sklearn 接口的方式提供访问。

这应该足够提供一个简短的概览。接下来我们进入建模部分！

Python 实践

正如本章引言所述，我们将使用 statsmodels 库进行建模。

需求

在本章中，我们将使用几个库，这些库可以从终端（或类似地从 Anaconda Navigator）快速安装：

pip install statsmodels pandas_datareader

我们将从 Python（或 IPython）终端执行这些命令，当然，我们也可以从 Jupyter notebook（或其他环境）执行。

让我们开始建模吧！

Python 中的建模

我们将使用 Yahoo 财经的股票代码数据集，通过 yfinance 库下载。我们首先加载数据集，进行一些快速探索，然后构建本章中提到的几个模型。

我们将加载一系列标准普尔存托凭证（SPDR S&P 500 ETF 信托基金）：

from datetime import datetimeimport yfinance as yfstart_date = datetime(2005, 1, 1)end_date = datetime(2021, 1, 1)df = yf.download( \'SPY\', start=start_date, end = end_date)

我们需要指定日期范围和股票代码。每日价格包括开盘价、收盘价等。我们将使用开盘价进行分析。

索引列已经是 pandas 的 DateTimeIndex，因此我们无需进行转换。现在，让我们绘制这个时间序列图！

import matplotlib.pyplot as pltplt.title(\'Opening Prices between {} and {}\'.format( start_date.date().isoformat(), end_date.date().isoformat()))df[\'Open\'].plot()plt.ylabel(\'Price\')plt.xlabel(\'Date\');

这将生成以下图表：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_04.png

图 5.6：标准普尔存托凭证价格随时间变化

由于这是每日数据，且每年有 253 或 252 个工作日，我决定将数据重采样为每周数据，使每年数据保持一致。

df1 = df.reset_index().resample(\'W\', on=\"Date\")[\'Open\'].mean()df1 = df1[df1.index.week < 53]

有些年份有 53 周。我们无法处理这种情况，因此将去掉第 53 周。现在我们有跨越 16 年的 52 周的每周数据。

最后一个修正：statsmodels 可以使用与 DateTimeIndex 相关联的频率信息；然而，这通常没有设置，df1.index.freq 为 None。所以，我们将自己设置：

df1 = df1.asfreq(\'W\').fillna(method=\'ffill\')

如果我们现在检查，df1.index.freq 为 。

设置频率可能会导致缺失值。因此，我们使用 fillna() 操作将缺失值替换为最后一个有效值。如果不这样做，某些模型将无法收敛，并且会返回 NaN（非数字）值，而不是预测值。

现在我们需要了解模型阶数的合理范围。我们将查看自相关和部分自相关函数来帮助确定：

import statsmodels.api as smfig, axs = plt.subplots(2)fig.tight_layout()sm.graphics.tsa.plot_pacf(df1, lags=20, ax=axs[0])sm.graphics.tsa.plot_acf(df1, lags=20, ax=axs[1])

这将生成以下图表：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_05.png

图 5.7：部分自相关与自相关

这些图表显示了时间序列在最多 20 个时间步的滞后期上的自相关性。R 或 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_053.png 值接近 0 表示滞后期的连续观察值之间没有相关性。相反，接近 1 或 -1 的相关性表示这些滞后期的观察值之间存在强烈的正相关或负相关。

自相关和部分自相关都返回置信区间。如果自相关值超出了置信区间（表示为阴影区域），则表明相关性显著。

我们可以看到滞后期为 1 的部分自相关非常高，而较高滞后的自相关则较低。所有滞后的自相关都显著且较高，但随着滞后期的增加，显著性逐渐降低。

让我们继续讨论自回归模型。从这里开始，我们将使用 statsmodels 的建模功能。这个接口非常方便，正如你将看到的那样。

我们不能直接使用自回归模型，因为它需要时间序列是平稳的，这意味着均值和方差在时间上是恒定的——没有季节性，没有趋势。

我们可以使用 statsmodels 工具来查看时间序列的季节性和趋势性：

from statsmodels.tsa.seasonal import seasonal_decomposeresult = seasonal_decompose(df, model=\'additive\', period=52)result.plot()

我们将周期设置为 1，因为每个数据点（行）对应一年。

让我们看看各个组成部分的样子：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_06.png

图 5.8：时间序列的季节性分解

第一个子图是原始时间序列。这个数据集包含季节性和趋势性，您可以在子图中看到这些成分被分离出来。

如前所述，我们需要一个平稳的序列来进行建模。为了建立平稳性，我们需要去除季节性和趋势性成分。我们也可以去除我们之前估算的季节性或趋势性成分。或者，我们可以使用 statsmodels 中的封装功能，或者在 ARIMA 中设置 d 参数。

我们可以使用扩展的迪基-富勒（Augmented Dickey-Fuller）和 KPSS 检验来检查平稳性：

from arch.unitroot import KPSS, ADFADF(df1)

我们本可以使用statsmodels.tsa.stattools.adfuller或statsmodels.tsa.stattools.kpss，但我们更喜欢 ARCH 库版本的方便性。我们将留给用户检查 KPSS 检验的输出。我们得到了如下输出：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_07.png

图 5.9：KPSS 检验平稳性的输出

给定 p 值为 0.997，我们可以拒绝单位根的原假设，并得出结论，我们的过程是弱平稳的。

那么我们如何找到合适的差分值呢？我们可以使用 pmdarima 库，它提供了一个专门用于此目的的函数：

from pmdarima.arima.utils import ndiffs# ADF Test:ndiffs(df1, test=\'adf\')

我们得到了一个值 1。对于 KPSS 和 PP 检验，我们会得到相同的值。这意味着我们可以从第一次差分开始工作。

让我们从自回归模型开始。

提醒一下，ARIMA 模型由参数 p、d、q 来定义，其中：

  • p 表示自回归模型：AR§

  • d 表示积分

  • q 表示移动平均：MA(q)

因此，ARIMA(p, d, 0) 就是带有差分阶数 d 的 AR§模型。

知道 statsmodels 会检查并警告平稳性假设是否成立，让人放心。让我们尝试运行以拟合以下 AR 模型：

mod = sm.tsa.arima.ARIMA(endog=df, order=(1, 0, 0))res = mod.fit()print(res.summary()) 
UserWarning: Non-stationary starting autoregressive parameters found. Using zeros as starting parameters. warn(\'Non-stationary starting autoregressive parameters\'

由于我们已经知道需要进行一次差分处理，我们可以将 d 设置为 1。再试一次。这次，我们将使用 STLForecast 包装器，它能去除季节性并重新加回季节性。这是必要的，因为 ARIMA 无法直接处理季节性：

from statsmodels.tsa.forecasting.stl import STLForecastmod = STLForecast( df1, sm.tsa.arima.ARIMA, model_kwargs=dict(order=(1, 1, 0), trend=\"t\"))res = mod.fit().model_resultprint(res.summary())

我们得到以下总结：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_08.png

图 5.10：我们的 ARIMA 模型概述

这个结果总结提供了所有关键统计数据。我们看到模型是 ARIMA(1, 1, 0)。对数似然值为 -1965。我们还看到了 BIC 和 AIC 值，这些值可以用于模型选择。

请注意，我们在这里需要设置 trend=\"t\"，这样模型就会包含常数项。如果不设置，我们将得到一个虚假的回归结果。

我们如何使用这个模型？让我们进行一些预测！

STEPS = 20forecasts_df = res.get_forecast(steps=STEPS).summary_frame()

这给我们提供了一个未来 20 步的预测。

让我们来可视化一下！

ax = df1.plot(figsize=(12, 6))plt.ylabel(\'SPY\')forecasts_df[\'mean\'].plot(style=\'k--\')ax.fill_between( forecasts_df.index, forecasts_df[\'mean_ci_lower\'], forecasts_df[\'mean_ci_upper\'], color=\'k\', alpha=0.1)

这是我们得到的结果：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_09.png

图 5.11：SPY 股票代码的价格预测

实线表示我们已知的数据，虚线代表我们预测的未来 20 年的结果。我们预测值周围的灰色区域是 95% 的置信区间。

这个看起来还不错。作为读者的练习，可以尝试使用不同的参数。值得更改的有趋势参数和模型的阶数。

对于移动平均，我们创建不同的模型来观察它们的预测差异！

首先，我们将生成预测：

forecasts = []qs = []for q in range(0, 30, 10): mod = STLForecast( df1, sm.tsa.arima.ARIMA, model_kwargs=dict(order=(0, 1, q), trend=\"t\") ) res = mod.fit() print(f\"aic ({q}): {res.aic}\") forecasts.append( res.get_forecast(steps=STEPS).summary_frame()[\'mean\'] ) qs.append(q)forecasts_df = pd.concat(forecasts, axis=1)forecasts_df.columns = qs

在循环中，我们正在迭代不同的 q 参数，选择 0、10 和 20。我们使用这些 q 值估计移动平均模型，并预测未来 20 年的价格。同时，我们还打印出每个 q 对应的 AIC 值。这是我们得到的结果：

aic (0): 3989.0104184919096aic (10): 3934.375909262983aic (20): 3935.3355340835

现在，让我们像之前一样绘制三个预测图：

ax = df1.plot()plt.ylabel(\'SPY\')forecasts_df.plot(ax=ax)

这是新的图表：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_10.png

图 5.12：使用不同 q 参数的预测

那么，这些模型中哪一个在统计学上更好呢？

让我们回到 AIC。AIC 值越低，模型的效果越好，因为它考虑了对数似然和参数的数量。

在这种情况下，q=10 的阶数给出了最低的 AIC 值，根据这一标准，我们应该选择 q=10。当然，我们只尝试了三个不同的值。接下来我将留给读者作为练习，找出一个更合理的 q 参数值。

请注意，pmdarima 库具有查找最优参数值的功能，而 SkTime 库提供了 ARIMA 模型最优阶数的自动发现实现：AutoARIMA。

接下来，我们使用指数平滑模型进行预测。

在循环中，我们正在迭代不同的 q 参数，选择 0、10 和 20。我们使用这些 q 值估计移动平均模型，并预测未来 20 年的价格。同时，我们还打印出每个 q 对应的 AIC 值。这是我们得到的结果：

mod = sm.tsa.ExponentialSmoothing( endog=df1, trend=\'add\' )res = mod.fit()

这将模型拟合到我们的数据上。

让我们为接下来的 20 年获取预测：

forecasts = pd.Series(res.forecast(steps=STEPS))

现在，让我们绘制预测图：

ax = df.plot(figsize=(12, 6))plt.ylabel(\'SPY\')forecasts.plot(style=\'k--\')

这是图表：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_11.png

图 5.13：指数平滑预测

直到现在，我们只是看了 20 步预测的图表。我们还没有对模型的表现进行深入分析。让我们看看误差！

为此，我们首先需要将数据集划分为训练集和测试集。我们可以进行 n 步预测并检查误差。我们只需将截至某一时间点的时间序列作为训练数据，之后的时间点作为测试数据，然后将预测值与实际数据点进行比较：

from statsmodels.tsa.forecasting.theta import ThetaModeltrain_length = int(len(df1) * 0.8)tm = ThetaModel(df1[:train_length], method=\"auto\",deseasonalize=True)res = tm.fit()forecasts = res.forecast(steps=len(df1)-train_length)ax = df1.plot(figsize=(12, 6))plt.ylabel(\'SPY\')forecasts.plot(style=\'k--\')

这是图表：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_05_12.png

图 5.14：Theta 模型预测

虚线是预测值。它似乎与时间序列的实际行为不太吻合。让我们使用我们在上一章《时间序列的机器学习入门》中讨论过的误差度量来量化这个问题：

from sklearn import metricsmetrics.mean_squared_error(forecasts, df1[train_length:], squared=False)

我们得到一个值为37.06611385754943。这是均方根误差（我们将squared参数设置为False）。

在预测比赛中，比如 Kaggle 网站上的比赛，最低误差获胜。在现实中，简洁性（简单性）也很重要；然而，我们通常仍然追求尽可能低的误差（无论采用哪种度量标准）。

还有很多其他模型可以探索和尝试，但现在是时候总结本章内容了。

总结

在这一章中，我们讨论了基于移动平均和自回归的时间序列预测。这个主题包含了一大套模型，在不同领域（如计量经济学和统计学）中非常流行。这些模型是时间序列建模的支柱，并提供了最先进的预测方法。

我们已经讨论了自回归和移动平均模型，以及结合这两者的其他模型，包括 ARMA、ARIMA、VAR、GARCH 等。在实践环节中，我们已经将一些模型应用于股票价格数据集。

第六章：时间序列的无监督方法

我们在前一章已经讨论了预测方法，接下来我们将在下一章中讨论时间序列的预测。这些预测模型的性能很容易受到数据中重大变化的影响。识别这些变化是无监督学习的领域。

在本章中，我们将描述使用时间序列数据进行无监督学习的具体挑战。无监督学习的核心是从时间序列中提取结构，最重要的是识别子序列之间的相似性。这正是异常检测的本质（同样适用于：离群点检测），我们希望识别出与其他序列明显不同的序列。

时间序列数据通常是非平稳的、非线性的，并且是动态演化的。处理时间序列的一个重要挑战是识别潜在过程中的变化。这被称为变点检测（CPD）或漂移检测。数据随着时间的推移而变化，识别这些变化的程度至关重要。这是值得深入探讨的，因为变点和异常点的存在是现实应用中的常见问题。

在本章中，我们将集中讨论异常检测和变点检测（CPD），而在第八章中，时间序列的在线学习，我们将更详细地探讨漂移检测。我们将从概述和定义开始，然后研究大科技公司在行业中的实践。

我们将讨论以下主题：

  • 时间序列的无监督方法

  • 异常检测

  • 变更检测

  • 聚类

  • Python 实践

我们将从无监督学习与时间序列的总体介绍开始。

时间序列的无监督方法

时间序列与其他类型数据的主要区别在于它对时间轴的依赖；在某一点 t[1] 上的相关结构可能与在点 t[2] 上的相同结构包含完全不同的信息。时间序列通常包含大量噪声且维度较高。

为了减少噪声并降低维度，可以应用降维、波形分析或信号处理技术，例如傅里叶分解。这些通常是异常检测或变点检测（CPD）技术的基础，我们将在本章讨论这些技术。我们将在第八章中讨论漂移检测，在线时间序列方法。

我们将详细讨论异常值和变化点，查看它们的表现形式可能会有所帮助。在 Ilona Otto 等人的文章《2050 年通过社会动态稳定地球气候的社会临界动态》中，他们分析了基于社会动态的温室气体排放变化是否以及如何将各国转变为碳中和社会。他们根据不同的情景预测了全球变暖，以下图显示了 2010 年及 2020 年代初期的临界点（图表改编自他们的文章）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_01.png

图 6.1：基于温室气体排放的可能变化点

全球气温在冰川期和温暖期之间循环变化，每个周期大约持续数万年。在过去几千年里，气候逐渐变冷，直到 1970 年代，关于可能导致下一次冰河时期的降温趋势有了广泛的猜测。然而，数据表明，自工业化开始以来，主要由燃烧化石燃料推动，全球气温已上升了约 1°C。

因此，工业化时期的开始可以视为全球气温的一个变化点，如下图所示（来源：维基媒体共享资源）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_02.png

图 6.2：全球气温变化点：工业时代的开始

在上面的图表中，工业革命开始时的变化点出现在现代气温上升的异常值之前。

对于人类来说，指出变化点或异常值相对容易，特别是在事后，历史数据完全可得的情况下。对于自动检测，有很多不同的方法可以找到显著的点。在实际应用中，重要的是要仔细平衡检测率和假阳性。

异常检测

在异常检测中，我们希望识别出与其余序列明显不同的序列。异常值或离群点有时可能是测量误差或噪声的结果，但它们也可能表明被观察系统的行为发生了变化或出现了异常行为，这可能需要采取紧急措施。

异常检测的一个重要应用是对潜在复杂、高维数据集的自动实时监控。

是时候尝试给出一个定义了（参考 D.M. Hawkins, 1980，《异常值识别》）：

定义：异常值是指与其他观测值差异极大的数据点，以至于它可能是由不同的机制产生的。

让我们从一个图表开始，这样我们可以直观地看到异常值在图形上可能的表现。这也为我们的讨论提供了背景。

异常值检测方法可以分为单变量方法和多变量方法。参数化异常值检测方法通过选择其分布参数（例如，算术平均值），对底层分布做出假设——通常是高斯分布。这些方法标记出偏离模型假设的异常值。

在最简单的情况下，我们可以定义一个异常值，如下所示，观察值 x[i] 相对于分布参数的 z-得分：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_001.png

z-得分衡量每个点与移动平均或样本均值之间的距离，https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_002.png，以移动或样本标准差的单位！。对于高于均值的值，它为正；对于低于均值的值，它为负。

在此公式中，https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_004.png 和 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_005.png 是时间序列的估计均值和标准差，x 是我们想要测试的点。最后，https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_006.png 是一个依赖于我们感兴趣的置信区间的阈值——通常，选择 2 或 1.96，对应于 95% 的置信区间。通过这种方式，异常值是指发生概率小于或等于 5% 的点。

z-得分假设数据是正态分布的；然而，以上异常值公式中使用的均值和标准差可以被其他去除此假设的度量替代。像中位数或四分位距（在第二章中有讨论，使用 Python 进行时间序列分析）等度量对分布更加稳健。

Hampel 滤波器（也叫 Hampel 标识符）是这一特殊情况，其中使用了中位数和中位绝对偏差（MAD）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_007.png

在这个方程中，样本均值被（样本）中位数取代，标准差被 MAD（绝对中位差）取代，MAD 的定义如下：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_008.png

中位数是按顺序排列的数值列表中的中间数。

在 Hampel 滤波器中，每个观察值 x 将与中位数进行比较。在正态分布的情况下，Hampel 滤波器等同于 z-得分，epsilon 可以像 z-得分一样选择。

在多变量情况下，异常值函数可以表示为到模型分布中某一点的距离（或者相反：相似度），例如重心、均值。例如，我们可以计算新观察值与均值之间的协方差。

尽管这些前述方法仅限于低维或单变量时间序列，但基于距离的方法可以处理更大的空间。基于距离的异常值检测方法有效地将点聚类为不同的组，其中小组会被标记为异常值。在这些方法中，距离度量的选择至关重要。

检测时间序列中异常值的挑战之一是：

  • 缺乏异常值的定义

  • 输入数据中的噪声

  • 时间序列的复杂性

  • 高度不平衡

我们通常并不知道异常值长什么样。在实际应用中，我们往往没有异常值的标签——这使得基于真实案例的基准测试变得不可能。至于复杂性，时间序列会随着时间变化，它们通常是非平稳的，变量之间的依赖关系可能是非线性的。最后，我们通常有比异常值更多的正常观测值。

部署大规模异常检测模型作为服务的一个要求是，它们应该能够实时检测到异常。

异常检测的应用包括此图中的内容：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_03.png

图 6.3：异常检测的应用

一些例子可能包括支付中的欺诈检测、网络安全（网络入侵）、医学监控或传感器网络。在医学监控中，我们希望实时监测生理变量，包括心率、脑电图和心电图，以便在急性紧急情况下发出警报。传感器网络中的异常警报有助于防止工业损害的发生。

该图表说明了根据数据集的可用知识，异常检测方法的主要类型：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_04.png

图 6.4：根据可用知识的异常检测方法

最早的异常检测例子是基于规则的系统。当模式可以清晰定义时，这种方法有效。当我们有一个标注好的异常集时，我们可以应用监督或半监督方法，如分类器或回归模型。然而，最常见的使用场景是异常没有标注，我们需要无监督方法来基于密度或分布检测异常点或异常区间。

观察大科技公司（如 Alphabet（谷歌）、Amazon、Facebook、Apple 和 Microsoft（GAFAM））在异常检测方面的做法是很有启发性的。我们逐一看看它们是如何处理异常检测的。

微软

在论文《微软的时间序列异常检测服务》（Hansheng Ren 等，2019）中，介绍了一种为微软生产数据的异常检测而部署的时间序列服务。其核心是谱残差（SR）和卷积神经网络（CNN），应用于单变量时间序列的无监督在线异常检测。

它们借用了来自计算机视觉中显著性图（saliency map）概念的 SR 方法。显著性图突出显示图像中对人类观察者来说突出的点。该算法对数据执行傅里叶变换，然后应用变换信号的对数幅度的 SR，最后通过逆傅里叶变换将频谱数据投影回时域。

作为扩展，他们基于人工数据使用 SR 方法训练了一个卷积神经网络（CNN）。他们展示了在公开可用数据上的基准测试，支持他们的主张：他们的方法是异常检测领域的最新技术。

他们进一步声称，在微软生产数据上，他们的检测准确性（F1 得分）提高了超过 20%。你可以在 alibi-detect 库中找到基本实现（\"谱残差\"方法）。

Google

在 Google Analytics 的常见问题中（support.google.com/analytics/answer/7507748?hl=en），谷歌提到了一个贝叶斯状态空间时间序列模型（“用贝叶斯结构时间序列预测当前状态”，作者 Steven L. Scott 和 Hal Varian，2013），用于变化点和异常检测。

谷歌发布了一个具有更具体时间序列功能的 R 包——CausalImpact。描述该包背后研究的论文于 2015 年发布（“通过贝叶斯结构时间序列模型推断因果影响”，作者 Kay H. Brodersen, Fabian Gallusser, Jim Koehler, Nicolas Remy, Steven L. Scott）。CausalImpact 基于结构化贝叶斯时间序列模型估计干预的因果效应。该方法已被移植到 Python（pycausalimpact 库）。我们将在第九章，时间序列的概率模型中实验使用贝叶斯结构时间序列（BSTS）进行因果影响分析。

Amazon

亚马逊通过其亚马逊云服务（AWS）平台提供大规模的机器学习解决方案，其中包含异常检测功能，作为其资源和应用监控解决方案 CloudWatch 的一部分。尽管其解决方案的具体原理尚不清楚，但经济学家 Corey Quinn 在一条推文中推测，他们的解决方案可能是指数平滑。作为其中的一部分，他们很可能将季节性分解作为算法的第一步。

他们还有第二项异常检测服务：Amazon Lookout for Metrics。关于该服务的具体工作原理也不清楚。该服务旨在监控业务指标，并且——根据文档——在亚马逊内部用于大规模监控。在此服务中，用户可以从不同细分的数据源中选择字段，例如，通过选择数据库列page_views和device_type，用户可以分别查看每种设备类型下页面浏览量的异常变化。

至于亚马逊在异常检测领域的研究，他们在声学场景与事件检测与分类研讨会（DCASE 2020）中的 117 个参赛作品中荣获前三名。他们在这次与时间序列异常检测类似的挑战中获得了最佳论文奖，论文题为\"基于组掩蔽自编码器的音频异常检测密度估计器\"（作者 Ritwik Giri 等，2020）。

Facebook

Facebook 的核心数据科学团队在 GitHub 上开源了他们用于时间序列预测和异常检测的实现。这个库叫做 Prophet。在 2017 年宣布这个库的博客文章中，他们表示，Prophet 是 Facebook 能够大规模创建预测的关键工具，并且在决策过程中被认为是重要的信息来源。

Sean J Taylor 和 Benjamin Letham（2017）在论文《大规模预测》中描述了他们在 Facebook 的设置，包括一个分析师参与的环节，并且能够自动标记预测结果进行人工审查和调整。异常检测建立在来自广义加法模型（GAM）预测的不确定性基础上。

Prophet 已在基准测试中与其他概率模型和非概率模型进行了比较，且很少展现出突出的成功。microprediction.com 的 Elo 评分表明，Prophet 在单变量预测方面表现不如指数移动平均和许多其他标准方法。

Twitter

Twitter 也发布了一个 R 包，名为 AnomalyDetection。该方法基于广义极端学生化偏差（ESD）测试，用于检测单变量近似正态分布的时间序列中的异常。该方法发表于 2017 年（“通过统计学习在云端进行自动异常检测”，Jordan Hochenbaum, Owen Vallis, Arun Kejariwal）。

对于他们的 ESD 测试的适配，季节性混合 ESD 方法在应用阈值之前，加入了基于 LOESS 的季节性趋势分解（STL），并对 z-score 应用了阈值（如上所述），或者对于异常值较多的数据集，基于中位数和 MAD 进行阈值化。Twitter 模型已被移植到 Python 中（sesd 库）。

实现

我们将以 Python 中现有的异常检测实现概述来结束。市面上有很多实现方法。它们的使用场景非常相似，但实现方式和用户群体各不相同。

以下是按 GitHub 上的星标数量排序的列表（截至 2021 年 5 月）：

Library 实现 维护者 星标 Prophet 关于预测趋势成分的估计的不确定性区间 Facebook 核心研究 12.7k PyOD 适用于多变量时间序列的 30 种检测算法——从经典的 LOF（SIGMOD 2000）到 COPOD（ICDM 2020） Yue Zhao 及其他人 4.5k alibi-detect 多种异常检测算法——特定于时间序列的有似然比、Prophet、谱残差、Seq2Seq、模型蒸馏 Seldon Technologies Ltd 683 Scikit-Lego 通过 PCA/UMAP 重构 Vincent D. Warmerdam 及其他人 499 Luminaire Luminaire 窗口密度模型 Zillow 371 Donut 用于季节性关键绩效指标（KPI）的变分自编码器 清华大学网络管理实验室 327 rrcf 用于流数据异常检测的鲁棒随机切割森林算法 密歇根大学实时水系统实验室 302 banpei 霍特林理论 Hirofumi Tsuruta 245 STUMPY 针对单变量和多变量时间序列的矩阵概况算法，如 STUMP、FLUSS 和 FLOSS（也可参考 matrixprofile-ts） TD Ameritrade 169 PySAD 超过十种流式异常检测算法 Selim Yilmaz, Selim 和 Suleyman Kozat 98

图 6.5：Python 中的异常检测方法

这些方法各自有其背景和形式化基础；然而，本章的范围并不涵盖对它们的详细描述。

这张图展示了三大仓库的星标历史（来自 star-history.t9t.io）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_05.png

图 6.6：Prophet、PyOD 和 alibi-detect 的星标历史

Prophet 和 PyOD 的受欢迎程度（GitHub 星标数）一直在持续增长。

最近，许多深度学习算法已经被应用于异常检测，既包括单变量时间序列，也包括多变量时间序列。

深度学习模型的特别之处在于，应用范围可以更加广泛：例如闭路电视中的视频监控异常检测。我们将在第十章，时间序列的深度学习中更详细地探讨深度学习架构。

变化点检测

时间序列的一个常见问题是观察系统行为的变化。一般来说，变化点表示在生成该序列的过程中，系统状态之间发生了突变和重大转变。例如，趋势可能会突然发生变化，而变化点可以指示趋势变化的位置。这在交易中的技术图表模式分析中非常常见。

这个列表展示了一些 变化点检测（CPD）的应用：

  • 语音识别：检测单词和句子的边界

  • 图像分析：对闭路电视视频监控进行监视

  • 健身：根据智能设备（如手表或手机）上的运动传感器数据，分割人的活动时间段

  • 金融：识别趋势模式的变化，可能表明从熊市到牛市，或反之的转变。

以股票市场为例，说明 CPD 的重要性。描述市场演变的时间序列数据，如股票价格，遵循趋势——它要么上涨，要么下跌，或者没有显著变化（停滞）。

当股票上涨时，投资者想要买入该股票。否则，当股票下跌时，投资者不希望持有该股票，而是希望将其卖出。不改变仓位会导致账面价值的损失——在最好的情况下，这会导致流动性问题。

对于投资者而言，了解市场从上涨到下跌，或从下跌到上涨的变化时机至关重要。识别这些变化可能决定是否盈利。

在预测中，诸如黑色星期五、圣诞节、选举、新闻发布或法规变动等特殊事件可能会对趋势或序列的水平造成短期（可能当时被视为异常）或长期的变化。这将不可避免地导致传统模型产生异常的预测。

CPD 算法面临的一个特别有趣的挑战是实时检测这些拐点。这意味着在拐点到来时立即检测到变化点（或者至少在下一个变化点发生之前）。

我们可以区分 CPD 的在线和离线方法，其中在线指的是实时处理，处理每一个新到的数据点。而离线算法则可以一次性处理整个时间序列。我们将在第八章，时间序列的在线学习中更多地讨论在线处理。

CPD 与分段、边缘检测、事件检测和异常检测相关，类似的技术可以应用于所有这些应用。CPD 可以看作与异常检测非常相似，因为识别变化点的一种方法是通过异常检测算法的异常分数。

从这个角度来看，变化点与高度异常点是相同的，任何超过某个阈值的点都对应一个变化。与异常检测相似，CPD 可以定义为在两个备择假设之间进行假设检验的问题，零假设为\"没有变化发生\"，而备择假设为\"发生了变化。\"

CPD 算法由三个组成部分构成：代价函数、搜索方法和约束条件。我们将逐一讲解这些内容。代价函数是可以应用于时间序列的一个子段（多变量或单变量）的距离函数。

代价函数的一个例子是 最小绝对偏差（LAD），它是一个估算分布中心点（均值、中位数和众数）变化的估计量，定义如下：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_009.png

在这个定义中，l 是时间序列 x 中一个子段的索引，而 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_010.png 是 x 的中心点。

搜索函数然后会遍历时间序列来检测变化点。这可以是近似的，例如基于窗口的检测、从下到上的方法或二分分割，或者是穷举的，如动态规划或 修剪精确线性时间（Pelt）方法。

Pelt（Gachomo Dorcas Wambui 等，2015）依赖于修剪启发式方法，计算成本与时间序列的点数成线性关系，https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_011.png。动态规划方法的计算成本要高得多，https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_012.png，其中 n 是预期变化点的最大数量。

最后，约束条件可以作为搜索算法中的惩罚项参与其中。这个惩罚项可以编码为一个成本预算，或者是我们预期找到的变化点数量的知识。

评估 CPD 算法的性能一直是个难题，因为缺乏基准数据集。直到最近（2020 年），来自阿兰·图灵研究所和爱丁堡大学的 Gerrit van den Burg 和 Christopher Williams 发布了一个基准数据集，包含来自世界银行、欧盟统计局、美国人口普查局、GapMinder 和维基百科等来源的 37 个时间序列。他们的基准数据集已发布在 GitHub 上，并且提到该数据集的变点注释集中在 2007-08 年的金融危机、英国的安全带法规、蒙特利尔议定书（用于调控氯氟烃排放）以及美国的自动电话呼叫监管等地方。

在同一篇论文（“An Evaluation of Change Point Detection Algorithm”）中，作者评估了各种 CPD 方法。他们指出，假设没有任何变点的“零”基准方法在 F1 度量和基于 Jaccard 指数的聚类重叠度量上优于许多其他方法。这是因为数据集中变点的比例很小，且这些方法返回了大量的假阳性结果。他们得出结论，二分法和贝叶斯在线 CPD 是时间序列中最有效的几种方法。

二分法（“On Tests for Detecting Change in Mean”，作者：Ashish Sen 和 Muni S. Srivastava，1975）属于基于窗口的 CPD 方法。二分法是一种贪心算法，它通过如下定义最小化代价之和：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_013.png

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_014.png是发现的变点，而*c()*是类似于 LAD 的代价函数，我们之前在本节中看到了。其基本思想是，当两个子序列高度不相似时，表示存在变点。

二分法是顺序执行的，这意味着首先在整个时间序列中检测变点，然后在变点前后两个子序列中再次检测。这也解释了其低复杂度！，其中 T 为时间序列的长度。这个计算成本使得它可以扩展到更大的数据集。

该表格概述了 CPD 方法的不同种类：

库 实现方式 维护者 Greykite 通过自适应 lasso 进行 CPD LinkedIn ruptures 离线 CPD：二分法、动态规划、Pelt、基于窗口的 Charles Truong 贝叶斯变点检测 贝叶斯 CPD Johannes Kulick banpei 奇异谱变换 Hirofumi Tsuruta changepy Pelt 算法 Rui Gil onlineRPCA 在线移动窗口稳健主成分分析 Wei Xiao

图 6.7：Python 中的 CPD 方法

我们省略了 Facebook 的 Prophet 库，因为它不是专门的 CPD 包。

下图展示了 CPD 方法随时间变化的受欢迎程度。

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_06.png

图 6.8：CPD 方法的历史演变

LinkedIn 的 Greykite 自发布以来在 GitHub 星标上迅速增长。同时，ruptures 也在流行度上大幅上升。

聚类

聚类分析或聚类是根据数据点或对象的相似性，在数据集中寻找有意义的组（簇）的过程。作为这种无监督数据挖掘技术的结果，我们希望每个簇中的点彼此相似，同时与其他簇中的点有所不同。

时间序列的聚类具有挑战性，因为每个数据点都是一个时间段（有序序列）。它已在多个领域得到应用，帮助发现模式，推动时间序列分析，从复杂数据集中提取洞察。

我们不打算深入讨论时间序列聚类，但下表提供了 Python 库在时间序列聚类中的概述：

库 实现 维护者 星标 tslearn 时间序列 K 均值、K-Shape 聚类、KernelKMeans Romain Tavenard 1.7k river DBStream、时间序列 K 均值、CluStream、DenStream、STREAMKMeans Albert Bifet 及其他人 1.7k

图 6.9：Python 中时间序列的聚类方法

你可以查看历史上顶级实现的 GitHub 星标数据：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_07.png

图 6.10：tslearn 和 river 的星标历史

两个库都在蓬勃发展。我们将在第八章，时间序列的在线学习中重新讨论 river。

Python 实践

首先让我们做一个异常检测的例子，然后再做一个 CPD 的例子。我们先看看下一节所需的库。

要求

在本章中，我们将使用多个库，这些库可以通过终端快速安装（或通过 anaconda navigator 进行安装）：

pip install ruptures alibi_detect

我们将从 Python（或 IPython）终端执行命令，但同样可以从 Jupyter 笔记本（或其他环境）中执行它们。

我们现在应该准备好深入实施 Python 中的无监督时间序列算法了。

异常检测

alibi-detect 提供了多个用于时间序列异常检测的基准数据集：

  • fetch_ecg—来自 BIDMC 充血性心力衰竭数据库的 ECG 数据集

  • fetch_nab—Numenta 异常基准

  • fetch_kdd—KDD Cup \'99 计算机网络入侵数据集

这些中的最后一个是通过 scikit-learn 加载的。

让我们加载计算机网络入侵的时间序列（KDD99）：

from alibi_detect.datasets import fetch_kddintrusions = fetch_kdd()

intrusions是一个字典，其中data键返回一个 494021x18 的矩阵。时间序列的 18 个维度是数据集的连续特征，主要是误差率和计数：

intrusions[\'feature_names\'][\'srv_count\', \'serror_rate\', \'srv_serror_rate\', \'rerror_rate\', \'srv_rerror_rate\', \'same_srv_rate\', \'diff_srv_rate\', \'srv_diff_host_rate\', \'dst_host_count\', \'dst_host_srv_count\', \'dst_host_same_srv_rate\', \'dst_host_diff_srv_rate\', \'dst_host_same_src_port_rate\', \'dst_host_srv_diff_host_rate\', \'dst_host_serror_rate\', \'dst_host_srv_serror_rate\', \'dst_host_rerror_rate\', \'dst_host_srv_rerror_rate\']

另一个键target包含异常的注释。

既然我们已经准备好了注释，我们本可以训练一个分类器，然而，我们将坚持使用无监督方法。此外，由于我们将使用的谱方法适用于单变量数据，我们只会从多变量数据集中提取一个维度，因此我们将完全忽略注释。

这是我们时间序列的快速图（我们将随意选择数据集的第一个维度）：

import pandas as pdpd.Series(intrusions[\'data\'][:, 0]).plot()

这是图表：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_08.png

图 6.11：时间序列图

我们将加载并运行实现微软提出方法的 SpectralResidual 模型：

from alibi_detect.od import SpectralResidualod = SpectralResidual( threshold=1., window_amp=20, window_local=20, n_est_points=10, n_grad_points=5)

然后，我们可以获取时间序列中每个点的异常得分：

scores = od.score(intrusions[\'data\'][:, 0])

让我们将得分绘制在时间序列之上！

import matplotlibax = pd.Series(intrusions[\'data\'][:, 0], name=\'data\').plot(legend=False, figsize=(12, 6))ax2 = ax.twinx()ax = pd.Series(scores, name=\'scores\').plot(ax=ax2, legend=False, color=\"r\", marker=matplotlib.markers.CARETDOWNBASE)ax.figure.legend(bbox_to_anchor=(1, 1), loc=\'upper left\');

我们在同一个图表中使用双 y 轴来绘制得分和数据。如下所示：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_09.png

图 6.12：带有异常值的时间序列

由于傅里叶滤波器移除了信号的周期性，一些点未被识别为离群点。

变化点检测

我们将首先使用 ruptures 库创建一个合成的多元时间序列。我们将维度数设置为 3，时间序列的长度设置为 500，并且我们的时间序列将有 3 个变化点，标准差为 5.0 的高斯噪声将被叠加：

import numpy as npimport matplotlib.pylab as pltimport ruptures as rptsignal, bkps = rpt.pw_constant( n_samples=500, n_features=3, n_bkps=3, noise_std=5.0, delta=(1, 20))

信号是一个 500x3 的 NumPy 数组。bkps是变化点的数组（123、251 和 378）。

我们可以使用一个实用函数绘制该时间序列，函数会突出显示由变化点分隔的子部分：

 rpt.display(signal, bkps)

这是我们带有三个变化点的时间序列图：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_10.png

图 6.13：带有变化点的时间序列

我们可以将二分段方法应用于此时间序列。ruptures 遵循 scikit-learn 的约定，因此，如果你之前使用过 scikit-learn，使用起来应该非常直观：

 algo = rpt.Binseg(model=\"l1\").fit(signal)my_bkps = algo.predict(n_bkps=3)

我们有几个二分段约束选项——可以选择l1、l2、rbf、linear、normal和ar。

我们可以用另一个实用函数绘制二分段的预测结果：

 rpt.show.display(signal, bkps, my_bkps, figsize=(10, 6))

这是我们从二分段模型中获得的变化点预测图：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_06_11.png

图 6.14：带有检测到的变化点的时间序列（Binary Segmentation）

让我们总结一下本章中的一些信息！

总结

在本章中，我们集中讨论了时间序列无监督方法的两个方面：

  • 异常检测

  • 变化点检测

异常检测（也叫：离群点检测）的本质是识别出明显与其他序列不同的部分。我们已经研究了不同的异常检测方法，以及一些大公司如何在大规模中处理它。

在处理时间序列时，重要的是要注意数据随时间变化，这可能导致模型失效（模型陈旧）。这就是所谓的变化点检测和漂移检测。

本章我们探讨了变化点检测。在第八章，时间序列的在线学习中，我们将更详细地探讨漂移检测。

第七章：时间序列的机器学习模型

近年来，机器学习取得了长足的进展，这一点在时间序列预测方法中得到了体现。我们在第四章，时间序列机器学习导论中介绍了一些最先进的时间序列机器学习方法。在本章中，我们将介绍更多的机器学习方法。

我们将介绍一些常用的基准方法，或者在性能、易用性或适用性方面表现突出的算法。我将介绍基于动态时间规整和梯度提升的 k 近邻算法作为基准，我们还将讨论其他方法，如 Silverkite 和梯度提升。最后，我们将进行一个应用练习，使用这些方法中的一些。

我们将涵盖以下主题：

  • 更多的时间序列机器学习方法

  • 带有动态时间规整的 k 近邻算法

  • Silverkite

  • 梯度提升

  • Python 练习

如果你想了解最先进的机器学习算法，请参阅第四章，时间序列机器学习导论。本章的算法讨论将假设你已经掌握了该章节的一些内容。我们将在接下来的部分中介绍的算法在预测和预报任务中都非常具有竞争力。

我们将在这里更详细地讨论这些算法。

更多的时间序列机器学习方法

我们将在本节中介绍的算法在预测和预报任务中都具有高度竞争性。如果你想了解最先进的机器学习算法，请参阅第四章，时间序列机器学习导论。

在上述章节中，我们简要讨论了其中的一些算法，但我们将在这里更详细地讲解它们，还会介绍一些我们之前未曾讨论过的算法，如 Silverkite、梯度提升和 k 近邻。

我们将专门设置一个实践部分，介绍一个 2021 年发布的库——Facebook 的 Kats。Kats 提供了许多高级功能，包括超参数调优和集成学习。在这些功能的基础上，它们实现了基于 TSFresh 库的特征提取，并包含多个模型，包括 Prophet、SARIMA 等。它们声称，与其他超参数调优算法相比，Kats 在时间序列上的超参数调优在基准测试中快了 6 到 20 倍。

该图展示了所选时间序列机器学习库的流行度概览：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_01.png

图 7.1：时间序列机器学习库的流行度

截至 2021 年中期，Kats 和 GreyKite 库最近才发布，尽管它们在 GitHub 上获得了越来越多的星标，但它们的受欢迎程度还未能与 TSFresh 相抗衡。我尽管 TSFresh 是一个特征生成库，而非预测库，但还是将其包含在内，因为我觉得它在本章所使用的其他库中非常重要。TSFresh 之后，SKTime 排名第二，并且在相对较短的时间内吸引了大量星标。

在本章的实际示例中，我们将使用一些这些库。

另一个重要的问题是验证，值得单独讨论这个问题。

验证

我们在第四章，时间序列的机器学习入门中已经讨论过验证。通常在机器学习任务中，我们使用 k 折交叉验证，其中数据的拆分是伪随机进行的，因此训练集和测试/验证集可以来自数据的任何部分，只要这些部分没有被用于训练（样本外数据）。

对于时间序列数据，这种验证方法可能会导致对模型性能的过度自信，因为现实中，时间序列往往随着趋势、季节性和时间序列特征的变化而变化。

因此，在时间序列中，验证通常采用所谓的前向验证（walk-forward validation）。这意味着我们在过去的数据上训练模型，然后在最新的数据片段上进行测试。这将消除过于乐观的偏差，并在模型部署后为我们提供更为现实的性能评估。

在训练、验证和测试数据集方面，这意味着我们将完全依赖训练和验证数据集来调整模型参数，并且我们将基于一个时间上更先进的数据集来评估测试，具体如下面的图示所示（来源：Greykite 库的 GitHub 仓库）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_02.png

图 7.2：前向验证

在前向验证中，我们先在数据的初始片段上训练，然后在训练集之后的某个时间段进行测试。接着，我们向前推进并重复这一过程。这样，我们有多个样本外的时间段，可以将这些时间段的结果进行整合。通过前向验证，我们不太可能遭遇过拟合的问题。

动态时间规整的 K 近邻算法

K 近邻是一个著名的机器学习方法（有时也被称为基于案例的推理）。在 kNN 中，我们可以使用距离度量来找到相似的数据点。然后，我们可以将这些最近邻的已知标签作为输出，并通过某种函数将它们整合在一起。

图 7.3 展示了 kNN 分类的基本思路（来源 – WikiMedia Commons: commons.wikimedia.org/wiki/File:KnnClassification.svg）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_03.png

图 7.3：用于分类的 K 最近邻

我们已经知道一些数据点。在前面的示意图中，这些数据点分别用方形和三角形表示，代表了两个不同类别的数据点。给定一个新的数据点，用圆圈表示，我们找到与之最接近的已知数据点。在这个例子中，我们发现新点与三角形相似，因此我们可以假设新点也属于三角形类。

尽管这种方法在概念上非常简单，但它通常作为一种强基线方法，或者有时甚至能与更复杂的机器学习算法竞争，即使我们只比较最接近的邻居（（𝑘=1））。

该算法中的重要超参数有：

  • 你希望根据其来生成输出的邻居数（k）

  • 集成函数（例如，平均值或最常见值）

  • 用于找到最近数据点的距离函数

我们在第四章《时间序列机器学习导论》中讨论过动态时间扭曲（Dynamic Time Warping），它是一种用于比较两条时间序列相似性（或等价地，距离）的方法。这些序列甚至可以有不同的长度。动态时间扭曲已经证明自己是一个异常强大的时间序列距离度量。

我们可以将 kNN 和动态时间扭曲结合起来，作为一种距离度量来寻找相似的时间序列，这种方法已经证明其在某些情况下很难被超越，尽管当前的技术已经超过了它。

Silverkite

Silverkite 算法与 LinkedIn 发布的 Greykite 库一起发布。它的设计目标是快速、准确且直观。该算法在 2021 年的文章中有描述（《生产系统的灵活预测模型》，Reza Hosseini 等人）。

根据 LinkedIn，它能够处理各种趋势和季节性因素，例如每小时、每日、每周、重复事件、假期以及短期效应。在 LinkedIn 内部，它既用于短期预测，例如 1 天的预测，也用于长期预测，例如 1 年后的预测。

在 LinkedIn 中的应用场景包括优化预算决策、设定业务指标目标，以及提供足够的基础设施来应对峰值流量。此外，一个应用场景是建模 COVID-19 大流行后的恢复情况。

时间序列被建模为趋势、变化点和季节性的加性组合，其中季节性包括假期/事件效应。趋势随后被建模如下：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_001.png

其中 K 是变化点的数量，t[i] 是第 i 个变化点的时间索引。因此，https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_002.png 是第 i 个变化点的指示函数。函数 f(t) 可以是线性、平方根、二次、任意组合或完全自定义。

Silverkite 还构建了假期的指示变量。假期可以通过名称或国家指定，甚至可以完全自定义。

变化点可以手动指定，或者可以通过回归模型自动检测候选变化点，随后使用自适应 Lasso 算法（Hui Zhou，2006）进行选择。

除了趋势、季节性和假期，Silverkite 还包括一个自回归项，该项是基于窗口平均值计算的，而不是独立地取滞后项（《为降水过程选择二元马尔科夫模型》，Reza Hosseini 等，2011）。

该自回归项使用 Pasty 库指定，使用类似以下字符串的公式迷你语言形式：

y ~ a + a:b + np.log(x)

在这个公式中，左侧的 y 定义为三个项的总和，a、a:b 和 np.log(x)。项 a:b 是两个因子 a 和 b 之间的交互作用。Pasty 中的模型模板本身具有高度的可定制性，因此该接口提供了很高的灵活性。

最后，Silverkite 提供了几种模型类型，如岭回归、弹性网和提升树，并支持损失函数，如均方误差（MSE）和分位数损失，用于稳健回归。

根据 LinkedIn 对多个数据集的基准测试，Silverkite 在预测误差方面优于 auto-Arima（pmdarima 库）和 Prophet。然而，Silverkite 的速度约为 Prophet 的四倍，我们将在 第九章，概率模型 中介绍 Prophet。

梯度提升

XGBoost（即极限梯度提升的缩写）是梯度提升的高效实现（Jerome Friedman，《贪婪函数逼近：一种梯度提升机器》，2001），用于分类和回归问题。梯度提升也被称为梯度提升机（GBM）或梯度提升回归树（GBRT）。一个特例是用于排序应用的 LambdaMART。除了 XGBoost，其他实现包括微软的轻量级梯度提升机（LightGBM）和 Yandex 的 Catboost。

梯度提升树是树的集成。这与类似于随机森林的袋装算法类似；然而，由于这是一个提升算法，每棵树都会被计算出来，逐步减少误差。每次新的迭代都会贪心地选择一棵树，并将其预测值基于权重项加到先前的预测值中。还有一个正则化项，用于惩罚复杂性并减少过拟合，类似于正则化贪婪森林（RGF）。

XGBoost算法由陈天奇和卡洛斯·格斯特林于 2016 年发布（“XGBoost：一种可扩展的树提升系统”），并推动了许多分类和回归基准的突破。它被用于许多 Kaggle 问题的获胜解决方案。事实上，在 2015 年，29 个挑战获胜解决方案中，有 17 个解决方案使用了 XGBoost。

它的设计非常具有可扩展性，并且扩展了梯度提升算法，用于加权分位数，并通过更智能的缓存模式、分片以及对稀疏性处理的改进，提升了可扩展性和并行化性能。

作为回归的一个特例，XGBoost 可以用于预测。在这种情况下，模型基于过去的值进行训练，以预测未来的值，这可以应用于单变量和多变量时间序列。

Python 练习

让我们将本章迄今为止学到的内容付诸实践。

至于依赖项，在本章中，我们将分别为每个部分安装依赖。可以通过终端、笔记本或 Anaconda Navigator 进行安装。

在接下来的几个部分中，我们将演示如何在预测中进行分类，因此这些方法中的一些可能无法进行比较。我们邀请读者使用每种方法进行预测和分类，并进行结果比较。

需要注意的是，Kats 和 Greykite（在写作时）都是非常新的库，因此它们的依赖项可能还会频繁变化。它们可能会锁定您的 NumPy 版本或其他常用库。因此，我建议您为每个部分分别在虚拟环境中安装这些库。

我们将在下一节中详细介绍这个设置过程。

虚拟环境

在 Python 虚拟环境中，安装到其中的所有库、二进制文件和脚本都是与其他虚拟环境中安装的内容以及系统中安装的内容隔离的。这意味着我们可以安装不同的库，如 Kats 和 Greykite，而无需担心它们之间的兼容性问题，或与我们计算机上安装的其他库之间的兼容性问题。

让我们通过一个简短的教程，介绍如何在使用 Anaconda 的 Jupyter 笔记本中使用虚拟环境（类似地，您也可以使用如 virtualenv 或 pipenv 等工具）。

在第一章，使用 Python 进行时间序列分析简介中，我们介绍了 Anaconda 的安装，因此我们会跳过这部分安装。请参考那一章，或者访问 conda.io 获取安装说明。

要创建一个虚拟环境，必须指定一个名称：

conda create --name myenv

这将创建一个同名的目录（myenv），其中所有库和脚本将被安装。

如果我们想使用这个环境，我们必须首先激活它，这意味着我们需要将PATH变量设置为包括我们新创建的目录：

conda activate myenv

现在我们可以使用像 pip 这样的工具，默认情况下它会使用与 conda 捆绑在一起的版本，或者直接使用 conda 命令来安装库。

我们可以在环境中安装 Jupyter 或 Jupyter Lab 然后启动它。这意味着我们的 Jupyter 环境将包含所有依赖项，因为我们已将它们单独安装。

让我们从一个带有动态时间规整（DTW）的 kNN 算法开始。正如我提到的，这个算法通常作为一个不错的基准进行比较。

使用动态时间规整（DTW）的 k-近邻算法（kNN）在 Python 中实现

在这一节中，我们将基于机器人在一段时间内的力和扭矩测量来分类故障。

我们将使用一个非常简单的分类器，kNN，并且或许我们应该提醒一下，这种方法涉及到逐点计算距离，这通常会成为计算瓶颈。

在这一节中，我们将把 TSFresh 的特征提取与 kNN 算法结合在一个管道中。时间序列管道确实能帮助我们简化过程，正如你在阅读代码片段时会发现的那样。

让我们安装 tsfresh 和 tslearn：

pip install tsfresh tslearn

我们将在 tslearn 中使用 kNN 分类器。我们甚至可以使用 scikit-learn 中的 kNN 分类器，它允许指定自定义的度量标准。

在这个示例中，我们将从 UCI 机器学习库下载一个机器人执行故障的数据集并将其存储到本地。该数据集包含故障检测后的机器人力和扭矩测量数据。对于每个样本，任务是分类机器人是否会报告故障：

from tsfresh.examples import load_robot_execution_failuresfrom tsfresh.examples.robot_execution_failures import download_robot_execution_failuresdownload_robot_execution_failures()df_ts, y = load_robot_execution_failures()

列包括时间和六个来自传感器的时间序列信号，分别是 F_x、F_y、F_z、T_x、T_y 和 T_z。目标变量 y，其值可以是 True 或 False，表示是否发生了故障。

始终检查两个类别的频率非常重要：

print(f\"{y.mean():.2f}\")

y 的均值是 0.24。

然后，我们可以使用 TSFresh 提取时间序列特征，正如在第三章，时间序列预处理中讨论的那样。我们可以填补缺失值，并根据与目标变量的相关性选择特征。在 TSFresh 中，统计测试的 p 值被用来计算特征的重要性：

from tsfresh import extract_featuresfrom tsfresh import select_featuresfrom tsfresh.utilities.dataframe_functions import imputeextracted_features = impute(extract_features(df_ts, column_id=\"id\", column_sort=\"time\"))features_filtered = select_features(extracted_features, y)

我们可以继续使用 features_filtered DataFrame，它包含我们的特征——来自传感器的信号和 TSFresh 特征。

让我们通过进行网格搜索来找到一个合适的邻居数值：

from sklearn.model_selection import TimeSeriesSplit, GridSearchCVfrom tslearn.neighbors import KNeighborsTimeSeriesClassifierknn = KNeighborsTimeSeriesClassifier()param_search = { \'metric\' : [\'dtw\'], \'n_neighbors\': [1, 2, 3]}tscv = TimeSeriesSplit(n_splits=2)gsearch = GridSearchCV( estimator=knn, cv=tscv, param_grid=param_search)gsearch.fit( features_filtered, y)

我们正在使用 scikit-learn 的 TimeSeriesSplit 来划分时间序列。这是为了网格搜索（GridSearch）。

或者，我们也可以仅仅基于索引进行拆分。

我们可以尝试许多参数，特别是在 kNN 分类器中的距离度量。如果你想尝试一下，请参阅 TSLEARN_VALID_METRICS 获取 tslearn 支持的度量标准的完整列表。

让我们预测一些 COVID 案例。在下一节中，我们将从 Silverkite 算法开始。Silverkite 是 LinkedIn 在 2021 年发布的 Greykite 库的一部分。

Silverkite

在写作时，Greykite 的版本为 0.1.1——它尚未完全稳定。它的依赖项可能与一些常用库的较新版本发生冲突，包括 Jupyter Notebooks。不过，如果你在虚拟环境或 Google Colab 中安装该库，不必担心。

只需安装该库及其所有依赖项：

pip install greykite

现在 Greykite 已经安装好了，我们可以使用它了。

我们将加载来自Our World in Data数据集的 COVID 病例数据，它可能是可用 COVID 数据中最好的来源之一：

import pandas as pdowid_covid = pd.read_csv(\"**https://covid.ourworldindata.org/data/owid-covid-data.csv**\")owid_covid[\"**date**\"] = pd.to_datetime(owid_covid[\"**date**\"])df = owid_covid[owid_covid.location == \"**France**\"].set_index(\"**date**\", drop=True).resample(\'**D**\').interpolate(method=\'**linear**\')

我们专注于法国的病例。

我们首先设置 Greykite 的元数据参数。然后，我们将此对象传递给预测器配置：

from greykite.framework.templates.autogen.forecast_config import ( ForecastConfig, MetadataParam)metadata = MetadataParam( time_col=\"date\", value_col=\"new_cases\", freq=\"D\")

我们的时间列是date，值列是new_cases。

现在我们将创建forecaster对象，它用于生成预测并存储结果：

import warningsfrom greykite.framework.templates.forecaster import Forecasterfrom greykite.framework.templates.model_templates import ModelTemplateEnumforecaster = Forecaster() warnings.filterwarnings(\"ignore\", category=UserWarning) result = forecaster.run_forecast_config( df=yahoo_df, config=ForecastConfig( model_template=ModelTemplateEnum.SILVERKITE_DAILY_90.name, forecast_horizon=90, coverage=0.95, metadata_param=metadata ) )

预测的时间跨度为 90 天；我们将预测未来 90 天。我们的预测区间为 95%。Silverkite 和 Prophet 都支持通过预测区间来量化不确定性。95%的覆盖率意味着 95%的实际值应该落在预测区间内。在 Greykite 中，_components.uncertainty模型提供了关于不确定性的额外配置选项。

我已经添加了一行代码，在训练过程中忽略UserWarning类型的警告，否则会有大约 500 行关于目标列中 0 值的警告。

让我们从结果对象中绘制原始时间序列图，我们可以将预测结果叠加在其上：

forecast = result.forecastforecast.plot().show(renderer=\"**colab**\")

如果你不是在 Google Colab 中，请不要使用renderer参数！

我们得到了以下图表：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_04.png

图 7.4：预测与实际时间序列（Silverkite）

预测结果存储在forecast对象的df属性中：

forecast.df.head().round(2)

这些是预测结果的上限和下限置信区间：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_05.png

图 7.5：预测与实际时间序列的表格（Silverkite）

我们可能需要获取一些关于模型的性能指标。我们可以像这样获取历史预测在保留测试集上的性能：

from collections import defaultdictbacktest = result.backtestbacktest_eval = defaultdict(list)for metric, value in backtest.train_evaluation.items(): backtest_eval[metric].append(value) backtest_eval[metric].append(backtest.test_evaluation[metric])metrics = pd.DataFrame(backtest_eval, index=[\"train\", \"test\"]).Tmetrics.head()

我们的性能指标如下所示：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_06.png

图 7.6：在保留数据上的性能指标（Silverkite）

我已将指标缩减为前五个。

我们可以方便地将模型应用到新的数据上，方法如下：

model = result.modelfuture_df = result.timeseries.make_future_dataframe( periods=4, include_history=False)model.predict(future_df)

预测结果如下所示：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_07.png

图 7.7：预测数据框（Silverkite）

请注意，您的结果可能会有所不同。

通过更改预测器运行配置中的model_template参数，我们可以使用其他预测模型。例如，我们可以将其设置为ModelTemplateEnum.PROPHET.name，以使用 Facebook 的 Prophet 模型。

这就结束了我们对 Silverkite 的介绍。接下来，我们将通过应用 XGBoost 的监督回归方法进行预测。让我们进行一些梯度提升吧！

梯度提升

我们也可以使用监督式机器学习进行时间序列预测。为此，我们可以使用日期和前期值来预测未来。

首先，我们需要安装 XGBoost：

pip install xgboost

在这个示例中，我们将使用 Yahoo 的每日收盘数据，和本章其他实践部分一样。

让我们一步一步地通过准备和建模过程。

我们首先需要对数据进行特征化。这里，我们通过提取日期特征来做到这一点，但请参见 kNN 部分，在那里我们使用了 TSFresh 的特征提取。你也许想通过结合这两种特征提取策略，或完全依赖 TSFresh，来改变这个示例。

我们将像之前一样重新加载来自Our World in Data数据集的 COVID 新病例数据：

import pandas as pdowid_covid = pd.read_csv(\"**https://covid.ourworldindata.org/data/owid-covid-data.csv**\")owid_covid[\"**date**\"] = pd.to_datetime(owid_covid[\"**date**\"])df = owid_covid[owid_covid.location == \"**France**\"].set_index(\"**date**\", drop=True).resample(\'**D**\').interpolate(method=\'**linear**\').reset_index()

对于特征提取，转换器非常方便。转换器基本上是一个类，包含fit()和transform()方法，可以使转换器适应数据集并相应地转换数据。以下是用于根据日期注释数据集的DateFeatures转换器的代码：

from sklearn.base import TransformerMixin, BaseEstimatorclass DateFeatures(TransformerMixin, BaseEstimator): features = [ \"hour\", \"year\", \"day\", \"weekday\", \"month\", \"quarter\", ] def __init__(self): super().__init__() def transform(self, df: pd.DataFrame): Xt = [] for col in df.columns: for feature in self.features: date_feature = getattr(  getattr( df[col], \"dt\"  ), feature ) date_feature.name = f\"{col}_{feature}\" Xt.append(date_feature) df2 = pd.concat(Xt, axis=1) return df2 def fit(self, df: pd.DataFrame, y=None, **fit_params): return self

这个转换器相对简单，它为日期列提取一系列特征，例如小时、年份、日期、星期几、月份、年度周数和季度。这些特征在机器学习上下文中，可能对描述或注释时间序列数据非常有用。

你可以在 GitHub 上找到这个示例的完整代码。我提供了一个额外的转换器，用于处理章节中未涉及的周期性特征。

我们按如下方式应用转换器到 DataFrame 的date列：

from sklearn.compose import ColumnTransformerfrom sklearn.pipeline import Pipeline, make_pipelinepreprocessor = ColumnTransformer( transformers=[( \"**date**\", make_pipeline( DateFeatures(), ColumnTransformer(transformers=[ (\"**cyclical**\", CyclicalFeatures(),  [\"**date_day**\", \"**date_weekday**\", \"**date_month**\"] ) ], remainder=\"passthrough\") ), [\"**date**\"], ),], remainder=\"passthrough\")

如果我们希望为预测提供额外的外生特征，可以设置remainder=\"passthrough\"参数。

我们可以定义一个包含这些预处理步骤和模型的管道，这样就可以进行拟合并应用于预测：

from xgboost import XGBRegressorpipeline = Pipeline( [ (\"**preprocessing**\", preprocessor), (\"xgb\", XGBRegressor(objective=\"**reg:squarederror**\", n_estimators=**1000**)) ])

预测器是一个 XGBoost 回归器。我在调整方面并没有做太多工作。唯一会改变的参数是估计器的数量。我们将使用 1,000 个树的集成大小（即树的数量）。

现在是时候将数据集拆分为训练集和测试集了。这包括两个问题：

  • 我们需要提前对特征和数值进行对齐

  • 我们需要根据截止时间将数据集拆分为两部分

首先设定基本参数。首先，我们希望基于时间范围预测未来。其次，我们需要决定用于训练和测试的数据点数量：

TRAIN_SIZE = int(len(df) * **0.9**)HORIZON = **1**TARGET_COL = \"**new_cases**\"

我们使用 90%的数据进行训练，并预测未来 90 天的情况：

X_train, X_test = df.iloc[HORIZON:TRAIN_SIZE], df.iloc[TRAIN_SIZE+HORIZON:]y_train = df.shift(periods=HORIZON).iloc[HORIZON:TRAIN_SIZE][TARGET_COL]y_test = df.shift(periods=HORIZON).iloc[TRAIN_SIZE+HORIZON:][TARGET_COL]

这既进行了对齐，也设定了预测范围。因此，我们有了用于测试和训练的数据集，两个数据集都包含了我们希望用 XGBoost 预测的特征和标签。

现在我们可以训练我们的 XGBoost 回归模型，根据我们通过转换器生成的特征和当前值，预测未来 HORIZON 内的值。

我们可以按如下方式拟合我们的管道：

FEATURE_COLS = [\"date\"]pipeline.fit(X_train[FEATURE_COLS], y_train)

我们可以看到以下管道参数：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_08.png

图 7.8：管道参数

如果我们从开始到结束创建一个日期系列，我们可以获得整个时间段内模型的预测：

MAX_HORIZON = **90**X_test_horizon = pd.Series(pd.date_range( start=df.date.**min()**, periods=**len**(df) + MAX_HORIZON, name=\"**date**\")).reset_index()

应用于X_test的管道的predict()方法为我们提供了预测结果：

forecasted = pd.concat( [pd.Series(pipeline.predict(X_test_horizon[FEATURE_COLS])), pd.Series(X_test_horizon.date)], axis=1)forecasted.columns = [TARGET_COL, \"**date**\"]

我们也可以对实际病例做同样的操作：

actual = pd.concat( [pd.Series(df[TARGET_COL]), pd.Series(df.date)], axis=1)actual.columns = [TARGET_COL, \"**date**\"]

现在，我们可以通过图表将预测结果与实际值y_test进行对比：

fig, ax = plt.subplots(figsize=(12, 6))forecasted.set_index(\"date\").plot(linestyle=\'--\', ax=ax)actual.set_index(\"date\").plot(linestyle=\'-.\', ax=ax)plt.legend([\"forecast\", \"actual\"])

这是我们得到的图表：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_09.png

图 7.9：预测与实际值对比（XGBoost）

我们可以通过以下方式提取测试期间的性能指标：

from sklearn.metrics import mean_squared_errortest_data = actual.merge(forecasted, on=\"**date**\", suffixes=(\"**_actual**\", \"**_predicted**\"))mse = mean_squared_error(test_data.new_cases_actual, test_data.new_cases_predicted, squared=False) # RMSE**print(\"The root mean squared error (RMSE) on test set: {:.2f}\".format(mse))**

我们应该看到类似这样的结果：

The root mean squared error (RMSE) on test set: 12753.41

接下来，我们将在 Kats 中创建一个用于时间序列预测的集成模型。

Kats 的集成方法

Kats 的安装应该很简单，只需两个步骤。首先，我们安装 Facebook 的 Prophet 库的旧版本 fbprophet：

conda install -c conda-forge fbprophet

现在我们使用 pip 安装 Kats：

pip install kats

或者，在 Colab 上，我们可以这样安装 Kats：

!MINIMAL=1 pip install kats!pip install \"numpy==1.20\"

我们将像之前一样加载 COVID 病例数据集。这里只展示最后一行：

df = owid_covid[owid_covid.location == \"**France**\"].set_index(\"**date**\", drop=True).resample(\'**D**\').interpolate(method=\'**linear**\').reset_index()

我们将配置集成模型，拟合它，然后进行预测。

首先是我们的集成模型的配置：

from kats.models.ensemble.ensemble import EnsembleParams, BaseModelParamsfrom kats.models.ensemble.kats_ensemble import KatsEnsemblefrom kats.models import linear_model, quadratic_modelmodel_params = EnsembleParams( [ BaseModelParams(\"linear\", linear_model.LinearModelParams()), BaseModelParams(\"quadratic\", quadratic_model.QuadraticModelParams()), ] )

这里，我们只包括了两种不同的模型，但我们本可以加入更多模型，也可以定义更好的参数。这只是一个示例；对于更现实的练习（我留给读者自己做），我建议加入 ARIMA 和 Theta 模型。我们需要为每个预测模型定义超参数。

我们还需要创建集成参数，定义如何计算集成聚合以及如何进行分解：

KatsEnsembleParam = { \"**models**\": model_params, \"**aggregation**\": \"**weightedavg**\", \"**seasonality_length**\": 30, \"**decomposition_method**\": \"**additive**\",}

要在 Kats 中使用时间序列，我们必须将数据从 DataFrame 或系列转换为 Kats 时间序列对象。我们可以如下转换我们的 COVID 病例数据：

from kats.consts import TimeSeriesDataTARGET_COL = \"new_cases\"df_ts = TimeSeriesData( value=df[TARGET_COL], time=df[\"date\"])

转换的关键是 Kats 能够推断索引的频率。这可以通过pd.infer_freq()进行测试。在我们的案例中，pd.infer_freq(df[\"date\"])应返回D，表示日频率。

现在我们可以创建我们的KatsEnsemble并进行拟合：

m = KatsEnsemble( data=df_ts, params=KatsEnsembleParam).fit()

我们可以使用predict()方法为每个模型获取单独的预测。如果我们想要获取集成输出，必须在predict()之后调用aggregate()：

m.predict(steps=90)m.aggregate()m.plot()plt.ylabel(TARGET_COL)

我们预测了 90 天的未来。这些预测结果作为模型的一部分被存储，因此我们不需要捕捉返回的预测结果。然后，我们可以将每个模型的预测结果进行聚合。再次说明，我们不需要获取返回的 DataFrame，因为它已经存储在模型对象内部（m.fcst_df）。

最终，我们使用 Kats 的便捷函数绘制了聚合后的 DataFrame：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_07_10.png

图 7.10：Kats 集成模型预测

由于我们可以通过更改基础模型参数和添加新模型来调整这个集成模型，因此它为我们提供了大量的改进空间。

是时候总结一下本章我们所学到的内容了。

总结

本章我们讨论了 Python 中流行的时间序列机器学习库。接着，我们讨论并尝试了带有动态时间规整的 k 最近邻算法，用于机器故障分类。我们还谈到了时间序列预测中的验证，并尝试了三种不同的方法来预测 COVID-19 疫情：Silverkite、XGBoost 梯度提升以及 Kats 中的集成模型。

第八章：时间序列的在线学习

在本章中，我们将深入探讨时间序列的在线学习和流数据。在线学习意味着随着新数据的到来，我们不断更新模型。在线学习算法的优势在于，它们能够处理高速度和可能的大规模流数据，并能够适应数据的新分布。

我们将讨论漂移，漂移非常重要，因为机器学习模型的性能可能会因为数据集的变化而受到强烈影响，直到模型变得过时（陈旧）。

我们将讨论什么是在线学习，数据如何变化（漂移），以及自适应学习算法如何结合漂移检测方法来适应这种变化，从而避免性能下降或代价高昂的再训练。

我们将涵盖以下主题：

  • 时间序列的在线学习

    • 在线算法

  • 漂移

    • 漂移检测方法

  • 自适应学习方法

  • Python 实践

我们将从在线学习的讨论开始。

时间序列的在线学习

学习有两种主要场景——在线学习和离线学习。在线学习意味着您在数据流入时逐步拟合模型（流数据）。另一方面，离线学习，即更常见的方式，意味着您有一个从一开始就已知的静态数据集，机器学习算法的参数一次性调整整个数据集（通常将整个数据集加载到内存中或分批处理）。

在线学习有三种主要的应用场景：

  • 大数据

  • 时间限制（例如，实时）

  • 动态环境

通常，在在线学习环境中，您有更多的数据，且非常适合大数据。在线学习可以应用于大数据集，因为在这些数据集上训练整个数据集在计算上不可行。

在线学习的另一个应用场景是在时间限制下进行推断和拟合（例如，实时应用），与离线算法相比，许多在线算法在资源消耗上非常高效。

在线学习的一个常见应用是在时间序列数据上，特别的挑战是时间序列观察的基础生成过程可能随时间变化。这被称为概念漂移。在离线设置中，参数是固定的，而在在线学习中，参数会根据新数据持续调整。因此，在线学习算法可以处理数据的变化，其中一些算法能够处理概念漂移。

下表总结了在线学习和离线学习的一些区别：

离线 在线 需要监控 是的，模型可能变得陈旧（模型将失去性能） 适应变化的数据 再训练成本 昂贵（从头开始） 便宜（增量式） 内存需求 可能需要较高的内存 低 应用领域 图像分类、语音识别等，假设数据是静态的 金融、电商、经济学和医疗保健等，数据是动态变化的 工具 tslearn, sktime, prophet Scikit-Multiflow, River

图 8.1：在线学习与离线学习方法的对比（时间序列）

还有许多不专门针对在线学习的工具，但支持在线学习，例如最受欢迎的深度学习库——PyTorch 和 TensorFlow，这些库的模型本身支持在线学习，数据加载器支持流式场景——通过迭代器，可以按需加载数据。

监督式机器学习问题的流式处理公式可以如下表示：

  1. 数据点 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_001.png 在时间t到达

  2. 在线算法预测标签

  3. 在下一个数据点到达之前，真实标签会被揭示

在批处理设置中，一组n个数据点 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_002.png 会在时间t同时到达，所有n个数据点将在真实标签被揭示之前，由在线模型进行预测，然后才会到达下一批数据点。

我们可以通过 Python 代码片段演示差异，展示在线与离线设置下机器学习的典型模式。你应该对离线学习熟悉，它的表现形式如下：特征X，目标向量y，和模型参数params：

from sklearn import linear_modeloffline_model = linear_model.LogisticRegression(params)offline_model.fit(X, Y)

这应该是从之前的章节中熟悉的内容，比如第七章，时间序列的机器学习模型。为了简化，我们省略了数据加载、预处理、交叉验证和参数调优等问题。

在线学习遵循以下模式：

from river import linear_modelonline_model = linear_model.LogisticRegression(params)for xi, yi in zip(X, y): online_model.learn_one(xi, yi)

这里，我们是逐个数据点喂给模型。再次强调，这只是简化版——我省略了设置参数、加载数据集等内容。

这些代码片段应该清楚地说明主要的区别：一次性在整个数据集上学习（离线）与逐个数据点学习（在线）。

我应该提到在线方法的评估方法：

  • 保留集

  • 预先评估

在保留集方法中，我们可以将当前模型应用于独立的测试集。这种方法在批处理和在线（流式）学习中都很流行，并能提供无偏的性能估计。

在预先评估中，我们在通过数据序列的过程中进行测试。每个新的数据点会先进行测试，然后再进行训练。

在线学习的一个有趣方面是模型选择，也就是如何在一组候选模型中选择最佳模型。我们在第四章，时间序列的机器学习模型中讨论了时间序列模型的模型选择。在在线环境中，模型选择有不同的选项。

在多臂强盗（也称为K 臂强盗）问题中，有限的资源必须在多个竞争选项之间分配，以最大化预期收益。每个选择（“臂”）都有其回报，可以随着时间的推移进行学习。随着时间的推移，我们可以调整对这些臂的偏好，并根据预期回报进行最优选择。同样，通过学习不同分类或回归模型的预期回报，基于多臂强盗的方法可以应用于模型选择。在实践部分，我们将讨论用于模型选择的多臂强盗算法。

在接下来的章节中，我们将更详细地探讨增量方法和漂移。

在线算法

当数据逐渐变得可用，或其大小超出系统内存限制时，增量式机器学习算法（无论是监督学习还是无监督学习）可以在数据的部分上更新参数，而不是从头开始学习。增量学习是指通过不断调整模型来适应新的输入数据。

一些机器学习方法天生支持增量学习。神经网络（如深度学习）、最近邻算法和进化方法（例如遗传算法）都是增量式的，因此可以应用于在线学习环境，其中它们会不断更新。

增量算法可能会随机访问以前的样本或原型（选定的样本）。这些算法，如基于最近邻算法的增量算法，称为具有部分记忆的增量算法。它们的变体适用于周期性漂移场景。

许多著名的机器学习算法都有增量式的变体，如自适应随机森林、自适应 XGBoost 分类器或增量支持向量机。

强化学习和主动学习可以看作是在线学习的类型，因为它们以在线或主动的方式进行工作。我们将在第十一章中讨论强化学习，时间序列的强化学习。

在在线学习中，更新是持续进行的。其核心是运行统计，因此展示如何在在线环境中增量计算均值和方差会很有帮助。

让我们来看一下在线算术平均值和在线方差的公式。至于在线均值，在时间点t更新均值 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_003.png 可以按以下方式进行：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_004.png

其中 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_005.png 是之前更新的次数——有时也写作 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_006.png。

在线方差 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_007.png 可以基于在线均值和运行中的平方和 https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_008.png 进行计算：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_009.pnghttps://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_010.png

离线算法的一个缺点是，它们有时更难实现，并且在掌握库、算法和方法时会有一定的学习曲线。

scikit-learn 是 Python 中机器学习的标准库，但它只有有限的增量算法，主要集中在批量学习模型上。相比之下，有许多专门用于在线学习的库，这些库具有自适应和增量算法，能够覆盖许多使用场景，如不平衡数据集。

来自新西兰怀卡托大学、巴黎电信学院（Télécom ParisTech）和巴黎综合理工学院（École Polytechnique）的研究工程师、学生和机器学习研究者们一直在开发River 库。River 是由两个库合并而成：Creme（作为增量的双关语）和 Scikit-Multiflow。River 包含了许多元方法和集成方法。作为点睛之笔，这些元方法或集成方法中的许多可以使用 scikit-learn 模型作为基础模型。

截至目前，River 库拥有 1700 个星标，并实现了许多无监督和有监督算法。虽然截至目前，River 的文档仍在完善中，但许多功能已经可用，我们将在本章结尾的实践部分中看到。

该图展示了 River 和 Scikit-Multiflow 随着时间的推移的受欢迎程度（以 GitHub 上的星标数为依据）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_01.png

图 8.2：River 和 Scikit-Multiflow 库的星标历史

我们可以看到，尽管 Scikit-Multiflow 稳步上升，但这一增长大多保持平稳。River 在 2019 年超过了 Scikit-Multiflow，并继续获得 GitHub 用户的大量星标。这些星标类似于社交媒体平台上的“点赞”。

本表展示了一些在线算法，其中部分适用于漂移场景：

算法 描述 非常快速决策树（VFDT） 基于少量示例进行分裂的决策树，也叫做 Hoeffding 树。对漂移较为敏感。 极快速决策树（EFDT） 当有足够信心时增量构建树，并在有更好分裂时替换当前分裂。假设数据分布是平稳的。 Learn++.NSE 用于非平稳环境的增量学习分类器集成。

图 8.3：在线机器学习算法——其中一些适用于漂移场景

最具代表性的在线算法是Hoeffding 树（Geoff Hulten、Laurie Spencer 和 Pedro Domingos，2001），也叫做非常快速决策树（VFDT）。它是最广泛使用的在线决策树归纳算法之一。

尽管一些在线学习算法相对高效，但其性能可能对数据点的顺序极为敏感，且它们可能永远无法摆脱由早期样本驱动的局部最小值。令人吸引的是，VFDT（非常快的决策树）提供了高分类精度，并且具有理论保证，随着时间的推移，它们将趋向于决策树的表现。事实上，VFDT 和传统训练树在树分裂上的差异的概率会随着样本数量的增加而指数下降。

Hoeffding 界限，由 Wassily Hoeffding 于 1963 年提出，表示以概率https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_011.png，随机变量Z的计算均值https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_012.png在n个样本中计算后，与真实均值https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_014.png的偏差小于https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_013.png：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_015.png

在这个方程中，R是随机变量Z的范围。这个界限与生成观测值的概率分布无关。

随着数据的输入，Hoeffding 树会不断添加新的分支，并且淘汰过时的分支。然而，问题在于，在概念漂移的情况下，某些节点可能不再满足 Hoeffding 边界。

在接下来的章节中，我们将讨论漂移，为什么你需要关注它，以及如何处理漂移。

漂移

数据质量的一个主要决定因素是漂移。漂移（也称为：数据集漂移）意味着数据中的模式随着时间的推移发生变化。漂移很重要，因为机器学习模型的表现可能会因为数据集的变化而受到不利影响。

漂移过渡可以突如其来，也可以逐渐发生、增量式发生，或者是周期性发生。下面是一个示例：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_02.png

图 8.4：四种类型的概念漂移过渡

当过渡突然而至时，它会从一个时间步骤跳跃到另一个时间步骤，没有明显的准备或警告。与此相对，它也可能是逐步的，首先是小幅变化，然后是更大的变化，再然后是更大的变化。

当过渡逐渐发生时，它可能表现为不同力量之间的来回波动，直到建立一个新的基线。另一种过渡类型是周期性的，当存在不同基线之间的规律性或重复性的变化时。

漂移有不同种类：

  • 协变量漂移

  • 先验概率漂移

  • 概念漂移

协变量漂移描述了自变量（特征）的变化。一个例子可能是法规干预，其中新的法律会震动市场格局，消费者行为也会与之前不同。例如，如果我们想预测吸烟行为下的慢性疾病风险，并且吸烟变得不再普遍，因为有了新的法律，这意味着我们的预测可能会变得不那么可靠。

概率漂移是目标变量的变化。例如，在欺诈检测中，欺诈发生的比例发生变化；在零售中，商品的平均价值增加。漂移的一个原因可能是季节性因素——例如，在冬季卖出更多的外套。

在概念漂移中，自变量与目标变量之间的关系发生了变化。这个术语所指的概念是自变量与因变量之间的关系。例如，如果我们想预测吸烟的数量，我们可以假设，在新法律出台后，我们的模型将变得无效。请注意，通常“概念漂移”一词在更广泛的意义上应用，指的是任何非平稳的变化。

协变量漂移：特征*P(x)*的变化。

标签漂移（或先验概率漂移）：目标变量*P(y)*的变化。

概念漂移：（在有监督的机器学习中）目标条件分布的变化——换句话说，自变量与因变量之间的关系发生了变化P(y|X)。

通常，在构建机器学习模型时，我们假设数据集不同部分中的数据点属于相同的分布。

尽管偶尔出现的异常现象，如异常事件，通常会被视为噪声并被忽略，但当分布发生变化时，模型通常需要基于新的样本从头开始重建，以捕捉最新的特征。这就是我们在第七章《时间序列的机器学习模型》中讨论的使用前向验证测试时间序列模型的原因。然而，从头开始训练可能非常耗时并且需要大量的计算资源。

漂移给机器学习模型带来了问题，因为模型可能会变得陈旧——它们随着时间推移变得不可靠，因为它们捕捉到的关系不再有效。这导致了这些模型的性能下降。因此，预测、分类、回归或异常检测的方法应能够及时检测并应对概念漂移，以便尽快更新模型。机器学习模型通常会定期重新训练，以避免性能下降。或者，也可以根据模型的性能监控或基于变化检测方法在需要时触发重新训练。

对于时间序列的应用，在许多领域，如金融、电商、经济学和医疗健康，时间序列的统计特性可能会发生变化，从而使得预测模型变得无用。令人困惑的是，尽管漂移问题的概念在文献中已有充分研究，但在使用时间序列方法应对这一问题方面，投入的努力却很少。

Gustavo Oliveira 等人于 2017 年提出（“《在概念漂移下的时间序列预测：一种基于 PSO 的方法》”）训练多个时间序列预测模型。在每个时间点，这些模型的参数会根据最新的性能加权变化（粒子群优化）。当最佳模型（最佳粒子）超出某个置信区间时，便触发了模型的重新训练。

下图展示了错误触发的重新训练和在线学习相结合的一种时间序列预测方法：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_03.png

图 8.5：时间序列预测的在线学习与重新训练（IDPSO-ELM-S）

你可以看到随着概念漂移的发生，误差率周期性地增加，并且根据漂移检测的概念，触发了重新训练。

许多在线模型已经专门适应了对概念漂移的鲁棒性或处理能力。在本节中，我们将讨论一些最流行或表现最佳的模型。我们还将讨论漂移检测方法。

漂移检测方法

有许多不同的方法可以显式地检测数据流中的漂移和分布变化。Page-Hinkley（Page，1954）和几何移动平均（Roberts，2000）是其中的先驱者。

漂移检测器通常通过性能指标来监控模型性能，但它们也可以基于输入特征，尽管这更多的是一种例外。基本思想是，当样本的类别分布发生变化时，模型不再与当前分布相对应，性能下降（误差率增加）。因此，模型性能的质量控制可以作为漂移检测的依据。

漂移检测方法可以至少分为三类（根据 João Gama 等人，2014 年）：

  • 统计过程控制

  • 顺序分析

  • 基于窗口的比较

统计过程控制方法考虑了模型预测的汇总统计量，如均值和标准差。例如，漂移检测方法（DDM；João Gama 等人，2004）会在误差率超过之前记录的最小误差率三倍标准差时发出警报。根据统计学习理论，在持续训练的模型中，随着样本数量的增加，错误应该减少，因此只有在发生漂移的情况下，这个阈值才会被超过。

顺序方法基于模型预测的阈值。例如，在线性四率（Wang，2015）方法中，列联表中的比率会递增更新。显著性是根据一个在开始时通过蒙特卡洛抽样估算的阈值来计算的。该方法比 DDM 更能处理类别不平衡。

列联表：比较变量频率分布的表格。具体来说，在机器学习分类中，表格显示了预测标签在测试集上的数量与实际标签的对比。在二分类的情况下，单元格显示真正例、假正例、假反例和真反例。

基于窗口的方法监控错误的分布。例如，ADWIN (自适应滑动窗口) 是由 Albert Bifet 和 Ricard Gavaldà于 2007 年提出的。时间窗口W内的预测错误被划分成更小的窗口，并将这些窗口内的平均误差率差异与 Hoeffding 界限进行比较。原始版本提出了一种变体，其时间复杂度为O(log W)，其中W是窗口的长度。

这里列出了一些漂移检测方法：

算法 描述 类型 自适应窗口法 (ADWIN) 基于阈值的自适应滑动窗口算法。 基于窗口 漂移检测方法 (DDM) 基于模型错误率应随时间减少的前提。 统计方法 早期漂移检测法 (EDDM) 基于两个错误之间的平均距离统计。与 DDM 类似，但更适用于渐进式漂移。 统计方法 霍夫丁漂移检测 (HDDM) 基于 Hoeffding 界限的非参数方法——移动平均检验或加权移动平均检验。 基于窗口 Kolmogorov-Smirnov 窗口法 (KSWIN) 在时间序列的窗口中进行 Kolmogorov-Smirnov 检验。 基于窗口 Page-Hinkley 用于高斯信号均值变化的统计检验。 序列化方法

图 8.6：漂移检测算法

Kolmogorov-Smirnov 是一个非参数检验，用于检验连续一维概率分布的相等性。

这些方法可以用于回归和分类（以及预测）场景中。它们可以用来触发模型的重新训练。例如，Hassan Mehmood 等人（2021）在检测到漂移时，会重新训练时间序列预测模型（其中包括 Facebook 的 Prophet 模型）。

漂移检测器都有其关于输入数据的假设。了解这些假设非常重要，我尝试在表格中概述了这些假设，以便你能选择适合你的数据集的检测器。

上述列出的漂移检测方法都有标签成本。由于它们都监控基本分类器或集成分类器的预测结果，因此要求在预测后立即获得类标签。在某些实际问题中，这一限制是不切实际的。这里没有列出其他一些方法，这些方法可以基于异常检测（或新颖性检测）、特征分布监控或模型依赖监控。我们在第六章，时间序列的无监督方法中看到了一些这些方法。

在下一节中，我们将介绍一些旨在抵抗漂移的方法。

自适应学习方法

自适应学习是指具有漂移调整的增量方法。这个概念指的是通过在线更新预测模型，以应对概念漂移。目标是通过考虑漂移，使得模型能够确保与当前数据分布的一致性。

集成方法可以与漂移检测器结合使用，以触发基础模型的重新训练。它们可以监控基础模型的表现（通常使用 ADWIN）——表现不佳的模型会被重新训练过的更精确模型替换。

作为一个例子，自适应 XGBoost 算法 (AXGB；Jacob Montiel 等人，2020 年) 是 XGBoost 在处理不断变化的数据流时的适应性改编，其中新的子树是从数据的小批量中创建的，随着新数据的到来。这种算法的最大集成大小是固定的，一旦达到该大小，集成就会在新数据上进行更新。

在 Scikit-Multiflow 和 River 库中，有几种方法将机器学习方法与漂移检测方法结合，这些方法调节适应性。这些方法中的许多都是由这两个库的维护者发布的。以下是其中一些方法的列表：

算法 描述 K-近邻 (KNN) 分类器与 ADWIN 漂移检测器 使用 ADWIN 漂移检测器的 KNN，决定保留或遗忘哪些样本。 自适应随机森林 每棵树都包含漂移检测器。它在检测到警告后开始在后台训练，如果发生漂移，则替换旧的树。 加法专家集成分类器 实现了修剪策略 —— 最旧或最弱的基础模型将被移除。 Hoeffding 自适应树 (HAT) 将 ADWIN 与 Hoeffding 树模型结合，检测漂移并进行学习。 非常快速的决策规则 类似于 VFDT，但使用规则集成而不是树。在 Scikit-Multiflow 中支持通过 ADWIN、DDM 和 EDDM 进行漂移检测。 Oza Bagging ADWIN 不同于有放回采样，每个样本都赋予一个权重。在 River 中，可以将其与 ADWIN 漂移检测器结合使用。 在线 CSB2 一种在线提升算法，折衷了 AdaBoost 和 AdaC2，并可选择使用漂移检测器。 在线提升 带有 ADWIN 漂移检测的 AdaBoost。

图 8.7：自适应学习算法

这些方法通过调节适应性或学习，以漂移检测的概念来应对漂移。

让我们试试这些方法中的一些！

Python 实践

本章的安装非常简单，因为在本章中，我们只会使用 River。我们可以从终端快速安装它（或通过 Anaconda Navigator 安装）：

pip install river

我们将在 Python（或 IPython）终端执行这些命令，但同样，我们也可以在 Jupyter Notebook（或其他环境）中执行它们。

漂移检测

让我们从尝试使用人工时间序列进行漂移检测开始。这个例子来自 River 库的测试。

我们将首先创建一个可以测试的人工时间序列：

import numpy as npnp.random.seed(12345)data_stream = np.concatenate( (np.random.randint(2, size=1000), np.random.randint(8, size=1000)))

这个时间序列由两组具有不同特征的序列组成。让我们看看漂移检测算法多快能够检测到这一点。

在这个数据集上运行漂移检测器意味着我们需要遍历数据集，并将值输入到漂移检测器中。我们将为此创建一个函数：

def perform_test(drift_detector, data_stream): detected_indices = [] for i, val in enumerate(data_stream): in_drift, in_warning = drift_detector.update(val) if in_drift: detected_indices.append(i) return detected_indices

现在我们可以在这个时间序列上尝试 ADWIN 漂移检测方法。让我们创建另一个方法，绘制漂移点与时间序列的叠加图：

import matplotlib.pyplot as pltdef show_drift(data_stream, indices): fig, ax = plt.subplots(figsize=(16, 6)) ax.plot(data_stream) ax.plot( indices, data_stream[indices], \"ro\", alpha=0.6, marker=r\'$\\circ$\', markersize=22, linewidth=4 )plt.tight_layout()

这是 ADWIN 漂移点的图示：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_04.png

图 8.9：我们人工数据集上的 ADWIN 漂移点

我鼓励你尝试一下，并尝试其他漂移检测方法。

接下来，我们将进行回归任务。

回归

我们将估计中等强度太阳耀斑的发生。

为此，我们将使用来自 UCI 机器学习库的太阳耀斑数据集。River 库附带了该数据集的压缩列分隔数据集，我们将加载它，指定列类型，并选择我们感兴趣的输出。

现在让我们绘制 ADWIN 结果：

from river import streamfrom river.datasets import baseclass SolarFlare(base.FileDataset): def __init__(self): super().__init__( n_samples=1066, n_features=10, n_outputs=1, task=base.MO_REG, filename=\"solar-flare.csv.zip\", ) def __iter__(self): return stream.iter_csv( self.path, target=\"m-class-flares\", converters={ \"zurich-class\": str, \"largest-spot-size\": str, \"spot-distribution\": str, \"activity\": int, \"evolution\": int, \"previous-24h-flare-activity\": int, \"hist-complex\": int, \"hist-complex-this-pass\": int, \"area\": int, \"largest-spot-area\": int, \"c-class-flares\": int, \"m-class-flares\": int, \"x-class-flares\": int, }, )

请注意，我们如何选择目标数量和转换器，这些转换器包含所有特征列的类型。

让我们看看这会是什么样子：

from pprint import pprintfrom river import datasetsfor x, y in SolarFlare(): pprint(x) pprint(y) break

我们看到数据集的第一个点（数据集的第一行）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_05.png

图 8.10：中等强度太阳耀斑数据集的第一个点

我们看到十个特征列作为字典，输出为浮动值。

让我们在 River 中构建我们的模型管道：

import numbersfrom river import composefrom river import preprocessingfrom river import treenum = compose.SelectType(numbers.Number) | preprocessing.MinMaxScaler()cat = compose.SelectType(str) | preprocessing.OneHotEncoder(sparse=False)model = tree.HoeffdingTreeRegressor()pipeline = (num + cat) | model

这样的管道非常易于阅读：数值特征进行最小-最大缩放，而字符串特征则进行独热编码。预处理后的特征输入到 Hoeffding 树模型中进行回归。

现在我们可以按预序列化的方式训练我们的模型，预测值并进行训练，就像之前讨论过的那样：

from river import evaluatefrom river import metricsmetric = metrics.MAE()evaluate.progressive_val_score(SolarFlare(), pipeline, metric)

我们使用平均绝对误差（MAE）作为我们的评估指标。

我们得到的 MAE 为 0.096979。

这个预序列化评估evaluate.progressive_val_score()等同于以下内容：

errors = []for x, y in SolarFlare(): y_pred = pipeline.predict_one(x) metric = metric.update(y, y_pred) errors.append(metric.get()) pipeline = pipeline.learn_one(x, y)

我已添加了两行代码来收集算法学习过程中的误差。

让我们绘制一下这个图：

fig, ax = plt.subplots(figsize=(16, 6))ax.plot( errors, \"ro\", alpha=0.6, markersize=2, linewidth=4)ax.set_xlabel(\"number of points\")ax.set_ylabel(\"MAE\")

这张图显示了这个误差如何随着算法遇到的点数而变化：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_06.png

图 8.11：根据点数计算的 MAE

我们可以看到，在 20 到 30 个点之后，度量稳定后，Hoeffding 树开始学习，误差不断下降，直到大约 800 个点，之后误差再次增加。这可能是行排序效应。

概念漂移的数据集是适应性模型的典型应用场景。让我们在一个有概念漂移的数据集上比较适应性和非适应性模型：

from river import ( synth, ensemble, tree, evaluate, metrics)models = [ tree.HoeffdingTreeRegressor(), tree.HoeffdingAdaptiveTreeRegressor(), ensemble.AdaptiveRandomForestRegressor(seed=42)]

我们将比较 Hoeffding 树回归器、适应性 Hoeffding 树回归器和适应性随机森林回归器。我们将采用每个模型的默认设置。

我们可以使用一个合成数据集来进行此测试。我们可以在数据流上训练前述的每个模型，并查看均方误差（MSE）指标：

for model in models: metric = metrics.MSE() dataset = synth.ConceptDriftStream( seed=42, position=500, width=40 ).take(1000) evaluate.progressive_val_score(dataset, model, metric) print(f\"{str(model.__class__).split(\'.\')[-1][:-2]}: {metric.get():e}\")

evaluate.progressive_val_score方法遍历数据集的每个点并更新度量。我们得到如下结果：

HoeffdingTreeRegressor: 8.427388e+42HoeffdingAdaptiveTreeRegressor: 8.203782e+42 AdaptiveRandomForestRegressor: 1.659533037987239+42

由于这些算法的性质，你的结果可能会有所不同。我们可以设置随机数生成器的种子来避免这种情况，然而，我认为强调这一点是值得的。

我们看到模型误差（MSE）以科学计数法表示，这有助于理解这些数字，因为它们相当大。你会看到误差分为两部分，首先是因子，然后是以 10 为指数的数量级。三个模型的数量级相同，然而，适应性随机森林回归器的误差大约只有其他两个模型的五分之一。

我们还可以通过可视化误差随时间的变化，观察模型的学习和适应过程：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_07.png

图 8.12：概念漂移数据流的模型表现（MSE）

River 中没有非适应性版本的随机森林算法，所以我们无法进行比较。我们无法得出适应性算法是否真正更有效的明确结论。

如果你想尝试不同的模型、元模型和预处理器，还有很多其他选项可以尝试。

模型选择

我们在本章早些时候提到过使用多臂赌博机进行模型选择，这里我们将通过一个实际示例来演示。这基于 River 文档中的内容。

让我们使用UCBRegressor来选择线性回归模型的最佳学习率。相同的模式可以更广泛地应用于选择任何一组（在线）回归模型。

首先，我们定义模型：

from river import composefrom river import linear_modelfrom river import preprocessingfrom river import optimmodels = [ compose.Pipeline( preprocessing.StandardScaler(), linear_model.LinearRegression(optimizer=optim.SGD(lr=lr)) ) for lr in [1e-4, 1e-3, 1e-2, 1e-1]]

我们在 TrumpApproval 数据集上构建并评估我们的模型：

from river import datasetsdataset = datasets.TrumpApproval()

我们将应用 UCB 赌博机算法，该算法计算回归模型的奖励：

from river.expert import UCBRegressorbandit = UCBRegressor(models=models, seed=1)

该赌博机提供了以在线方式训练模型的方法：

for x, y in dataset: bandit = bandit.learn_one(x=x, y=y)

我们可以检查每个“臂”被拉动的次数（百分比）。

for model, pct in zip(bandit.models, bandit.percentage_pulled): lr = model[\"LinearRegression\"].optimizer.learning_rate print(f\"{lr:.1e} — {pct:.2%}\")

四个模型的百分比如下：

1.0e-04 — 2.45%1.0e-03 — 2.45%1.0e-02 — 92.25%1.0e-01 — 2.85%

我们还可以查看每个模型的平均奖励：

for model, avg in zip(bandit.models, bandit.average_reward): lr = model[\"LinearRegression\"].optimizer.learning_rate print(f\"{lr:.1e} — {avg:.2f}\")

奖励如下：

1.0e-04 — 0.001.0e-03 — 0.001.0e-02 — 0.741.0e-01 — 0.05

我们还可以绘制奖励随时间的变化，并根据模型表现进行更新：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_08.png

图 8.13：奖励随时间的变化

你可以看到，随着我们逐步处理数据并更新模型，奖励逐渐变得明显。模型奖励在大约 100 个时间步长时明显分开，而在大约 1000 个时间步长时，似乎已经收敛。

我们还可以绘制不同模型在每个步骤中被选择的百分比（这是基于奖励的）：

https://github.com/OpenDocCN/freelearn-ml-zh/raw/master/docs/ml-ts-py/img/B17577_08_09.png

图 8.14：随着时间推移选择的模型比例

这个分布大致跟随随时间变化的奖励分布。这是可以预期的，因为模型选择依赖于奖励（以及一个调节探索的随机数）。

我们还可以选择最佳模型（即具有最高平均奖励的模型）。

best_model = bandit.best_model

赌博机选择的学习率是：

best_model[\"LinearRegression\"].intercept_lr.learning_rate

学习率为 0.01。

总结

在本章中，我们讨论了在线学习。我们谈到了在线学习方法的一些优点：

  • 它们高效，能够处理高速吞吐量。

  • 它们可以处理非常大的数据集。

  • 它们能够适应数据分布的变化。

概念漂移是数据与目标之间关系的变化。我们已经讨论了漂移的重要性，机器学习模型的性能可能会受到数据集变化的强烈影响，甚至导致模型过时（陈旧）。

漂移检测器不直接监控数据本身，而是用于监控模型性能。漂移检测器可以使流学习方法对概念漂移具有鲁棒性，在 River 中，许多自适应模型使用漂移检测器进行部分重置或调整学习参数。自适应模型是将漂移检测方法结合使用的算法，旨在避免性能退化或避免高成本的重新训练。我们概述了几种自适应学习算法。

在 Python 实践中，我们尝试了 River 库中的一些算法，包括漂移检测、回归和使用多臂赌博机方法的模型选择。

网络标签：

相关问题