TowardsDataScience 博客中文翻译 2019（二百零一）

原文：TowardsDataScience Blog

协议：CC BY-NC-SA 4.0

用时尚愚弄面部检测

原文：https://towardsdatascience.com/fooling-facial-detection-with-fashion-d668ed919eb?source=collection_archive---------12-----------------------

Photo by Pavel Anoshin on Unsplash

Accompanying GitHub repository: [https://github.com/BruceMacD/Adversarial-Faces](https://github.com/BruceMacD/Adversarial-Faces)

面部识别的使用正在增加。随着最近关于面部识别伦理的辩论，我一直在考虑对面部识别的潜在敌对攻击。从 T2 机场到社交媒体，面部识别正被广泛应用。选择不接受面部扫描似乎几乎是不可能的。

对面部检测的理想攻击应该是一件在不知情的人看来不显眼的衣服。受 Hyperface 项目的启发，我决定研究并实现一个可穿戴的对抗实例。在这篇文章中，我将详细介绍创建一个敌对图像来欺骗选定类型的面部检测的过程，以及我如何在一个面具上实现一个实际的例子。

面部检测与面部识别

An illustration of facial detection (left) vs. facial recognition (right)

在深入研究这个项目之前，首先需要注意的是面部检测和面部识别之间的区别。面部检测指的是检测图像中何时出现面部的能力。面部识别依靠面部检测来确定图像中是否存在人脸，但它更进一步，试图确定这是谁的脸。

对于这个项目，我选择把重点放在面部检测上。主要是因为它更容易测试。为了正确地测试面部识别，访问面部识别数据库将是理想的。

面部检测模型

下一步是选择哪种面部检测模型来建立对抗的例子。目前有许多不同的面部检测模型在使用。Vikas Gupta 在“Learn OpenCV”上有一篇关于面部检测模型及其实现的很好的入门文章，其中有深入的解释。我将在这里简单回顾一下。

[## 人脸检测- OpenCV、Dlib 和深度学习|学习 OpenCV

在本教程中，我们将讨论 OpenCV 和 Dlib 中的各种人脸检测方法，并对这些方法进行比较

www.learnopencv.com](https://www.learnopencv.com/face-detection-opencv-dlib-and-deep-learning-c-python/)

• 深度神经网络(DNNs): DNNs 可以使用输入数据集进行训练，以检测许多不同方向的人脸。基于 DNN 的面部检测的一种流行方法是单镜头多盒检测器。dnn 是精确和通用的。
• **卷积神经网络(CNN)😗*CNN 是一种深度神经网络，旨在为图像的不同部分分配重要性。它是健壮的，但是在 CPU 上相当慢。
• **哈尔级联分类器:**哈尔级联使用具有大量标记的阳性和阴性示例图像的数据集来训练。Haar 级联分类器的主要缺点是它们只能正面识别人脸。它们不再被广泛使用，因为神经网络更加通用。
• **梯度方向直方图(HOG)😗*HOG 是一种面部检测方法，它将经过处理的输入图像分成具有梯度方向的单元，然后将结果输入到支持向量机中。猪检测是快速和轻量级的，但对一些不寻常的面部角度不起作用。

A modeled example of a face as a histogram of oriented gradients from dlib

很快脱颖而出成为最简单的攻击候选模型的是梯度方向直方图。最值得注意的是，猪的预期输入可以很容易地被可视化并反馈到面部检测模型中。将面部可视化为定向梯度的直方图还具有对于人类观察者来说不是明显的面部的优点。

Python 中的梯度方向直方图人脸检测

**Note:** Expanded code samples with the functionality to display results are available on the accompanying [GitHub repository](https://github.com/BruceMacD/Adversarial-Faces).

为了测试这些例子，我需要一个简单的基于 HOG 的面部检测实现。幸运的是，dlib 库在其frontier _ face _ detector中内置了一个猪面部检测器。

import dlib
import cv2cv2.imread("path/to/input_img.png")
frontal_face_detector = dlib.get_frontal_face_detector()
upscaling_factor = 1
detected_faces = frontal_face_detector(img, upscaling_factor)

用输入图像和放大因子运行正面人脸检测器。放大因子 1 表示输入图像将被放大一次。放大可以创建更大的图像，从而更容易检测人脸。正面人脸检测的结果是一个边界框列表，每个检测到的人脸对应一个边界框。

The result of using the dlib to detect a face in our visualized HOG

通过猪的预期输入的可视化，您可以看到它被检测为一张脸。太好了！我们有了对抗性进攻的基础。

使用随机优化创建对抗性设计

现在我知道猪的预期输入的可视化将被检测为正面脸的假阳性，我需要创建一个在看起来不显眼的面具上打印的设计。然而，仍然有许多影响设计的因素，我不知道如何优化。人脸的位置、方向和大小都会影响图像中检测到的人脸数量。我可以简单地尝试不同的设计，直到我找到一个好的，但是让一个学习模型为我做艰苦的工作似乎更有趣，也不那么乏味。

我考虑了几种不同的模型来寻找最佳输入。我研究了强化学习、生成对抗网络和 Q 学习。最终，我决定使用随机优化的模拟退火，因为它最适合我的问题，即找到与 dlib 检测到的大多数人脸相对应的输入。

我使用 PIL(Python 图像库)和 mlrose(用于随机优化的 Python 库)来生成图像并找到最佳状态。用 mlrose 优化需要一个初始状态和一个适应度函数。在我的情况下，找到这个最佳状态是一个非常昂贵的计算，因为生成的状态需要作为图像保存到磁盘，以便找到检测到的人脸数量。

# indexes:
# 0 % 4 = pos_x
# 1 % 4 = pos_y
# 2 % 4 = rotation
# 3 % 4 = scale
initial_state = np.array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0])

从初始状态开始，mlrose 需要一个 1D 数组(据我所知)。这意味着我不得不使用一个简单的解决方案，给不同的数组位置赋予不同的意义(参见索引解释)。我选择优化 6 个面的输入，因为我总是可以复制设计来增加它的大小。

def detected_max(state):
    # converts the 1D state array into images
    get_img_from_state(state)
    return len(detect_faces(cv2.imread(OUTPUT)))

我的适应度函数只是由状态到图像的转换组成，然后检测图像中的人脸数量。找到的人脸数量越多，拟合度越好。我还试着根据输入猪脸图像的大小将适应度函数修改得更高。这可能更好，因为在现实生活中，更大的人脸更有可能被检测到。然而，我发现在视觉效果相似的情况下，考虑面部大小会导致计算时间延长。

fitness = mlrose.CustomFitness(detected_max)
problem = mlrose.DiscreteOpt(length=24, fitness_fn=fitness,
                             maximize=True, max_val=scale_factor)
schedule = mlrose.ExpDecay()
best_state, max_faces = mlrose.simulated_annealing(problem, schedule=schedule, max_attempts=10, max_iters=1000,
                                          init_state=initial_state, random_state=1)

print('Optimal state found: ', best_state)
print('Max fitness found: ', max_faces)
# save the optimal found
get_img_from_state(best_state)
print("Number of faces in output: ", len(detect_faces(cv2.imread(OUTPUT))))

有了适应度和初始状态集，为模拟退火配置 mlrose 就很简单了。我只是分配我们的输入，并让它运行，直到找到一个最佳结果。我运行了几次，找到了一个视觉上有趣的结果。

An interesting output from the simulated annealing

最后，随着这个有趣的输出，我添加了一些收尾工作来模糊它的面部设计。我决定我更有资格用手来做这件事，因为我的意图是愚弄人类。

The final design that obscures the facial structure

测试掩模上的设计

Faces detected on a prototype mask covered in the adversarial face design

随着最终设计的完成，我创造了一些模拟面具设计来测试他们是如何被猪面部检测评估的。最初的结果似乎很有希望。上述设计始终返回 4–5 个错误检测的人脸。

足球转会 Python 美丽的汤 15 年分析

原文：https://towardsdatascience.com/football-wages-15-year-analysis-with-pythons-beautiful-soup-f3a7b8882524?source=collection_archive---------10-----------------------

一个多页网页抓取教程使用美丽的汤

Image Courtesy of Nathan Dumlao via Unsplash

介绍

Web 抓取能够自动收集丰富的数据集。如果您可以在 web 浏览器中查看一些数据，您将能够通过程序访问和检索这些数据。如果你可以通过一个程序访问它，这些数据就可以以任何方式存储、清理和使用。

Web scarping 提供了一些优于应用程序编程接口(API)的优势。网络抓取是免费的，没有速率限制，让你接触到所有你想获取的数据，最重要的是，你想从中提取数据的网站可能没有 API。这就是网络抓取进入画面的时候。

本教程将探讨如何编写一个简单的 Web scraper，它可以收集过去 15 年英超联赛中每个球员的平均转会收入数据。

入门指南

首先，需要导航至 transfermarkt 。下一步是查看这个页面的底层 HTML，在 Chrome 中右键单击并选择“inspect”。

这里特别有用的是，你可以将鼠标悬停在 Elements 标签中的 HTML 标签上，Chrome 会在网页本身的表示上描绘一个透明框。这可以快速帮助我们找到我们正在搜索的内容。

或者，您也可以右键单击网页上的任何元素，然后单击“检查元素”。这将立即在 Elements 选项卡中突出显示相应的 HTML 代码。

在这里，我想提炼一下:收入，每个俱乐部的收入，每个球员的收入。为此，我通过将鼠标悬停在带有“ transferbilanz ”类的“div”标签上来找到包装器。

美丽的汤库

为了解析来自网页的收入信息，我使用了漂亮的 Soup 库。美丽的汤很容易安装使用画中画。

信息所在页面的 url 放在名为 URL 的变量的引号中。requests.get()方法用于“获取”URL 页面，结果保存在一个响应变量中(记住安装请求库并将其导入脚本)。

为了确认连接已经建立，我使用了**。status_code** 属性，并检查服务器是否返回了超文本传输协议(HTTP)状态代码 200。这意味着“成功”，页面已被检索。

包含当前页面的原始 HTML 内容的 response.text 被保存到变量“financial_data”中。现在，我创建一个 soup BeautifulSoup 对象，将原始 HTML 数据 financial_data 作为第一个参数传递，将“html.parser”作为第二个参数传递。这个解析器内置在 Python 中，不需要安装，用于解析 HTML 标签。

一旦我知道成功连接到网页，我就可以利用 URL 的结构。

在下面显示的 URL 中，每个足球赛季唯一变化的部分是年份。当我第一次检查连接时，年份是 2018 年(下面用红色和粗体格式突出显示)。为了检索过去 15 个赛季的相同信息，我只需将赛季编号更改为 2004！

我创建了一个简单的 for 循环来遍历从 2004 年到 2018 年的一系列年份，以检索每一年的收入数字。(记住最后一个数字 2019 被排除在名单之外

在这个 for 循环中，我需要的三条信息可以在 span 标记中找到，带有如下所示的“greentext”类。

为了检索这些数据，我使用汤找到了包含所有这些信息的第一个包装器。查找方法并将结果保存到名为 grouped_data 的变量中。

然后，我对 grouped_data 使用 find_all 方法来查找所有带有“greentext”类的 span 标签。由于 find_all 方法返回一个类似列表的对象，我可以索引这个列表来检索我需要的值，即列表中第一个元素的[0]指向该季节的“收入”(参见下面的 github 要点)。

我需要删除数字中的空格、欧元符号和句点，以便将值转换为浮点数，以便稍后进行基于数字的分析。我还将返回的每个值乘以 0.89，因为这是欧元和英镑之间的当前(2019 年 13 月 6 日)汇率，我希望数据框架中的结果以英镑为单位。

对于每一年(我的 for 循环的每一次迭代)，我将数据附加到适当标记的列表中，即 income_per_player_list。

一旦所有迭代完成，我使用 pd 创建一个数据帧。DataFrame 方法(需要将 pandas 模块导入到脚本中)，并写入 CSV 文件。

然后我读入数据帧；

Income 2004 _ 18 _ df = PD . read _ CSV(’ Income _ euro _ to _ pounds . CSV ')

这使我们能够开始得出结论:

(使用 Python 的 Seaborn 库生成的条形图)

我们可以看到，除了 2018 年，每个赛季每个球员的转会收入都在稳步攀升。有趣的是，这一变化反映了 2016 年签署的英超电视转播权协议。现在游戏中的钱比以往任何时候都多。这使得支付给球员的转会费更高，足球俱乐部每个赛季每个球员的平均收入也相应增加！

天空电视台和英国电信体育公司为 2016-17 赛季三个赛季的英超联赛电视转播权支付了创纪录的 51.36 亿英镑。

结论

本示例教程演示了如何执行多页面网络抓取，并将数据转换为适合分析的 Pandas DataFrame。

本质上，我们已经扔掉了我们的网络浏览器，使用 Python 程序在网上冲浪来提取我们需要的信息！

“足球天气”——探究天气对 NFL QB 表现的影响

原文：https://towardsdatascience.com/football-weather-diving-into-the-effects-of-weather-on-nfl-qb-performance-f0edb420623d?source=collection_archive---------7-----------------------

作为一个 NFL 的狂热粉丝，也许是一个更大的梦幻足球粉丝，我的梦幻足球爱好中一直困扰我的一个挫折来源是我不确定恶劣天气预报会对比赛和比赛中球员的统计产量产生多大影响。我经常发现自己对雨和风的预报反应过度，不管天气如何，坐冷板凳只会错过一场重要的比赛。这仅仅是小样本/坏运气，还是说“坏”天气实际上对一个专业运动员来说不是一个有意义的损害？如果 10 英里每小时的大风和小雨没什么大不了的，那么下雪呢？时速 15 英里的风怎么样？20?温度呢？在开始产生有意义的影响之前，温度需要达到哪个极端？这些是我在这篇博文中试图回答的问题，或者至少是阐明了一些问题。

找到 NFL 比赛的历史天气数据出奇地困难，但幸运的是，通过一些相当简单的 python 代码，我能够从profootballarchives.com收集到 1985 年至 2016 年间几乎每场比赛的基本天气信息。我也能够找到一些 NFL 的统计数据，尽管不幸的是这并不完全是按游戏的比分格式。从 Kaggle NFL 数据集，我决定使用 QB 游戏日志数据。经过一点争论，我能够将数据集合并到一个熊猫数据框架中。

在进行了相当多的清理、重命名列、对变量进行分类、宁滨一些连续数据以及转换/重铸数据类型之后，dataframe 已经准备好进行一些分析了。首先，由于四分卫通过空中传球获得报酬，我认为查看不同风速的数据将是一个好的开始。我将数据分为 4 类:0-10 英里，10-15 英里，15-20 英里和 20+英里。让我们来看看风对 QB 性能的几个测量的影响:

随着每个子图的风速类别从左到右增加，QB 性能有一个非常明显的下降趋势。随着游戏平均风速的增加，生产指标(码数、TD 传球数和幻想点数)以及效率指标(每次尝试的码数、QB 评分和完成百分比)都会受到影响。

接下来，我分析了降水的影响，分为 5 类:无(或室内/圆顶)、小雨、小雪、雨/风暴和雪/风暴。

奇怪的是，小雨似乎比暴雨有更深远的负面影响。然而，小雪、雨/风暴和雪的误差线相当大，这表明该数据可能不可靠。我对显著性进行了一些统计测试，得到的 p 值确实证明了这一点，小雪的 p 值为 0.75，暴雨为 0.26，雪为 0.09。这只是意味着结果可能是由于随机变化或混杂变量(如风速/温度——稍后会有更多介绍)造成的。

似乎温度本身并没有太大的影响。虽然似乎有理由假设极热或极冷的温度会对四分卫的表现产生现实世界的影响，特别是极寒，但这些类别的大误差棒再次表明样本量不足以达到统计意义。只是为了好玩，我决定绘制一张在最冷的温度下打最多比赛的 NFL 球队的热图。

好了，我们已经看到了孤立的天气类别如何影响 QB 统计，但是在现实世界中，这些天气类别很少孤立出现。在不同的温度范围内，有风和降水的游戏会有什么影响？使用 plotly express，只需一行代码，我们就可以将所有这些变量可视化:

这里发生了很多事情:每个点或泡泡代表一个 QB 游戏。温度用颜色表示，得分的幻想点数对应于 y 轴上的高度，风速用 x 轴表示。此外，分散点的大小代表该游戏的最大阵风，这就是为什么在平均风速较高的游戏中，较大的气泡大多出现在每个图形的右侧。最后，每个支线剧情代表该游戏的沉淀类别。该图中的两个主要要点是幻想点数和风速/阵风之间的负相关性，跨所有降水类别，以及在有降水的游戏中得分的总体下降趋势和较低的总体幻想点数。

上面的图显示了相同的数据(减去温度)，但将所有的游戏分组在一起，这里气泡的颜色代表降水类别。

现在我们已经看到了一些风、温度和不同类型降水的总体影响的图形表示，我真的想深入分析一下在好天气和坏天气下玩的游戏中 QB 的表现。要成为好天气，它必须非常原始，没有任何天气因素会影响 QB 的表现。我选取了整体数据的一个子集，其中温度适中(在 60 到 80 度之间)，天空晴朗(或在圆顶中)，没有降水，平均风速低于 5 英里/小时。对于恶劣天气，我选择了有雨/风暴或雪(小雪和小雨除外)、平均风速超过 12 英里/小时或极端极地温度(低于 10 度或高于 95 度)的子集。由此产生的子集每个都包括 2000 多一点的 QB 游戏，这是一个不错的样本量。让我们看看恶劣天气对 QB 性能的影响有多大:

Negative % difference for all categories indicates worse QB stats in bad vs great weather, as expected

这里的主要要点是，在恶劣天气下，生产指标通常会有较大的百分比下降。效率指标虽然也普遍较差，但没有那么明显。这可能是因为在天气不好的比赛中，教练会减少传球次数。也许在未来，我会看一看运行数据，并在恶劣天气下寻找相应的量/冲尝试的增加。

最后，我想看看几个玩家，比较他们在好天气和坏天气游戏中的表现。我从经常在恶劣天气的主场比赛的球队中挑选了 3 名球员，他们是新英格兰爱国者队的汤姆·布拉迪、匹兹堡钢人队的本·罗特利斯伯格和绿湾包装工队的亚伦·罗杰斯。为了代表主场比赛在圆顶/好天气下进行的球队，我选择了印第安纳波利斯小马队的培顿·曼宁和圣路易斯公羊队/亚利桑那红雀队的库尔特·华纳。

不出所料，在《好天气 vs 坏天气》中，所有 5 个 QB 平均都有更多的幻想点数。恶劣的天气适应了 QBs 每个得分之间的伟大天气游戏多 2 至 3 个幻想点。同样不足为奇的是，库尔特·华纳的股价大幅下跌。虽然培顿·曼宁的微小差异乍一看有点令人惊讶，但曼宁以他的准备和一丝不苟的职业道德而闻名，这表现在他的表现一致性上，不管不利天气带来的任何挑战。

最后，作为一名铁杆新奥尔良圣徒队球迷，如果我不看一下德鲁·布里斯的统计数据，那就是我的失职。在梦幻足球世界中，人们相当普遍地声称，德鲁·布里斯在联赛中有一些更激烈的主场分裂，这意味着当他离开梅赛德斯·奔驰超级穹顶的友好“穹顶球场优势”时，他的表现会明显下降。这种定性准确吗？让我们看看这些数字是否支持这一声誉:

Say it ain’t so, Drew!!

唉，我不得不痛苦地承认，专家们是对的…

总之，这里是我从探索 NFL 天气数据中收集到的主要信息，因为它与 QB 的表现有关:

• 风/阵风的影响是影响 QB 性能的最重要的天气变量。从 0 英里/小时到大约 15 英里/小时的平均风速的影响相当轻微，但在这个阈值以上，幻想点的差异大约少 12%。在风速超过 20 英里/小时的游戏中，负面影响更加显著，每场游戏的幻想点数减少约 17%。
• 降水，就其本身而言，基本没有定论。尽管小雨游戏的平均幻想点数最低，但有雨/风暴和雪的游戏的数据在统计上并不显著。
• 温度作为一个孤立的天气变量，对 QB 的表现没有太大的影响。
• 总体来说，QBs 在天气好的时候比天气不好的时候多获得 8-10%的幻想分数。在相对百分比的基础上，在圆顶或天气好的城市打主场比赛的 QB 的一小部分精选样本，在恶劣天气下的表现并不比往往更经常在恶劣天气下比赛的 QB 差多少。

感谢阅读！为了更深入地探究这篇博文背后的数据和分析，这里是我的 github repo 的链接。

足球:为什么赢家赢了，输家输了

原文：https://towardsdatascience.com/football-why-winners-win-and-losers-lose-56665e5be90e?source=collection_archive---------36-----------------------

Photo by Vienna Reyes on Unsplash

探索欧洲足球 5 年

介绍

在这本笔记本中，我们将探索现代足球指标(xG、xGA 和 xPTS)及其对体育分析的影响。

• 预期目标(xG) —根据几个变量来衡量射门质量，如助攻类型、射门角度和与球门的距离，是否是头球射门以及是否被定义为大概率射门。
• 预期助攻(xGA) —衡量给定传球成为进球助攻的可能性。它考虑了几个因素，包括传球类型、传球终点和传球长度。
• 【xPTS】—衡量某场比赛给团队带来积分的可能性。

这些指标让我们更深入地观察足球统计数据，了解球员和球队的总体表现，并认识到运气和技巧在其中的作用。声明:它们都很重要。

这个 Kaggle 内核描述了这个笔记本的数据收集过程:网页抓取足球统计

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
import collections
import warnings

from IPython.core.display import display, HTML

*# import plotly* 
import plotly
import plotly.figure_factory as ff
import plotly.graph_objs as go
import plotly.offline as py
from plotly.offline import iplot, init_notebook_mode
import plotly.tools as tls

*# configure things*
warnings.filterwarnings('ignore')

pd.options.display.float_format = '**{:,.2f}**'.format  
pd.options.display.max_columns = 999

py.init_notebook_mode(connected=True)

%load_ext autoreload
%autoreload 2

%matplotlib inline
sns.set()

*# !pip install plotly --upgrade*

导入数据和可视化 EDA

df = pd.read_csv('../input/understat.com.csv')
df = df.rename(index=int, columns={'Unnamed: 0': 'league', 'Unnamed: 1': 'year'}) 
df.head()

在下一个可视化中，我们将检查在过去 5 年中每个联盟有多少支球队排名前 4。它可以给我们一些关于不同国家顶尖球队稳定性的信息。

f = plt.figure(figsize=(25,12))
ax = f.add_subplot(2,3,1)
plt.xticks(rotation=45)
sns.barplot(x='team', y='pts', hue='year', data=df[(df['league'] == 'Bundesliga') & (df['position'] <= 4)], ax=ax)
ax = f.add_subplot(2,3,2)
plt.xticks(rotation=45)
sns.barplot(x='team', y='pts', hue='year', data=df[(df['league'] == 'EPL') & (df['position'] <= 4)], ax=ax)
ax = f.add_subplot(2,3,3)
plt.xticks(rotation=45)
sns.barplot(x='team', y='pts', hue='year', data=df[(df['league'] == 'La_liga') & (df['position'] <= 4)], ax=ax)
ax = f.add_subplot(2,3,4)
plt.xticks(rotation=45)
sns.barplot(x='team', y='pts', hue='year', data=df[(df['league'] == 'Serie_A') & (df['position'] <= 4)], ax=ax)
ax = f.add_subplot(2,3,5)
plt.xticks(rotation=45)
sns.barplot(x='team', y='pts', hue='year', data=df[(df['league'] == 'Ligue_1') & (df['position'] <= 4)], ax=ax)
ax = f.add_subplot(2,3,6)
plt.xticks(rotation=45)
sns.barplot(x='team', y='pts', hue='year', data=df[(df['league'] == 'RFPL') & (df['position'] <= 4)], ax=ax)

正如我们从这些条形图中看到的，有些球队在过去 5 年中只进入过前 4 名，这意味着这并不常见，这意味着如果我们深入挖掘，我们会发现有一种运气因素可能对这些球队有利。这只是一个理论，所以让我们更接近那些离群值。

在过去的 5 个赛季中只有一次进入前 4 名的球队是:

• 来自德甲的沃尔夫斯堡(2014)和沙尔克 04 (2017)
• 来自 EPL 的莱斯特(2015)
• 来自西甲的比利亚雷亚尔(2015)和塞维利亚(2016)
• 来自意甲的拉齐奥(2014)和佛罗伦萨(2014)
• 来自法甲的里尔(2018)和圣艾蒂安(2018)
• 来自 RFPL 的罗斯托夫足球俱乐部(2015)和莫斯科迪纳摩足球俱乐部(2014)

让我们拯救这些队伍。

*# Removing unnecessary for our analysis columns* 
df_xg = df[['league', 'year', 'position', 'team', 'scored', 'xG', 'xG_diff', 'missed', 'xGA', 'xGA_diff', 'pts', 'xpts', 'xpts_diff']]

outlier_teams = ['Wolfsburg', 'Schalke 04', 'Leicester', 'Villareal', 'Sevilla', 'Lazio', 'Fiorentina', 'Lille', 'Saint-Etienne', 'FC Rostov', 'Dinamo Moscow']*# Checking if getting the first place requires fenomenal execution*
first_place = df_xg[df_xg['position'] == 1]

*# Get list of leagues*
leagues = df['league'].drop_duplicates()
leagues = leagues.tolist()

*# Get list of years*
years = df['year'].drop_duplicates()
years = years.tolist()

了解赢家如何获胜

在本节中，我们将尝试找到一些模式，这些模式可以帮助我们理解胜利之汤的一些成分:d .从德甲开始。

德甲联赛

first_place[first_place['league'] == 'Bundesliga']

pts = go.Bar(x = years, y = first_place['pts'][first_place['league'] == 'Bundesliga'], name = 'PTS')
xpts = go.Bar(x = years, y = first_place['xpts'][first_place['league'] == 'Bundesliga'], name = 'Expected PTS')

data = [pts, xpts]

layout = go.Layout(
    barmode='group',
    title="Comparing Actual and Expected Points for Winner Team in Bundesliga",
    xaxis={'title': 'Year'},
    yaxis={'title': "Points",
    }
)

fig = go.Figure(data=data, layout=layout)
py.iplot(fig)

通过看表格和柱状图我们看到，拜仁每年都获得了比他们应该得到的更多的分数，他们的得分比预期的多，失误比预期的少(除了 2018 年，它没有打破他们赢得赛季的计划，但它给出了一些暗示，拜仁今年打得更差，尽管竞争对手没有利用这一点)。

*# and from this table we see that Bayern dominates here totally, even when they do not play well*
df_xg[(df_xg['position'] <= 2) & (df_xg['league'] == 'Bundesliga')].sort_values(by=['year','xpts'], ascending=False)

西甲联赛

first_place[first_place['league'] == 'La_liga']

pts = go.Bar(x = years, y = first_place['pts'][first_place['league'] == 'La_liga'], name = 'PTS')
xpts = go.Bar(x = years, y = first_place['xpts'][first_place['league'] == 'La_liga'], name = 'Expected PTS')

data = [pts, xpts]

layout = go.Layout(
    barmode='group',
    title="Comparing Actual and Expected Points for Winner Team in La Liga",
    xaxis={'title': 'Year'},
    yaxis={'title': "Points",
    }
)

fig = go.Figure(data=data, layout=layout)
py.iplot(fig)

从上面的图表中我们可以看到，在 2014 年和 2015 年，巴塞罗那创造了足够多的赢得冠军的时刻，而不是依靠个人技能或运气，从这些数字中我们可以实际上说球队正在那里比赛。

2016 年，马德里和巴塞罗那之间有很多竞争，最终，马德里在一场特定的比赛中变得更幸运/更有勇气(或者巴塞罗那变得不幸运/没有球)，这是冠军的代价。我确信如果我们深入挖掘那个赛季，我们可以找到那场特殊的比赛。

2017 年和 2018 年，巴塞罗那的成功主要归功于莱昂内尔·梅西的行动，他在正常球员不会这样做的情况下得分或助攻。是什么导致 xPTS 差异如此之大。让我觉得(在这个赛季皇家马德里在转会市场上非常活跃的背景下)结果可能会很糟糕。只是基于数字和看巴萨比赛的主观看法。真希望我是错的。

*# comparing with runner-up*
df_xg[(df_xg['position'] <= 2) & (df_xg['league'] == 'La_liga')].sort_values(by=['year','xpts'], ascending=False)

EPL

first_place[first_place['league'] == 'EPL']

pts = go.Bar(x = years, y = first_place['pts'][first_place['league'] == 'EPL'], name = 'PTS')
xpts = go.Bar(x = years, y = first_place['xpts'][first_place['league'] == 'EPL'], name = 'Expected PTS')

data = [pts, xpts]

layout = go.Layout(
    barmode='group',
    title="Comparing Actual and Expected Points for Winner Team in EPL",
    xaxis={'title': 'Year'},
    yaxis={'title': "Points",
    }
)

fig = go.Figure(data=data, layout=layout)
py.iplot(fig)

在 EPL，我们看到一个明显的趋势，告诉你:“要想赢，你必须比统计数据更好”。有趣的案例是莱斯特在 2015 年的胜利故事:他们比他们应该得到的多了 12 分，同时阿森纳比预期少了 6 分！这就是我们热爱足球的原因，因为这样莫名其妙的事情发生了。我不是说这完全是运气，但它在这里发挥了作用。

另一件有趣的事是 2018 年的曼城——他们超级稳定！他们只比预期多进了一个球，少失了两个球，多得了 7 分，而利物浦打得很好，运气稍微好一点，但尽管比预期多了 13 分，还是没能赢。

Pep 正在完成制造毁灭机器。曼城根据技术创造和转换他们的时刻，而不是依靠运气——这使得他们在下一个赛季非常危险。

*# comparing with runner-ups* df_xg[(df_xg['position'] <= 2) & (df_xg['league'] == 'EPL')].sort_values(by=['year','xpts'], ascending=False)

法甲联赛

first_place[first_place['league'] == 'Ligue_1']

pts = go.Bar(x = years, y = first_place['pts'][first_place['league'] == 'Ligue_1'], name = 'PTS')
xpts = go.Bar(x = years, y = first_place['xpts'][first_place['league'] == 'Ligue_1'], name = 'Expected PTS')

data = [pts, xpts]

layout = go.Layout(
    barmode='group',
    title="Comparing Actual and Expected Points for Winner Team in Ligue 1",
    xaxis={'title': 'Year'},
    yaxis={'title': "Points",
    }
)

fig = go.Figure(data=data, layout=layout)
py.iplot(fig)

在法甲联赛中，我们继续看到这样的趋势:“要想赢，你必须付出 110%的努力，因为 100%是不够的”。在这里，巴黎圣日耳曼占据绝对优势。只有在 2016 年，我们在摩纳哥的比赛中遇到了一个异数，比预期多进了 30 个球！！！而且比预期多考了差不多 17 分！运气？相当不错的一部分。那一年 PSG 不错，但是摩纳哥非同一般。同样，我们不能说这是纯粹的运气或纯粹的技巧，而是两者在正确的时间和地点的完美结合。

*# comparing with runner-ups* df_xg[(df_xg['position'] <= 2) & (df_xg['league'] == 'Ligue_1')].sort_values(by=['year','xpts'], ascending=False)

意甲联赛

first_place[first_place['league'] == 'Serie_A']

pts = go.Bar(x = years, y = first_place['pts'][first_place['league'] == 'Serie_A'], name = 'PTS')
xpts = go.Bar(x = years, y = first_place['xpts'][first_place['league'] == 'Serie_A'], name = 'Expecetd PTS')

data = [pts, xpts]

layout = go.Layout(
    barmode='group',
    title="Comparing Actual and Expected Points for Winner Team in Serie A",
    xaxis={'title': 'Year'},
    yaxis={'title': "Points",
    }
)

fig = go.Figure(data=data, layout=layout)
py.iplot(fig)

在意大利甲级联赛中，尤文图斯已经连续 8 年称霸，尽管在冠军联赛中没有取得任何重大成功。我认为，通过检查这张图表和数字，我们可以了解到尤文在国内没有足够强的竞争，并获得了很多“幸运”的分数，这也源于多种因素，我们可以看到那不勒斯两次以 xPTS 超过尤文图斯，但这是现实生活，例如，在 2017 年，尤文疯狂地打进了额外的 26 个进球(或不知从哪里创造的进球)，而那不勒斯比预期多错过了 3 个进球(由于守门员的错误或可能是一些球队在 1 或 2 场特定比赛中的出色表现)。就像皇马成为冠军时西甲联赛的情况一样，我确信我们可以找到那一年的一两场关键比赛。

足球中细节很重要。你看，这里一个错误，那里一个木制品，你就失去了冠军。

*# comparing to runner-ups* df_xg[(df_xg['position'] <= 2) & (df_xg['league'] == 'Serie_A')].sort_values(by=['year','xpts'], ascending=False)

RFPL

first_place[first_place['league'] == 'RFPL']

pts = go.Bar(x = years, y = first_place['pts'][first_place['league'] == 'RFPL'], name = 'PTS')
xpts = go.Bar(x = years, y = first_place['xpts'][first_place['league'] == 'RFPL'], name = 'Expected PTS')

data = [pts, xpts]

layout = go.Layout(
    barmode='group',
    title="Comparing Actual and Expected Points for Winner Team in RFPL",
    xaxis={'title': 'Year'},
    yaxis={'title': "Points",
    }
)

fig = go.Figure(data=data, layout=layout)
py.iplot(fig)

我不关注俄罗斯超级联赛，所以仅仅通过冷眼看数据，我们就可以看到同样的模式，即从 2015 年到 2017 年，莫斯科中央陆军的得分超过了你应得的分数，这也是一个有趣的情况。在这些年里，这些人很优秀，但他们只转化了一次优势，其他两次——如果你不转化，你会受到惩罚，或者你的主要竞争对手只是转化得更好。

足球是没有公平可言的:d .虽然，我相信随着 VAR 的出现，接下来的几个赛季数字会变得更加稳定。因为那些额外的目标和分数的原因之一是仲裁人的错误。

*# comparing to runner-ups* df_xg[(df_xg['position'] <= 2) & (df_xg['league'] == 'RFPL')].sort_values(by=['year','xpts'], ascending=False)

统计概述

由于有 6 个联盟有不同的球队和统计数据，我决定集中在一个，在开始时，测试不同的方法，然后在其他 5 个上复制最终的分析模型。因为我主要看西甲联赛，所以我会从这场比赛开始，因为我对它了解最多。

*# Creating separate DataFrames per each league*
laliga = df_xg[df_xg['league'] == 'La_liga']
laliga.reset_index(inplace=True)
epl = df_xg[df_xg['league'] == 'EPL']
epl.reset_index(inplace=True)
bundesliga = df_xg[df_xg['league'] == 'Bundesliga']
bundesliga.reset_index(inplace=True)
seriea = df_xg[df_xg['league'] == 'Serie_A']
seriea.reset_index(inplace=True)
ligue1 = df_xg[df_xg['league'] == 'Ligue_1']
ligue1.reset_index(inplace=True)
rfpl = df_xg[df_xg['league'] == 'RFPL']
rfpl.reset_index(inplace=True)laliga.describe()

def print_records_antirecords(df):
  print('Presenting some records and antirecords: **n**')
  for col **in** df.describe().columns:
    if col **not** **in** ['index', 'year', 'position']:
      team_min = df['team'].loc[df[col] == df.describe().loc['min',col]].values[0]
      year_min = df['year'].loc[df[col] == df.describe().loc['min',col]].values[0]
      team_max = df['team'].loc[df[col] == df.describe().loc['max',col]].values[0]
      year_max = df['year'].loc[df[col] == df.describe().loc['max',col]].values[0]
      val_min = df.describe().loc['min',col]
      val_max = df.describe().loc['max',col]
      print('The lowest value of **{0}** had **{1}** in **{2}** and it is equal to **{3:.2f}**'.format(col.upper(), team_min, year_min, val_min))
      print('The highest value of **{0}** had **{1}** in **{2}** and it is equal to **{3:.2f}**'.format(col.upper(), team_max, year_max, val_max))
      print('='*100)

*# replace laliga with any league you want*
print_records_antirecords(laliga)***Presenting some records and antirecords:*** 

***The lowest value of SCORED had Cordoba in 2014 and it is equal to 22.00
The highest value of SCORED had Real Madrid in 2014 and it is equal to 118.00
================================================================
The lowest value of XG had Eibar in 2014 and it is equal to 29.56
The highest value of XG had Barcelona in 2015 and it is equal to 113.60
================================================================
The lowest value of XG_DIFF had Barcelona in 2016 and it is equal to -22.45
The highest value of XG_DIFF had Las Palmas in 2017 and it is equal to 13.88
================================================================
The lowest value of MISSED had Atletico Madrid in 2015 and it is equal to 18.00
The highest value of MISSED had Osasuna in 2016 and it is equal to 94.00
================================================================
The lowest value of XGA had Atletico Madrid in 2015 and it is equal to 27.80
The highest value of XGA had Levante in 2018 and it is equal to 78.86
================================================================
The lowest value of XGA_DIFF had Osasuna in 2016 and it is equal to -29.18
The highest value of XGA_DIFF had Valencia in 2015 and it is equal to 13.69
================================================================
The lowest value of PTS had Cordoba in 2014 and it is equal to 20.00
The highest value of PTS had Barcelona in 2014 and it is equal to 94.00
================================================================
The lowest value of XPTS had Granada in 2016 and it is equal to 26.50
The highest value of XPTS had Barcelona in 2015 and it is equal to 94.38
================================================================
The lowest value of XPTS_DIFF had Atletico Madrid in 2017 and it is equal to -17.40
The highest value of XPTS_DIFF had Deportivo La Coruna in 2017 and it is equal to 20.16***

现在让我们在图表上看看这些数字。

trace0 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2014], 
    y = laliga['xG_diff'][laliga['year'] == 2014],
    name = '2014',
    mode = 'lines+markers'
)

trace1 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2015], 
    y = laliga['xG_diff'][laliga['year'] == 2015],
    name='2015',
    mode = 'lines+markers'
)

trace2 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2016], 
    y = laliga['xG_diff'][laliga['year'] == 2016],
    name='2016',
    mode = 'lines+markers'
)

trace3 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2017], 
    y = laliga['xG_diff'][laliga['year'] == 2017],
    name='2017',
    mode = 'lines+markers'
)

trace4 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2018], 
    y = laliga['xG_diff'][laliga['year'] == 2018],
    name='2018',
    mode = 'lines+markers'
)

data = [trace0, trace1, trace2, trace3, trace4]

layout = go.Layout(
    title="Comparing xG gap between positions",
    xaxis={'title': 'Year'},
    yaxis={'title': "xG difference",
    }
)

fig = go.Figure(data=data, layout=layout)
py.iplot(fig)

trace0 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2014], 
    y = laliga['xGA_diff'][laliga['year'] == 2014],
    name = '2014',
    mode = 'lines+markers'
)

trace1 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2015], 
    y = laliga['xGA_diff'][laliga['year'] == 2015],
    name='2015',
    mode = 'lines+markers'
)

trace2 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2016], 
    y = laliga['xGA_diff'][laliga['year'] == 2016],
    name='2016',
    mode = 'lines+markers'
)

trace3 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2017], 
    y = laliga['xGA_diff'][laliga['year'] == 2017],
    name='2017',
    mode = 'lines+markers'
)

trace4 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2018], 
    y = laliga['xGA_diff'][laliga['year'] == 2018],
    name='2018',
    mode = 'lines+markers'
)

data = [trace0, trace1, trace2, trace3, trace4]

layout = go.Layout(
    title="Comparing xGA gap between positions",
    xaxis={'title': 'Year'},
    yaxis={'title': "xGA difference",
    }
)

fig = go.Figure(data=data, layout=layout)
py.iplot(fig)

trace0 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2014], 
    y = laliga['xpts_diff'][laliga['year'] == 2014],
    name = '2014',
    mode = 'lines+markers'
)

trace1 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2015], 
    y = laliga['xpts_diff'][laliga['year'] == 2015],
    name='2015',
    mode = 'lines+markers'
)

trace2 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2016], 
    y = laliga['xpts_diff'][laliga['year'] == 2016],
    name='2016',
    mode = 'lines+markers'
)

trace3 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2017], 
    y = laliga['xpts_diff'][laliga['year'] == 2017],
    name='2017',
    mode = 'lines+markers'
)

trace4 = go.Scatter(
    x = laliga['position'][laliga['year'] == 2018], 
    y = laliga['xpts_diff'][laliga['year'] == 2018],
    name='2018',
    mode = 'lines+markers'
)

data = [trace0, trace1, trace2, trace3, trace4]

layout = go.Layout(
    title="Comparing xPTS gap between positions",
    xaxis={'title': 'Position'},
    yaxis={'title': "xPTS difference",
    }
)

fig = go.Figure(data=data, layout=layout)
py.iplot(fig)

从上面的图表中，我们可以清楚地看到，顶级球队得分更多，失球更少，获得的分数比预期的更多。这就是为什么这些球队都是顶级球队。和与外人完全相反的情况。中游球队的平均水平。完全符合逻辑，没有什么深刻的见解。

*# Check mean differences*
def get_diff_means(df):  
  dm = df.groupby('year')[['xG_diff', 'xGA_diff', 'xpts_diff']].mean()

  return dm

means = get_diff_means(laliga)
means

*# Check median differences*
def get_diff_medians(df):  
  dm = df.groupby('year')[['xG_diff', 'xGA_diff', 'xpts_diff']].median()

  return dm

medians = get_diff_medians(laliga)
medians

离群点检测

z 分数

z 得分是数据点相对于平均值的标准偏差数。假设|z-score| > 3 是异常值，我们可以使用它来查找数据集中的异常值。

*# Getting outliers for xG using zscore*
from scipy.stats import zscore
*# laliga[(np.abs(zscore(laliga[['xG_diff']])) > 2.0).all(axis=1)]*
df_xg[(np.abs(zscore(df_xg[['xG_diff']])) > 3.0).all(axis=1)]

*# outliers for xGA*
*# laliga[(np.abs(zscore(laliga[['xGA_diff']])) > 2.0).all(axis=1)]*
df_xg[(np.abs(zscore(df_xg[['xGA_diff']])) > 3.0).all(axis=1)]

*# Outliers for xPTS*
*# laliga[(np.abs(zscore(laliga[['xpts_diff']])) > 2.0).all(axis=1)]*
df_xg[(np.abs(zscore(df_xg[['xpts_diff']])) > 3.0).all(axis=1)]

zscore 总共检测到 12 个异常值。2016 年可怜的奥萨苏纳——近 30 个不该进的球。

正如我们从这个数据中看到的，处于异常空间顶端并不能让你赢得这个赛季。但是如果你错过了你的机会，或者在不该进球的地方进了球，而且进得太多——你应该降级。失败和平庸比胜利要容易得多。

四分位数间距(IQR)

IQR —是一组数据的第一个四分位数和第三个四分位数之间的差值。这是描述一组数据分布的一种方式。

一个常用的规则是，如果一个数据点高于第三个四分位数或低于第一个四分位数超过 1.5 个⋅ IQR，则该数据点是异常值。换句话说，低异常值低于 Q1 1.5⋅iqr，高异常值高于 Q3 + 1.5 ⋅ IQR。

我们去看看。

*# Trying different method of outliers detection* 
df_xg.describe()

*# using Interquartile Range Method to identify outliers*
*# xG_diff*
iqr_xG = (df_xg.describe().loc['75%','xG_diff'] - df_xg.describe().loc['25%','xG_diff']) * 1.5
upper_xG = df_xg.describe().loc['75%','xG_diff'] + iqr_xG
lower_xG = df_xg.describe().loc['25%','xG_diff'] - iqr_xG

print('IQR for xG_diff: **{:.2f}**'.format(iqr_xG))
print('Upper border for xG_diff: **{:.2f}**'.format(upper_xG))
print('Lower border for xG_diff: **{:.2f}**'.format(lower_xG))

outliers_xG = df_xg[(df_xg['xG_diff'] > upper_xG) | (df_xg['xG_diff'] < lower_xG)]
print('='*50)

*# xGA_diff*
iqr_xGA = (df_xg.describe().loc['75%','xGA_diff'] - df_xg.describe().loc['25%','xGA_diff']) * 1.5
upper_xGA = df_xg.describe().loc['75%','xGA_diff'] + iqr_xGA
lower_xGA = df_xg.describe().loc['25%','xGA_diff'] - iqr_xGA

print('IQR for xGA_diff: **{:.2f}**'.format(iqr_xGA))
print('Upper border for xGA_diff: **{:.2f}**'.format(upper_xGA))
print('Lower border for xGA_diff: **{:.2f}**'.format(lower_xGA))

outliers_xGA = df_xg[(df_xg['xGA_diff'] > upper_xGA) | (df_xg['xGA_diff'] < lower_xGA)]
print('='*50)

*# xpts_diff*
iqr_xpts = (df_xg.describe().loc['75%','xpts_diff'] - df_xg.describe().loc['25%','xpts_diff']) * 1.5
upper_xpts = df_xg.describe().loc['75%','xpts_diff'] + iqr_xpts
lower_xpts = df_xg.describe().loc['25%','xpts_diff'] - iqr_xpts

print('IQR for xPTS_diff: **{:.2f}**'.format(iqr_xpts))
print('Upper border for xPTS_diff: **{:.2f}**'.format(upper_xpts))
print('Lower border for xPTS_diff: **{:.2f}**'.format(lower_xpts))

outliers_xpts = df_xg[(df_xg['xpts_diff'] > upper_xpts) | (df_xg['xpts_diff'] < lower_xpts)]
print('='*50)

outliers_full = pd.concat([outliers_xG, outliers_xGA, outliers_xpts])
outliers_full = outliers_full.drop_duplicates()***IQR for xG_diff: 13.16
Upper border for xG_diff: 16.65
Lower border for xG_diff: -18.43
==================================================
IQR for xGA_diff: 13.95
Upper border for xGA_diff: 17.15
Lower border for xGA_diff: -20.05
==================================================
IQR for xPTS_diff: 13.93
Upper border for xPTS_diff: 18.73
Lower border for xPTS_diff: -18.41
==================================================****# Adding ratings bottom to up to find looser in each league (different amount of teams in every league so I can't do just n-20)*
max_position = df_xg.groupby('league')['position'].max()
df_xg['position_reverse'] = np.nan
outliers_full['position_reverse'] = np.nan

for i, row **in** df_xg.iterrows():
  df_xg.at[i, 'position_reverse'] = np.abs(row['position'] - max_position[row['league']])+1

for i, row **in** outliers_full.iterrows():
  outliers_full.at[i, 'position_reverse'] = np.abs(row['position'] - max_position[row['league']])+1total_count = df_xg[(df_xg['position'] <= 4) | (df_xg['position_reverse'] <= 3)].count()[0]
outlier_count = outliers_full[(outliers_full['position'] <= 4) | (outliers_full['position_reverse'] <= 3)].count()[0]
outlier_prob = outlier_count / total_count
print('Probability of outlier in top or bottom of the final table: **{:.2%}**'.format(outlier_prob))***Probability of outlier in top or bottom of the final table: 8.10%***

因此，我们可以说，很有可能每年在 6 个联赛中的一个联赛中，会有一支球队凭借他们出色的技术和运气获得冠军联赛或欧罗巴联赛的门票，或者有一支球队进入了乙级，因为他们不能转换他们的时刻。

*# 1-3 outliers among all leagues in a year*
data = pd.DataFrame(outliers_full.groupby('league')['year'].count()).reset_index()
data = data.rename(index=int, columns={'year': 'outliers'})
sns.barplot(x='league', y='outliers', data=data)
*# no outliers in Bundesliga*

我们表现出色和表现不佳的赢家和输家。

top_bottom = outliers_full[(outliers_full['position'] <= 4) | (outliers_full['position_reverse'] <= 3)].sort_values(by='league')
top_bottom

*# Let's get back to our list of teams that suddenly got into top. Was that because of unbeliavable mix of luck and skill?*
ot = [x for x  **in** outlier_teams if x **in** top_bottom['team'].drop_duplicates().tolist()]
ot
*# The answer is absolutely no. They just played well during 1 season. Sometimes that happen.*[]

结论

足球是一项低分的比赛，一个进球可以改变比赛的整个画面，甚至结束比赛结果。这就是为什么长期分析能让你对形势有更好的了解。

随着 xG 指标(以及由此衍生的其他指标)的引入，我们现在可以真正评估团队的长期表现，并了解顶级团队、中产阶级团队和绝对局外人之间的差异。

xG 为围绕足球的讨论带来了新的论点，这使得足球变得更加有趣。而同时，游戏又不失这种不确定的因素，以及疯狂的事情发生的可能性。其实现在，这些疯狂的事情有了解释的机会。

最终，我们发现，在某个联赛中，发生怪事的几率几乎是 100%。这将是一个多么史诗般的时刻，只是时间的问题。

带互动图表的原创作品可以在这里找到。

原载于 2019 年 9 月 11 日【http://sergilehkyi.com】。

为了成功采用人工智能，我们需要更多的女性领导者

原文：https://towardsdatascience.com/for-the-successful-future-of-ai-women-have-to-take-the-lead-180f09be4e50?source=collection_archive---------14-----------------------

Why we need more women in AI

为什么女性拥有正确的品质来更成功地领导人工智能项目，并创造一个协作和包容的环境来构建现实世界的人工智能产品。

“21 世纪最激动人心的突破不会因为技术而出现，而是因为对人类意义的不断扩展的概念”

约翰·奈斯比特

第一部分:男女之间的差异

一个关于女性和男性如何以不同方式工作的故事

在我们深入探讨为什么更多女性应该领导人工智能团队之前，我想分享一个我从卢塞恩应用科学与艺术大学三年级学生 Tania Biland 那里听到的有趣故事。

塔尼亚讲述的故事:

上学期，为了给瑞士或德国品牌开发安全技术解决方案，我们班分成了三个不同的小组:

第一组: 仅限女性(我的组)

第二组: 仅限男性

第三组: 四女一男

经过 4 周的工作后，每个团队都要展示他们的作品。

第 1 组 ，仅由女性组成，在黑暗中为女性制定了安全解决方案。由于评审团只有男性，我们决定用角色、音乐和视频来讲述一个故事，以便让他们感受女性在日常生活中所经历的事情。我们还强调了这样一个事实，即每个人的生活中都有母亲、姐妹或妻子，他们可能不希望她/他们受苦。最终，我们的解决方案是相当简单的技术上的 :利用灯光提供安全但在情感上与观众 相连。

组 2 ，大多由男性组成，呈现了一个更具 高科技解决方案使用 AIGPS视频会议 。他们以 事实和数字 为论据，指出自己的 竞争优势 。**

*****第 3 组，*4 女 1 男，结局好像还没完。小组中唯一的男性不同意被女性领导，因此他们花太多时间讨论小组动态，而不是工作。

这些小组不仅有不同的结果，而且处理问题的方式也不同。我所在的小组(小组 1)决定从定义彼此的工作偏好和风格开始，以便分配一些责任，并尽可能保持层级的扁平化。

另一方面，另外两个小组选出了一名组长。事实证明，这些“领导者”更像是独裁者，这导致了严重的冲突，团队花了几个小时讨论和争论，而我们的团队只是在工作和生产。

科学告诉我们性别差异

关于性别差异和对行为的影响的科学图景仍在发展，还没有为不同的行为提出一套清晰的科学解释。

通过汇总大部分研究，影响行为的主要因素有两个:

1. 男女之间潜在的生理差异
2. ****社会规范和压力形成不同的行为

*在上面的故事中，正如 Tania 所讲述的，女性以一种协作的领导风格(【adhocracy culture】*)开发了解决方案，以一种近乎扁平的层级结构根据任务调整领导地位。他们通过让所有利益相关者(在这种情况下，母亲和妻子=用户)参与进来，对他们的问题表示同情，从而得出他们的论点。他们看到了更广阔的前景，也构建了一个实际上已经完成的更简单的解决方案。

通过这个故事，我能够将为什么大多数人工智能项目最终都没有从原型阶段发展到现实应用的原因联系起来。

第二部分:让人工智能成功

为什么不采用 AI 产品？

根据我的经验，大多数人工智能和机器学习(ML)解决方案没有从原型阶段走向现实世界有三个主要原因:

1. ****信任缺失:AI 或者 ML 产品最大的困难之一就是信任缺失。数百万美元已经花在原型上，但在现实世界的发射中几乎没有成功。本质上，做生意和为客户提供价值的最基本的价值观之一是信任，当涉及到道德问题和相关的信任问题时，人工智能是争论最激烈的技术。信任来自于在整个开发阶段涉及不同的选项和各方，这在原型阶段是做不到的。
2. ****发布的复杂性:构建一个原型很容易，但在进入现实世界时，需要考虑数十个其他外部实体。除了技术挑战之外，还有其他需要与原型集成的领域(比如营销、设计和销售)。
3. 人工智能产品通常不会考虑到所有的利益相关者:我听说过一个故事，在家庭暴力的情况下，Alexa 和 Google Home 被男性用来锁定他们的配偶。他们把音乐开得很大声，或者把他们锁在门外。在一个大多数男性工程师构建这些产品的环境中，可能没有人会考虑这种情况。此外，有许多例子表明人工智能和数据传感器可能带有偏见、性别歧视和种族主义[1]。

有趣的是，这三点都与技术挑战无关，所有这些都可以通过创建正确的团队来克服。

如何让 AI 更成功的被采用？

为了解决上述挑战并构建更成功的人工智能产品，我们需要专注于一种更加协作和社区驱动的方法**。**

这考虑了不同利益攸关方的意见，特别是那些代表性不足的利益攸关方的意见。以下是实现这一目标的步骤:

第一步。让不同的团体参与进来。来自人才金字塔中间的女性

在技术领域，大多数公司都专注于雇佣人才金字塔顶端的人，主要由于历史原因，这里的女性较少。例如，大多数计算机科学班的女生不到 10%。然而，许多有才华的女性隐藏在金字塔的中间，通过在线课程自学却缺乏机会和鼓励。

Talent Pyramid

举个例子，我和极客女孩胡萝卜的主席聊过，这是一个在技术领域促进女性的组织。他们正在组织一个人工智能研讨会，有超过 125 名女性申请，但她们只有 25 个席位**，所以自然地，她们不得不留下100 多名有才华的女性。**

想象一下，如果我们能让其他 100 名女性中的大多数参与进来，而不仅仅是高层。这将给更多女性提供在人工智能等新技术领域工作的机会。

第二步。与不同的利益相关者建立一个自下而上的公共协作团队

接下来，我们需要男性和女性以及不同利益相关者之间的更多合作，以便在真实市场中成功推出产品。这可以通过形成包容性的项目社区来实现，这些社区基于共同的价值观、信仰和更大的愿景来构建人工智能产品。

为了证明这一点，在过去的六个月里，我们召集了 50 多名男女学生来构建一个 ML 模型。在很短的时间内，成员们开始合作并互相帮助来建立模型。形成了四个小组，其中一个由两个女人领导，两个男人支持(数据标签)。其他组都是男性。在 4 个月的时间里，由两名女性和两名男性组成的小组建立了最精确的模型。从一开始，女性就比男性更愿意合作。然而，更有趣的是，我发现，由于团队中的其他女性，团队中的男性最终也表现得更加合作。这太迷人了！！

第三步。为协作创建正确的组织结构

如果我们可以创造不需要授权的组织结构和实践，因为从设计上来说，每个人都是强大的，没有人是无能为力的，那会怎么样？我已经看到，这可以通过连接内在和外在动机(与金钱无关)并创建一个没有竞争力的激励结构来实现。

在我的例子中，我建立了一个社区，其中导师位于金字塔的顶端，接着是社区经理，然后是致力于建立模型的工程师，最后是数据标记者。每个团队的成员都在努力向上爬，以达到下一个级别，这创造了一个外在的动机。然而，同一级别人员的货币补偿是相同的。这促进了合作。

在这种情况下，领导者的角色不是做老板，而是培养合作领导能力。这样的组织结构将减少控制人的需要，并提供共同学习和成长的机会[2]。

第三部分:连接第一部分和第二部分。

为什么女性应该领导人工智能团队

在一开始的故事中，女性小组遵循了一种更加协作的领导风格**，表现出更多的顾客同理心和协作意愿。**

考虑到 solar 项目中的有限实验，我们发现使用社区来构建产品的方法也有助于促进社区成员之间的协作和信任。

虽然上面提到的品质都不能一概而论，但下图旨在总结许多女性非常适合协作型领导的一些原因。

总之，我认为:

我们应该更全面地思考，尽最大努力创造合适的环境，让我们超越性别、种族和文化背景，专注于我们人类如何合作建设更美好的未来。

最后，我要感谢所有帮助我撰写这篇文章的人(包括男性和女性)——他们创造了体验和内容，也提炼了文本。

如果你想获得我们的人工智能挑战的更新，获得专家采访和提高你的人工智能技能的实用技巧，订阅我们的每月简讯。

我们也在脸书、 LinkedIn 和 Twitter 上。

预测 KPI:RMSE、梅伊、MAPE 和 Bias

原文：https://towardsdatascience.com/forecast-kpi-rmse-mae-mape-bias-cdc5703d242d?source=collection_archive---------0-----------------------

选择正确的预测指标并不简单。让我们回顾一下 RMSE、梅、MAPE 和偏见的利弊。剧透:MAPE 是最差的。不要用。

本文摘自我的著作 供应链预测的数据科学 。你可以在这里 阅读我的其他文章 。我也活跃在LinkedIn。

Credit

衡量预测的准确性(或误差)并不是一件容易的事情，因为没有放之四海而皆准的指标。只有实验才能告诉你什么样的关键绩效指标(KPI)最适合你。正如你将看到的，每个指标将避免一些陷阱，但会倾向于其他。

我们要做的第一个区分是预测的精确度和它的偏差之间的区别:

• 偏差代表历史平均误差。基本上，平均而言，你的预测是过高(即你超过需求)还是过低(即你低于需求)？这将给出错误的总体方向。
• 精度衡量预测值和实际值之间的差距。预测的精度给出了误差大小的概念，但不是它们的总体方向。

当然，正如你在下图中看到的，我们想要的是一个既精确又公正的预测。

预测 KPI

错误

让我们首先将误差定义为预测减去需求。

请注意，如果根据此定义预测超出需求，误差将为正。如果预测低于需求，那么误差将为负。

偏见

偏差定义为平均误差:

其中 n 是既有预测又有需求的历史周期数。

由于一个项目上的正误差可以抵消另一个项目上的负误差，预测模型可以实现非常低的偏差，同时又不精确。显然，单凭偏差不足以评估你的预测精度。但是高度偏差的预测已经表明模型中有问题。

multidimensional assessment of philosophy of education 教育哲学的多维评价

平均绝对百分比误差(MAPE) 是衡量预测准确度最常用的 KPI 之一。

MAPE 是个体绝对误差的总和除以需求(每个时期分开)。它是百分比误差的平均值。

MAPE 是一个非常奇怪的预测 KPI。

尽管这是一个不太准确的指标，但它在企业经理中是众所周知的。正如你在公式中看到的，MAPE 将每个误差单独除以需求，所以它是偏斜的:低需求时期的高误差将显著影响 MAPE。因此，优化 MAPE 将导致一个奇怪的预测，很可能会低于需求。避开就好。

平均绝对误差

平均绝对误差(MAE) 是衡量预测准确性的一个非常好的 KPI。顾名思义，就是绝对误差的均值。

此 KPI 的首要问题之一是它没有根据平均需求进行调整。如果有人告诉你，某个特定项目的 MAE 是 10，你无法知道这是好是坏。如果你的平均需求是 1000，这当然是惊人的。尽管如此，如果平均需求是 1，这是一个非常差的精度。为了解决这个问题，通常将 MAE 除以平均需求得到一个百分比:

MAPE/梅混淆——似乎许多从业者使用梅公式，并将其称为 MAPE。这可能会造成很多混乱。当与某人讨论预测误差时，我总是建议你明确地展示你如何计算预测误差，以确保比较苹果和苹果。

均方根误差

均方根误差(RMSE) 是一个奇怪的 KPI，但是非常有用，我们将在后面讨论。它被定义为平均平方误差的平方根。

就像 MAE 一样，RMSE 也不能满足需求。我们可以这样定义 RMSE%,

实际上，许多算法(尤其是机器学习)都是基于均方误差(MSE) ，这与 RMSE 直接相关。

许多算法使用 MSE，因为它比 RMSE 计算更快，更容易操作。但是它没有被缩放到原始误差(因为误差是平方的)，导致我们无法将 KPI 与原始需求缩放相关联。因此，我们不会用它来评估我们的统计预测模型。

误差加权的问题

与梅相比，RMSE 对每个错误都一视同仁。它更重视最重要的错误。这意味着一个大的错误就足以得到一个非常糟糕的 RMSE。

让我们用一个虚拟需求时间序列来做一个例子。

假设我们想要比较两个略有不同的预测。唯一的区别是对最新需求观察的预测:预测#1 比它少 7 个单位，而预测#2 只比少 6 个单位。

如果我们查看这两个预测的 KPI，我们会发现:

有趣的是，通过将最后一个周期的误差改变一个单位，我们将总 RMSE 降低了 6.9% (2.86 到 2.66)。尽管如此，MAE 只减少了 3.6% (2.33 到 2.25)，所以对 MAE 的影响低了近两倍。显然，RMSE 强调最重要的错误，而梅对每个错误都给予同等的重视。你可以自己尝试一下，减少一个最准确时期的误差，观察对梅和 RMSE 的影响。

剧透:对 RMSE 几乎没有影响。

稍后你会看到，RMSE 还有其他一些非常有趣的特性。

你想预测什么？

我们仔细检查了这些 KPI 的定义(bias、MAPE、MAE、RMSE)，但是仍然不清楚使用一个 KPI 代替另一个 KPI 会对我们的模型产生什么影响。有人可能会认为，用 RMSE 代替梅或梅代替 MAPE 不会改变任何事情。但是没有什么是不真实的。

让我们做一个简单的例子来说明这一点。

让我们想象一下，一个产品的周需求量很低，而且相当平稳，但不时会有一个大订单(可能是由于促销或客户批量订购)。这是我们目前观察到的每周需求:

现在让我们设想我们对这种产品提出三种不同的预测。第一个预测 2 件/天，第二个预测 4 件，最后一个预测 6 件。让我们把我们观察到的需求和这些预测画出来。

让我们看看这些预测在历史时期的偏差、MAPE、梅伊和 RMSE 方面的表现:

这意味着，就 MAPE 而言，预测#1 是历史期间最好的，就平均月平均日而言，预测#2 是最好的。预测#3 在 RMSE 和偏差方面是最好的(但在梅和 MAPE 方面是最差的)。现在让我们来揭示这些预测是如何做出的:

• 预测 1 只是一个很低的量。
• 预测 2 是需求中位数:4。
• 预测 3 是平均需求。

中位数与平均数—数学优化

在进一步讨论不同的预测 KPI 之前，让我们花一些时间来理解为什么中位数的预测会得到好的 MAE 和平均数的预测会得到好的 RMSE 。

前面有一点数学。如果你对这些等式不清楚，这不是问题，不要气馁。跳过它们，直接跳到 RMSE 和梅埃塔勒的结论。

RMSE

让我们从 RMSE 开始:

为了简化下面的代数，让我们使用一个简化版本:均方误差(MSE):

如果您将 MSE 设置为预测模型的目标，它会将其最小化。人们可以通过将一个数学函数的导数设为零来最小化它。让我们试试这个。

结论为了优化预测的 MSE，该模型必须以总预测等于总需求为目标。也就是说，优化 MSE 的目的是产生平均正确的预测，因此是无偏的。

梅

现在让我们为梅做同样的事情。

或者，

和

这意味着

结论为了优化 MAE(即将其导数设置为 0)，预测需要比需求高很多倍，因为它比需求低。换句话说，我们正在寻找一个将数据集分成两个相等部分的值。这是中位数的确切定义。

MAPE

不幸的是，MAPE 的衍生物不会显示一些优雅和简单的属性。我们可以简单地说，MAPE 提出了一个非常低的预测，因为当需求较低时，它为预测误差分配了较高的权重。

结论

正如我们上面看到的，在任何模型中，RMSE 的优化将寻求平均正确。相比之下，MAE 的优化将试图经常超过需求和低于需求，这意味着以需求中值为目标。我们必须明白，一个显著的区别在于梅和 RMSE 的数学根源。一个瞄准中位数，一个瞄准平均数。

梅和 RMSE——选择哪一个？

以需求的中间值或平均值为目标更糟糕吗？嗯，答案不是非黑即白的。每种技术都有一些好处和一些风险，我们将在下一页讨论。只有实验才能揭示哪种方法最适合数据集。你甚至可以选择同时使用 RMSE 和梅。

让我们花些时间来讨论选择 RMSE 或 MAE 对偏差、对异常值的敏感性和间歇性需求的影响。

偏见

对于许多产品，你会发现中间值与平均需求并不相同。需求很可能会在这里或那里出现一些峰值，这将导致一个偏斜的分布。这些不对称的需求分布在供应链中普遍存在，因为峰值可能是由于定期促销或客户批量订购造成的。这将导致需求中位数低于平均需求，如下所示。

这意味着最小化 MAE 的预测将导致偏差。相比之下，最小化 RMSE 的预测不会导致偏差(因为它的目标是平均值)。这绝对是 MAE 的主要弱点。

对异常值的敏感度

正如我们所讨论的，RMSE 更重视最高的错误。这是有代价的:对异常值的敏感性。让我们想象一个需求模式如下的物品。

中位数是 8.5，平均数是 9.5。我们已经观察到，如果我们做出最小化 MAE 的预测，我们将预测中位数(8.5)，并且我们将平均低于需求 1 个单位(偏差= -1)。然后，您可能希望最小化 RMSE，并预测平均值(9.5)来避免这种情况。

然而，现在让我们假设我们有一个新的需求观察值 100。

中位数还是 8.5(没变！)，但是现在平均是 18.1。在这种情况下，您可能不想预测平均值，而可能会恢复到预测中值。

一般来说，中位数比平均数对异常值更稳健。这在供应链环境中相当重要，因为由于编码错误或需求高峰(营销、促销、现货交易)，我们可能会面临许多异常值。

对异常值的鲁棒性总是一件好事吗？号码

间歇需求

事实上，不幸的是，中位数对异常值的稳健性会对具有间歇性需求的项目产生非常恼人的影响。

让我们想象一下，我们把一个产品卖给一个客户。这是一种高利润的产品。不幸的是，我们独一无二的客户似乎每三周就下一次订单，没有任何规律可循。客户总是以 100 件为一批订购产品。然后，我们的平均周需求为 33 件，需求中位数为… 0。

我们必须填充该产品的每周预测。假设我们做了第一次预测，目标是平均需求(33 件)。从长期来看，我们将获得 6，667 的总平方误差(RMSE 为 47)和 133 的总绝对误差(平均误差为 44)。

现在，如果我们预测需求中位数(0)，我们会得到 100 的总绝对误差(MAE 为 33)和 10.000 的总平方误差(RMSE 为 58)。

正如我们所看到的，MAE 对于间歇性需求来说是一个不好的 KPI。一旦有超过一半的期间没有需求，最佳预测就是… 0！

结论

梅保护离群值，而 RMSE 向我们保证得到一个无偏见的预测。您应该使用哪个指标？不幸的是，没有明确的答案。作为一名供应链数据科学家，您应该进行试验:如果使用 MAE 作为 KPI 会导致高偏差，您可能会想要使用 RMSE。如果数据集包含许多异常值，导致预测不准确，您可能需要使用 MAE。

请注意，您可以选择用一个或多个 KPI(通常是 MAE 和 bias)报告预测误差，并使用另一个 KPI(RMSE？)来优化你的模型。

针对低需求项目使用的最后一个技巧是将需求汇总到更高的时间范围。例如，如果每周的需求很低，您可以测试月度预测甚至季度预测。您总是可以通过简单的划分将预测分解回原始时段。这种技术可以让您将 MAE 用作 KPI，同时平滑需求峰值。

关于作者

在 LinkedIn 上关注我！

[## Nicolas vande put——顾问，创始人——供应链| LinkedIn

查看 Nicolas Vandeput 在世界上最大的职业社区 LinkedIn 上的个人资料。尼古拉斯有 7 份工作列在…

www.linkedin.com](https://www.linkedin.com/in/vandeputnicolas/)

icolas Vandeput 是一名供应链数据科学家，擅长需求预测和库存优化。他在 2016 年创立了自己的咨询公司 SupChains ，并在 2018 年共同创立了 SKU 科学——一个快速、简单、实惠的需求预测平台。尼古拉斯对教育充满热情，他既是一个狂热的学习者，也喜欢在大学教学:自 2014 年以来，他一直在比利时布鲁塞尔为硕士学生教授预测和库存优化。自 2020 年以来，他还在法国巴黎的 CentraleSupelec 教授这两个科目。他于 2018 年出版了 供应链预测的数据科学(2021 年第 2 版)和 2020 年出版了 库存优化:模型与模拟 。

先知中使用附加回归器的预测模型调整

原文：https://towardsdatascience.com/forecast-model-tuning-with-additional-regressors-in-prophet-ffcbf1777dda?source=collection_archive---------3-----------------------

关于如何在先知模型中添加回归因子以提高预测精度的说明

Use case — bike rental (picture source: Pixabay)

我将把我的实验结果分享给预言家附加回归器。我的目标是检查额外的回归将如何影响先知计算的预测。

使用来自 Kaggle 的数据集— 华盛顿州的自行车共享数据集。数据附带了每天自行车租金和天气情况的数字。我创建并比较了三个模型:

1. 带有自行车租赁日期和数量的时间序列先知模型
2. 带有附加回归器的模型——天气温度
3. 带有附加回归器 s 的模型——天气温度和状态(下雨、晴天等)。)

我们应该看到回归的效果，并比较这三个模型。

预测是为未来 10 天计算的。数据集中的最后一天是 2012 年 12 月 31 日，这意味着预测从 2013 年 1 月 1 日开始。

A .无附加回归器的模型

从 2013 年 1 月 1 日开始的自行车租赁预测值:

Bike rentals forecast for the next ten days

B .带附加回归器的模型-天气温度

回归器的值必须在过去和将来都是已知的，这是它如何帮助先知调整预测的。未来值必须是预定义的和已知的(例如，在特定日期发生的特定事件)，或者应该在其他地方预测。就我们的情况而言，我们使用的是天气温度预报，可以从天气频道获得。

回归器的值必须与时间序列数据位于同一数据框中。我复制了一个要定义为索引列的日期列:

df = pd.read_csv('weather_day.csv')
df = df[['dteday', 'cnt', 'temp']].dropna()d_df['date_index'] = d_df['dteday']
d_df['date_index'] = pd.to_datetime(d_df['date_index'])
d_df = d_df.set_index('date_index')

我们需要构建未来十天的数据框架——从 2013 年 1 月 1 日到 2001 年 10 天创建熊猫数据框架，并使用温度预测(标准化值)初始化每个元素:

t = 13
min_t = -8
max_t = 39
n_t = (t - min_t)/(max_t - min_t)
print(n_t)future_range = pd.date_range('2013-01-01', periods=10, freq='D')
future_temp_df = pd.DataFrame({ 'future_date': future_range, 'future_temp' : 0})future_temp_df['future_date'] = pd.to_datetime(future_temp_df['future_date'])
future_temp_df = future_temp_df.set_index('future_date')future_temp_df.at['2013-01-01', 'future_temp'] = 0.319148
future_temp_df.at['2013-01-02', 'future_temp'] = 0.255319
future_temp_df.at['2013-01-03', 'future_temp'] = 0.234042
future_temp_df.at['2013-01-04', 'future_temp'] = 0.319148
future_temp_df.at['2013-01-05', 'future_temp'] = 0.340425
future_temp_df.at['2013-01-06', 'future_temp'] = 0.404255
future_temp_df.at['2013-01-07', 'future_temp'] = 0.361702
future_temp_df.at['2013-01-08', 'future_temp'] = 0.404255
future_temp_df.at['2013-01-09', 'future_temp'] = 0.425531
future_temp_df.at['2013-01-10', 'future_temp'] = 0.446808future_temp_df.tail(10)

在下面的代码中，我添加了天气温度的回归器。如果日期属于训练集，则从训练集返回温度，否则从未来预测数据框(上面构建的数据框)返回温度。未来十天的期限已经设定。先知模型是用拟合函数构建的，调用预测函数来计算预测:

def weather_temp(ds):
    date = (pd.to_datetime(ds)).date()

    if d_df[date:].empty:
        return future_temp_df[date:]['future_temp'].values[0]
    else:
        return (d_df[date:]['temp']).values[0]

    return 0m = Prophet()
m.add_regressor('temp')
m.fit(d_df)future = m.make_future_dataframe(periods=10)
future['temp'] = future['ds'].apply(weather_temp)forecast = m.predict(future)
forecast[['ds', 'yhat', 'yhat_lower', 'yhat_upper']].tail(15)

从 2013 年 1 月 1 日开始，自行车租赁的预测值，以及附加的温度回归器。天气预报显示气温良好，预计 1 月份天气温暖，这有助于将自行车租赁数量调整到更高水平(如果天气温度良好，自然会有更多的租赁):

Bike rentals forecast for the next ten days, with temperature regressor

C .具有两个附加回归因子的模型-天气温度和条件

数据集包含相当多的描述天气的属性，使用更多的这些属性将有助于计算更准确的自行车租赁预测。我将展示增加一个回归变量如何改变预测。

def weather_temp(ds):
    date = (pd.to_datetime(ds)).date()

    if d_df[date:].empty:
        return future_temp_df[date:]['future_temp'].values[0]
    else:
        return (d_df[date:]['temp']).values[0]

    return 0def weather_condition(ds):
    date = (pd.to_datetime(ds)).date()

    if d_df[date:].empty:
        return future_temp_df[date:]['future_weathersit'].values[0]
    else:
        return (d_df[date:]['weathersit']).values[0]

    return 0m = Prophet()
m.add_regressor('temp')
m.add_regressor('weathersit')
m.fit(d_df)future = m.make_future_dataframe(periods=10)
future['temp'] = future['ds'].apply(weather_temp)
future['weathersit'] = future['ds'].apply(weather_condition)forecast = m.predict(future)
forecast[['ds', 'yhat', 'yhat_lower', 'yhat_upper']].tail(15)

使用上述代码添加了额外的回归变量—天气条件(数据集中的 check weathersit 属性)以及天气温度回归变量。出于测试的目的，我将未来十天的天气条件设置为 4(这意味着天气不好)。即使一月份的天气温度在上升(这对自行车租赁来说是好事)，整体天气还是很糟糕(这应该会减少自行车租赁的数量)。

使用第二个回归变量(指向预期的坏天气)，Prophet 返回较小的预期自行车租赁数量:

Bike rentals forecast for the next ten days, with temperature regressor and weather condition

总结:在 Prophet 中，附加回归量功能对于精确预测计算非常重要。它有助于调整预测的构建方式，并使预测过程更加透明。回归量必须是过去已知的变量，并且是已知的(或对未来单独预测的)。

资源:

• GitHub 回购与源代码
• 我以前的一篇关于同一主题的文章— 用烧瓶为先知模型服务—预测未来

用机器学习预测收入惊喜

原文：https://towardsdatascience.com/forecasting-earning-surprises-with-machine-learning-68b2f2318936?source=collection_archive---------16-----------------------

如何预测哪些公司将超过或低于分析师的收益预期

上市公司发布季度收益报告，当结果偏离分析师的估计时，会引起重大的价格波动。这是因为根据有效市场假说，资产价格完全反映了所有可获得的信息，因此会成为共识估计的一个因素。在本文中，我们将了解如何使用机器学习来预测一家公司是否会超出或低于其预期。

数据

我们考虑来自汤姆森路透 I/B/E/S 估计数据库的每股收益分析师估计，并从 Sentieo 下载。该数据库收集并汇编了分析师对 20 多项指标的估计。对于每家公司，我们给出了估计值的平均值、估计数、低值、高值和实际值，如下所示:

不幸的是，对于这个数据库，我们每个公司只有 70 个数据点，这不足以根据以前公布的结果和他们的 Beat/Miss vs 估计来预测一个公司的收益，但我们可以重新构建问题以增加数据点的数量。

不要问我们自己一个公司是否会超出或错过估计，我们可以问估计值是否会高于或低于实际值。

然后，我们将对这些值进行归一化，以便对它们进行聚合。在这种情况下，我们将为我们的模型考虑的特征是:

• #估计值
• 低/平均百分比
• 高/平均百分比
• 实际/平均百分比

然后，我们决定按部门汇总估计值，以检验分析师准确预测收益的能力与公司性质相关的假设。在这项研究中，我们将重点关注医疗保健股。

然后，我们对以下 117 家公司进行了 6000 次评估:

AAC ‘，’ ABT ‘，’ ABBV ‘，’ ACHC ‘，’ XLRN ‘，’ ACOR ‘，’ AERI ‘，’ AGIO ‘，’ AIMT ‘，’ AKCA ‘，’ AKBA’
，’ AKRX ‘，’ ALXN ‘，’ ALGN ‘，’ AGN ‘，’ ALNY ‘，’ AMRN ‘，’ AMGN’
，’ FOLD ‘，’ ARRY ‘，’ ASND ‘，’ AZN ‘，’ ATRC ‘，’ AVNS’
，’ BHC ‘，’ BAX ‘，’ BDX ‘，’ BCRX ‘，’ BMRN ‘，’ TECH ‘，’ BEAT ‘，’ BLUE’ 【T7 ‘，’ BSX ‘，’ BMY ‘，’

处理数据

我们将数据上传到 AuDaS ，这是一个由 Mind Foundry 为分析师打造的数据科学和教育平台。为了提高模型的准确性，我们创建了一个新列来表示实际值是高于(1)还是低于(-1)实际值(相对于 a %)。

我们还可以通过自动生成的直方图来可视化数据，并查看节拍/缺失如何跨其他特征分布。

构建我们的机器学习模型

由于我们预测节拍/缺失列，因此我们将构建一个分类器，同时排除实际/平均列。我们还将坚持由 AuDaS 自动选择的推荐训练配置，以防止过度拟合(10 倍交叉验证和 10%坚持验证)。

然后，AuDaS 将开始使用 Mind Foundry 的专有优化器 OPTaaS 搜索最佳机器学习模型，该优化器已迅速成为量化基金行业最受欢迎的快速全局优化和超参数调整工具。

在不到一分钟的时间里，AuDaS 已经尝试了 37 种不同的机器学习模型，找到的最佳解决方案是一个简单的梯度提升分类器。将鼠标悬停在模型上可显示其参数值:

该模型的相对特征相关性表明低/均值、高/均值比率包含最多的信息。

测试模型

然后，我们可以在 10%的保留率上测试模型:

AuDaS 实现了 69.4%的分类准确度，并且最终的模型健康建议是好的。

结论和延伸

凭借相对较少的功能，AuDaS 能够建立一个准确的预测模型，支持投资分析师对 IBES 估值的评估。这可以使他们预测重大的价格变动。

这项研究的延伸将是使用 AuDaS 的聚类功能对多个部门的估计值进行分组。这将使我们能够检验这样一个假设，即公司治理等其他因素将比行业更能影响预测的成功/失败比率。

如果您有兴趣了解我们的 Quant 和基本面对冲基金客户如何使用澳元进行风险分析和投资管理，请不要犹豫，通过电子邮件或 LinkedIn 联系我。您还可以阅读以下更多案例研究:

[## 集群风险管理

如何揭开表面下的结构

towardsdatascience.com](/risk-management-with-clustering-fd594e064806) [## 机器学习的价值投资

你最喜欢的持有期不一定是永远…

towardsdatascience.com](/value-investing-with-machine-learning-e41867156108) [## 用数据科学充实投资分析师

基本面投资如何受益于机器学习

towardsdatascience.com](/augmenting-investment-analysts-with-data-science-98297cb1ccb0)

注来自《走向数据科学》的编辑: 虽然我们允许独立作者根据我们的 规则和指导方针 发表文章，但我们不认可每个作者的贡献。你不应该在没有寻求专业建议的情况下依赖一个作者的作品。详见我们的 读者术语 。

更新:我开了一家科技公司。你可以在这里找到更多的

团队和资源

Mind Foundry 是牛津大学的一个分支，由斯蒂芬·罗伯茨和迈克尔·奥斯本教授创立，他们在数据分析领域已经工作了 35 年。Mind Foundry 团队由 30 多名世界级的机器学习研究人员和精英软件工程师组成，其中许多人曾是牛津大学的博士后。此外，Mind Foundry 通过其分拆地位，拥有超过 30 名牛津大学机器学习博士的特权。Mind Foundry 是牛津大学的投资组合公司，其投资者包括牛津科学创新、牛津技术与创新基金、牛津大学创新基金和 Parkwalk Advisors 。

利用历史数据预测加密货币的未来价格

原文：https://towardsdatascience.com/forecasting-future-prices-of-cryptocurrency-using-historical-data-83604e72bc68?source=collection_archive---------10-----------------------

加密货币(cryptocurrency)这个词已经席卷了金融界，然而目前还缺乏对数字资产数据的正式和公开的研究。这篇博文旨在建立在这篇 文章 中进行的探索性数据分析的基础上。我们将在这里使用可用的数据集来预测不同加密货币的未来价格。

在 EDA 之后，预测是下一步，你需要预测价格的未来值。对于有兴趣投资这些加密货币的人来说，这可能具有巨大的商业价值。

我们将尝试通过使用 closing_price 特性来预测比特币的未来价格。

用什么型号？

为了进行预测，我们需要一个机器学习模型。大多数人想预测值的时候都会想到多元线性回归。但是对于时间序列数据，这不是一个好主意。不选择时间序列数据回归的主要原因是我们对预测未来感兴趣，这将是线性回归的外推(预测数据范围之外)。我们知道，在线性回归中，任何形式的外推都是不可取的。

对于时间序列数据，最好使用自回归综合移动平均或 ARIMA 模型。

ARIMA

ARIMA 实际上是一类模型，它基于给定的时间序列的过去值来“解释”给定的时间序列，即其自身的滞后和滞后预测误差，因此该方程可用于预测未来值。任何“非季节性”的时间序列显示模式，不是一个随机的白噪声可以用 ARIMA 模型建模。如下所述，假设检验表明价格不是季节性的，因此我们可以使用 ARIMA 模型。在我们的 EDA 阶段，我们还发现任何加密货币都没有系统模式。

Non-Seasonal Trend of Omisego’s Prices

要创建 ARIMA 模型，我们需要 3 个参数:

• p；自回归项的阶
• 问；移动平均项的阶
• d；使时间序列平稳所需的差分次数

确定平稳性

我们将需要时间序列数据是平稳的，因为术语自动回归在 ARIMA 意味着它是一个线性回归模型，使用它自己的滞后作为预测器。如你所知，当预测因子不相关且相互独立时，线性回归模型效果最佳。

为了实现平稳性，我们需要执行差分:从当前值中减去前一个值。有时，根据系列的复杂程度，可能需要一个以上的差异。在我们这样做之前，我们需要检查我们的数据是否已经是稳定的。为此，我们将使用 增强的 Dickey Fuller 测试。这是一个假设检验**，假设零假设为*“自回归模型中存在单位根”*。单位根检验时间序列变量是否是非平稳的，是否有单位根。如果一个数列有一个单位根，它就显示出一个系统的模式。在这种情况下，另一个假设是“没有单位根存在”，这意味着时间序列是平稳的。**

result = adfuller(df["Close]).dropna())

测试返回一个为 0.002114 的 P 值。这远低于我们设定的 0.05 的显著性水平，因此我们拒绝零假设，并得出数据是平稳的结论。

由于数据是静态的，我们不需要执行任何差分。因此，参数 d 的值将为 0。

确定 AR 项的阶数

p 是自回归 (AR)项的顺序。它是指用作预测值的 Y 的滞后数。我们可以通过查看偏相关 函数 (PACF)图来找出所需的 AR 项数。在排除中间滞后的贡献后，部分自相关可以被想象为序列和它的滞后之间的相关性。**

我们可以观察到，PACF 滞后 1 是唯一一个非常显著的值，与其他值相比，它远远高于显著性线。因此，我们可以安全地将 p 设置为 1。

确定 MA 术语的顺序

q 是移动平均线 (MA)项的顺序。它是指应该进入 ARIMA 模型的滞后预测误差的数量。MA 术语在技术上是滞后预测的误差。**

就像我们看 PACF 图中 AR 项的数量一样，你也可以看 ACF 图中 MA 项的数量。ACF 告诉我们需要多少 MA 项来消除平稳序列中的任何自相关。

计算显著线以上的滞后数。所以，我们暂定 q 为 8。

构建模型

既然我们已经有了所有需要的参数，现在我们将构建模型。

model = ARIMA(df["Close"], order=(1,0,8))
model_fit = model.fit(disp=0)
print(model_fit.summary())

ARIMA MODEL FOR p=1, d=0, q=8

我们可以注意到，MA 项的 P > |z|中的 P 值远小于 0.05。因此，我们可以理想地去掉一些 MA 项并重建模型。由于三项显著高于显著性水平，我们将 q 设为 3。

model = ARIMA(df["Close"], order=(1,0,3))

ARIMA Model for revised parameters

评估模型

既然模型已经被拟合，我们现在将测试我们的 ARIMA 模型。为此，我们现在将数据集分为两部分:

• 训练数据(原始数据集的 80%)
• 测试数据(原始数据集的 20%)

X = bit_df["Close"].values
train_size = int(len(X) * 0.80)
predictions = model_fit.predict(train_size,len(X)-1)
test = X[train_size:len(X)]
error = mean_squared_error(test, predictions)

model_fit.predict()获取测试数据的开始和结束索引。然后，通过计算均方误差，将结果与测试数据集的目标值进行比较。

均方差的值出来是~ 575.24 。当我们看到数据跨越了 5 年的每日条目时，这个误差并不大。然而，仅仅是均方误差值并不能让我们清楚地知道我们的模型有多精确。因此，我们将可视化测试数据的预测值和实际值的比较。

Actual closing values VS Predicted closing values

查看图表，我们可以看到收盘价的预测值和实际值非常接近。

那都是乡亲们！祝你有愉快的一天。

预测:如何检测异常值？

原文：https://towardsdatascience.com/forecasting-how-to-detect-outliers-cb65faafcd97?source=collection_archive---------5-----------------------

对于需求计划者来说，异常值可能是一个痛苦。让我们来看看是什么导致了它们，以及标记它们的 3 种技术。

下面这篇文章摘自我的《供应链预测的数据科学》一书，这里有。你可以在这里* 找到我的其他出版物 。我也活跃在LinkedIn*。

“今天我不想进一步尝试去定义这种材料(……)，也许我永远也不能明白地做到这一点。但我一看就知道了。”
波特·斯图尔特

1964 年，波特·斯图尔特是美国最高法院法官。他不是在讨论局外人，而是在讨论电影《情人》是否淫秽。

当您进行预测时，您会注意到您的数据集会有异常值。尽管当我看到它时我就知道它可能是唯一实用的定义，但这些离群值对供应链构成了真正的威胁。这些高点(或低点)将导致你的预测或安全库存的过度反应，最终导致(最好的情况下)人工修正或(最坏的情况下)死库存、损失和严重的牛鞭效应。实际上，当你查看博客、书籍、文章或预测软件时，异常值检测的问题经常被回避。这是一个遗憾。离群点检测是一件严肃的事情。

这些异常值在现代供应链中无时无刻不在出现。它们主要是由于两个主要原因:

错误&错误这些都是明显的离群值。如果您发现了此类错误或编码错误，就需要流程改进来防止这些错误再次发生。
异常需求尽管一些需求观察是真实的，但这并不意味着它们不是异常并且不应该被清理或平滑。这种非凡的销售在供应链中其实并不罕见。想想促销、营销、奇怪的客户行为或减少库存。通常情况下，您可能不想在预测中考虑去年为处理几乎过时的旧库存而实现的 80%的特殊销售额。

如果你能发现异常值并消除它们，你就能做出更好的预测。我见过许多例子，仅仅由于异常值清理，预测误差就减少了百分之几。实际上，数据集越大，自动检测和清理就越重要。让我们看看我们如何能做到这一点。

在这篇文章中，我们将讨论三个半想法来找出这些异常值，并把它们放回到一个合理的水平。

想法 1——Winsorization

正如我们所说，异常值是异常高或异常低的值。基于这个简单的定义，检测异常值的第一个想法是简单地删除数据集的前 x 个最高点和最低点。让我们看看这将如何在下表中的两个(虚拟)数据集上工作。

**

第一种方法只是将历史需求的上/下 x%值降低到第 x 个百分点的极限。

第 x 百分位是一个值，低于该值，一个组中的 x%的观察值将下降。例如，一个产品 99%的需求观察值将低于它的第 99 个百分位数。

这种简单地将需求缩减到特定百分比的技术被称为 winsorization 。这个名字来自 20 世纪上半叶的统计学家查尔斯·p·温索尔。

如果我们查看上面两个虚拟数据集的第 1 和第 99 个百分位数，我们会得到以下结果:

在该表中，我们看到在两个数据集中，所有低值都将增加到 4.4。你可以在下图中看到，这剪切了我们数据集的一部分。对于没有异常值的数据集，高值将降至 16.6(见图 10.1)，对于有异常值的数据集，高值将降至 70.9(见图 10.2)。

你可能已经注意到，winsorization 没有给我们 4 或 5 这样的整数结果，但是我们得到了 4.4。实际上，由于我们没有一个精确的值可以将数据集削减 99%，所以我们基于两个最接近的值进行线性近似。这就是我们如何得到这些数字而不是整数的。

**

那么，我们对这种技术满意吗？不，我们不是。
-我们在没有异常值的数据集上发现了假异常值。
-在有离群值的数据集上，我们没有充分减少离群值(它从 100 降到了 70.9)。

当然，可以简单地提议将 winsorization 的上限从 99%降低到 95%，以进一步减少数据集#2 上的离群值。然而，不幸的是，这也会影响数据集#1。这不是一个好的解决办法。也可以提议取消这个下限，这样我们就不会将需求增加到 4.4。但是，如果我们有需求缺失的时期呢？如果有的话，我们不应该也清理一下吗？

自己动手
Excel 使用公式=PERCENTILE 可以快速得到 Excel 中某一范围单元格的不同百分位数。INC(范围，极限)。当然，您必须对上限(值大约为 0.95-0.99)使用一次该公式，对下限(值大约为 0.01-0.05)使用一次。多亏了 NumPy，我们可以轻松地用 Python 编写数据集。借助函数 np.percentile(array，percentile) ，我们可以计算数组的不同百分位数。

*import numpy as np
higher_limit = np.percentile(array, 99)
lower_limit = np.percentile(array, 1)*

请注意， percentile 函数采用表示为 0 到 100 之间的值的百分点，而不是像 Excel 中那样的比率(即 0 到 1 之间的值)。

由于函数 np.clip(array，min，max) ，我们可以简单地将数组剪切到这些下限和上限:

*array = np.clip(array,a_min=lower_limit,a_max=higher_limit)*

想法 2 标准偏差

正如我们刚刚看到的，winsorization 并不是排除异常值的完美方法，因为它会取出数据集的高值和低值，即使它们并不异常。

另一种方法是查看历史平均值周围的需求变化，并排除异常地远离该平均值的值。
让我们将需求标准差定义为:

其中 n 是我们拥有的需求观察的数量。

如果我们假设我们的数据正态分布在历史平均值附近，我们可以计算需求在两个阈值之间的概率。这里涉及的精确数学超出了本文的范围，不幸的是，通常情况下，正态性假设并没有得到严格的尊重。这两个阈值将以需求平均值(μ)为中心，在两个方向上的分布为 x 乘以标准差(σ)。需求越混乱(即σ越重要)，阈值就越宽。

**

例如，我们有 98%的概率处于范围:需求平均值+/- 2.33 倍标准差(如上图所示)。因此，如果我们想要移除高值和低值的前 1%,我们会将需求限制在μ +/-2.33 σ。

请注意，这意味着我们有 99%的概率比μ + 2.33 σ低*。并且有 99%的概率比μ-2.33σ高 T25。***

如果我们将此应用于我们的示例数据集(参见前两个表)，我们将得到这些限制:

让我们看看这些新的普通限制与 winsorization 限制相比表现如何。**

这已经比我们用 winsorization 得到的结果好太多了:
——在没有离群值的数据集上(见图 10.4)，我们不改变任何需求观察(完美！—正如我们所愿)。
-在有异常值的数据集上，我们不改变低需求点，只改变实际的异常值(见图 10.5)。

尽管如此，即使我们将离群值减少到比 winsorization (70.9)更易管理的数量(47.9)，这可能还不够。

那么，我们现在幸福吗？
还不完全是。

你可能还记得，我们假设误差在历史平均值左右。这对于需求平淡的产品来说没什么问题，但是当你有一个有趋势或季节性的产品时，实际的限制就会出现。例如，在下面的季节表中，最高点(或最低点)不再是您想要移除的异常值。

在下图中，您可以看到 winsorization 和 normalization 如何影响季节性需求。

这根本没有意义:这两种技术都将季节峰值标记为异常值，并且它们跳过了真正的异常值，即 Y2 M11。

我们将用下一项技术来解决这个问题。

自己动手
Excel 借助公式 =STDEV，你可以计算一系列单元格的标准差。p(范围)。和往常一样，你可以通过 =AVERAGE(range) 来计算平均值。一旦你有了这两个，你就可以通过 =NORM 来计算上限和下限。INV(百分位数、平均值、标准偏差)。通常，您会希望上百分位在 0.99 左右，下百分位在 0.01 左右。
Python 你可以通过 np.std(array) 计算一个类似数组的(例如，一个列表，一个数据帧等)的标准差。)或直接通过方法获取数据帧。std()* 。因此，如果您有一个*数据帧 df，您可以简单地键入:

***m = df.mean()
s = df.std()***

然后我们将再次使用 SciPy 库来计算正态概率。然后我们将使用。在我们的数据帧上剪辑方法来限制它。**

***from scipy.stats import norm
#Print the probabilities of each demand observation
print(norm.cdf(df.values, m, s).round(2))limit_high = norm.ppf(0.99,m,s)
limit_low = norm.ppf(0.01,m,s)
df = df.clip(lower=limit_low, upper=limit_high)***

想法 3 误差标准偏差

我们必须标记异常值的第二个想法是将每个观察值与需求的平均值进行比较。我们看到，如果我们有一个趋势或季节性，这是没有意义的，因为观察值和历史平均值之间的差异是不相关的。

好了，让我们回到离群值的定义:离群值是你没有预料到的值，就像《巨蟒之灾》里的西班牙宗教裁判所显示的那样。也就是说离群值是一个远离你的预测(即你的预测)的值。因此，为了发现异常值，我们将分析预测误差，并查看哪些时段异常错误。为此，我们将使用之前使用的标准差方法。

让我们拿回上面做的季节性例子。我们将历史需求与我们对它的简单(但季节性的)预测进行比较。

如果我们计算这种预测的误差(它只是历史需求的平均值)，我们将获得 0.4 的平均误差和 3.2 的标准偏差(这当然受到 Y2 M11 误差的严重影响)。如果我们在这个平均值附近取 99%的置信区间，我们会将预测误差缩小到-0.4+/-2.33 x 3.2 =-8.7。您可以在下图中看到预测的这些限制是如何完美地符合季节性需求的。

我们现在可以纠正 Y2 M11 的异常值。该期间的需求为 19，但预测为 5。可接受的最大值为 5 + 7 = 12。这意味着我们可以用这个新值(12)替换 Y2 M11 (19)的异常值。**

结论如上图所示，标准化和 winsorization 无法为季节性需求带来任何有意义的结果。

**这个方法的微调，你应该取多少个标准差作为极限？当然，这是留给你去试验的……

自己动手
Python 如果您有一个 pandas DataFrame ，其中一列为预测，另一列为需求(我们的指数平滑模型的典型输出)，我们可以使用以下代码:

***df[“Error”] = df[“Forecast”] — df[“Demand”]
m = df[“Error”].mean()
s = df[“Error”].std()
from scipy.stats import normlimit_high = norm.ppf(0.99,m,s)+df[“Forecast”]
limit_low = norm.ppf(0.01,m,s)+df[“Forecast”]
df[“Updated”] = df[“Demand”].clip(lower=limit_low,upper=limit_high)
print(df)***

多走一步！

如果你回想一下我们分析预测误差并设定可接受误差阈值的想法，我们实际上还有一个小问题。我们计算的阈值基于数据集，包括离群值。这种异常值驱使误差向上变化，使得可接受的阈值被偏置和高估。为了纠正这一点，可以不将异常值缩小到基于原始需求数据集计算的阈值，而是缩小到基于没有该特定异常值的数据集计算的限制。食谱如下:**

1. 根据历史需求填充第一个预测。
2. 计算误差、误差均值和误差标准差
3. 计算可接受的下限和上限阈值(基于误差平均值和标准偏差)。
4. 如前所述，识别异常值。
5. 重新计算误差平均值和标准偏差，但排除异常值。
6. 基于这些新值更新可接受的下限和上限阈值。
7. 基于新的阈值更新异常值。

如果我们回到上面的季节性例子，我们最初的预测误差平均值为 0.4，标准偏差为 3.22。如果我们去掉点 Y2 M11，我们得到的误差平均值为-0.1，标准偏差为 2.3。这意味着现在预测的阈值是-5.3，5.2。然后，我们在 Y2 M11 中的异常值将被更新为 10(而不是我们之前技术中的 12)。

自己动手
我们将从之前的想法中收回代码，并添加一个新步骤来更新误差均值和标准偏差值。

***df[“Error”] = df[“Forecast”] — df[“Demand”]
m = df[“Error”].mean()
s = df[“Error”].std()from scipy.stats import norm
prob = norm.cdf(df[“Error”], m, s)
outliers = (prob > 0.99) | (prob < 0.01)m2 = df[“Error”][~outliers].mean()
s2 = df[“Error”][~outliers].std()limit_high = norm.ppf(0.99,m2,s2)+df[“Forecast”]
limit_low = norm.ppf(0.01,m2,s2)+df[“Forecast”]
df[“Updated”] = df[“Demand”].clip(lower=limit_low,upper=limit_high)
print(df)***

关于作者

*** [## Nicolas vande put——顾问，创始人——供应链| LinkedIn

查看 Nicolas Vandeput 在世界上最大的职业社区 LinkedIn 上的个人资料。尼古拉斯有 7 份工作列在…

www.linkedin.com](https://www.linkedin.com/in/vandeputnicolas/)

N icolas Vandeput 是一名供应链数据科学家，擅长需求预测和库存优化。他在 2016 年创立了他的咨询公司 SupChains ，并在 2018 年共同创立了 SKU 科学——一个快速、简单、实惠的需求预测平台。尼古拉斯对教育充满热情，他既是一个狂热的学习者，也喜欢在大学教学:自 2014 年以来，他一直在比利时布鲁塞尔为硕士学生教授预测和库存优化。自 2020 年以来，他还在法国巴黎的 CentraleSupelec 教授这两个科目。2018 年出版了 【供应链预测的数据科学】(2021 年第 2 版)，2020 年出版了 库存优化:模型与模拟 。

***

用 Python 预测脸书预言家

原文：https://towardsdatascience.com/forecasting-in-python-with-facebook-prophet-29810eb57e66?source=collection_archive---------2-----------------------

如何使用领域知识调整和优化 Prophet，以便更好地控制您的预测。

我写了一本关于脸书先知的书，已经由 Packt 出版社出版了！这本书在亚马逊上有售。

这本书涵盖了使用 Prophet 的每个细节，从安装到模型评估和调整。十几个数据集已经可用，并用于演示 Prophet 功能，从简单到高级，代码完全可用。如果你喜欢这篇中帖，请考虑在这里订购:【https://amzn.to/373oIcf】T4！在超过 250 页的篇幅中，它涵盖的内容远远超过了媒体所能教授的内容！

非常感谢你支持我的书！

被困在付费墙后面？点击此处阅读完整故事，并附上好友链接！

我是 Greg Rafferty，湾区的数据科学家。这个项目的代码可以在我的 GitHub 上找到。

在这篇文章中，我将解释如何使用脸书的先知进行预测，并演示一些通过使用领域知识来处理趋势不一致的高级技术。网上有很多 Prophet 教程，但是没有一个深入到调整 Prophet 模型，或者整合分析师知识来帮助模型导航数据。我打算在这篇文章中做到这两点。

https://www.instagram.com/p/BaKEnIPFUq-/

在之前一个关于 Tableau 中预测的故事中，我使用了一个 ARIMA 算法的修改版来预测美国商业航班的乘客数量。ARIMA 方法在处理静态数据和预测短期框架时表现不错，但脸书的工程师已经为 ARIMA 无法处理的情况开发了一个工具。Prophet 的后端是 STAN，一种概率编码语言。这使得 Prophet 拥有贝叶斯统计提供的许多优势，包括季节性、包含领域知识和置信区间，以添加数据驱动的风险估计。

我将从三个数据来源来说明如何使用 Prophet，以及它的一些优点。如果你想跟着做，你首先需要安装 Prophet脸书的文档提供了简单的说明。我在本文中使用的笔记本提供了构建所讨论模型的完整代码。

飞机乘客

让我们从简单的开始。同样的航空乘客数据来自我的上一篇文章。Prophet 要求时间序列数据至少有两列:ds是时间戳，y是值。加载数据后，我们需要将其格式化为:

passengers = pd.read_csv('data/AirPassengers.csv')df = pd.DataFrame()
df['ds'] = pd.to_datetime(passengers['Month'])
df['y'] = passengers['#Passengers']

只需几行代码，Prophet 就可以制作一个和我之前制作的 ARIMA 模型一样复杂的预测模型。在这里，我打电话给 Prophet 做一个 6 年的预测(频率是每月，周期是 12 个月/年乘以 6 年):

prophet = Prophet()
prophet.fit(df)
future = prophet.make_future_dataframe(periods=12 * 6, freq='M')
forecast = prophet.predict(future)
fig = prophet.plot(forecast)
a = add_changepoints_to_plot(fig.gca(), prophet, forecast)

Number of passengers (in the thousands) on commercial airlines in the US

Prophet 将原始数据作为黑点，蓝线是预测模型。浅蓝色区域是置信区间。使用add_changepoints_to_plot功能添加红线；垂直的虚线是预测趋势变化的变化点，红色实线是去除所有季节性因素后的趋势。我将在本文中使用这种绘图格式。

这个简单的例子结束后，让我们继续讨论更复杂的数据。

Divvy 自行车共享

Divvy 是芝加哥的一家自行车共享服务公司。我以前做过一个项目，在那里我分析了他们的数据，并把它与从地下天气搜集的天气信息联系起来。我知道这个数据表现出很强的季节性，所以我认为这将是一个伟大的先知的能力展示。

Divvy 数据是基于每次乘坐的级别，因此为了对 Prophet 的数据进行格式化，我聚合到了每日级别，并为每天的“事件”列的模式(即天气条件:'not_clear', 'rain or snow', ‘clear', ‘cloudy', ‘tstorms', ‘unknown')、乘坐次数和平均温度创建了列。

格式化后，让我们看看每天的乘车次数:

因此，这些数据显然有季节性，而且趋势似乎是随着时间的推移而增加。有了这个数据集，我想演示如何添加额外的回归量，在这个例子中是天气和温度。让我们来看看温度:

它看起来很像之前的图表，但没有增加的趋势。这种相似性是有道理的，因为当天气晴朗温暖时，骑自行车的人会骑得更频繁，所以这两个图应该一前一后地上升和下降。

为了通过添加另一个回归量来创建预测，该额外的回归量必须具有预测期的数据。出于这个原因，我将 Divvy 数据缩短为一年，这样我就可以用天气信息来预测那一年。您可以看到，我还添加了美国的先知默认假期:

prophet = Prophet()
prophet.add_country_holidays(country_name='US')
prophet.fit(df[d['date'] < pd.to_datetime('2017-01-01')])
future = prophet.make_future_dataframe(periods=365, freq='d')
forecast = prophet.predict(future)
fig = prophet.plot(forecast)
a = add_changepoints_to_plot(fig.gca(), prophet, forecast)
plt.show()
fig2 = prophet.plot_components(forecast)
plt.show()

上述代码块创建了趋势图，如前面的航空乘客部分所述:

Divvy trend plot

和组件图:

Divvy component plot

组件图由 3 部分组成:趋势、假日和季节性。事实上，这三个部分的总和构成了模型的整体。如果你减去所有其他成分，趋势就是数据显示的。假日图显示了模型中包含的所有假日的影响。在 Prophet 中实现的假日可以被认为是不自然的事件，当趋势偏离基线，但一旦事件结束又会返回。正如我们将在下面探讨的，额外的回归变量就像假日一样，它们导致趋势偏离基线，只是趋势在事件发生后会保持变化。在这种情况下，假期都会导致乘客减少，如果我们意识到很多乘客都是通勤者，这也是有道理的。每周季节性因素显示，乘客量在一周内相当稳定，但在周末会急剧下降。这是支持大多数乘客都是通勤者这一理论的证据。我要注意的最后一点是，每年的季节性曲线是相当波动的。这些图是用傅立叶变换创建的，本质上是叠加的正弦波。显然，这种情况下的默认有太多的自由度。为了平滑曲线，我接下来将创建一个 Prophet 模型，关闭年度季节性，并添加一个额外的回归变量来解释它，但自由度更少。我还将继续在此模型中添加这些天气回归因素:

prophet = Prophet(growth='linear',
                  yearly_seasonality=False,
                  weekly_seasonality=True,
                  daily_seasonality=False,
                  holidays=None,
                  seasonality_mode='multiplicative',
                  seasonality_prior_scale=10,
                  holidays_prior_scale=10,
                  changepoint_prior_scale=.05,
                  mcmc_samples=0
                 ).add_seasonality(name='yearly',
                                    period=365.25,
                                    fourier_order=3,
                                    prior_scale=10,
                                    mode='additive')prophet.add_country_holidays(country_name='US')
prophet.add_regressor('temp')
prophet.add_regressor('cloudy')
prophet.add_regressor('not clear')
prophet.add_regressor('rain or snow')
prophet.fit(df[df['ds'] < pd.to_datetime('2017')])
future = prophet.make_future_dataframe(periods=365, freq='D')
future['temp'] = df['temp']
future['cloudy'] = df['cloudy']
future['not clear'] = df['not clear']
future['rain or snow'] = df['rain or snow']
forecast = prophet.predict(future)
fig = prophet.plot(forecast)
a = add_changepoints_to_plot(fig.gca(), prophet, forecast)
plt.show()
fig2 = prophet.plot_components(forecast)
plt.show()

趋势图看起来非常相似，所以我只分享组件图:

Divvy component plot with smooth annual seasonality and weather regressors

去年的趋势是向上的，而不是向下的！这是因为去年的数据显示平均气温较低，这使得乘客人数比预期的要少。我们还看到年度曲线变得平滑，还有一个额外的图:extra_regressors_multiplicative图。这显示了天气的影响。我们所看到的是意料之中的:乘客量在夏季增加，冬季减少，这种变化很大程度上是由天气造成的。我想再看一样东西，只是为了演示一下。我再次运行了上述模型，但这次只包括了rain or snow的回归变量。这是组件图:

Divvy component plot of just the effect of rain or snow

这表明，下雨或下雪时，每天的乘车次数会比其他时候少 1400 次。很酷，对吧！？

最后，我想按小时聚合这个数据集，以创建另一个分量图，即每日季节性。这个图看起来是这样的:

Divvy component plot for daily seasonality

正如 T2·里维斯指出的，凌晨 4 点是最不适合醒来的时间。显然，芝加哥的自行车骑手同意这一点。不过，当地在上午 8 点过后出现了一个高峰:早上的通勤者；下午 6 点左右是全球高峰:傍晚时分。我还看到午夜过后有一个小高峰:我喜欢认为这是人们从酒吧回家的路。Divvy 数据到此为止！让我们继续看 Instagram。

照片墙

脸书开发了 Prophet 来分析自己的数据。因此，在合适的数据集上测试 Prophet 似乎是公平的。我在 Instagram 上搜索了几个展现出我想探究的有趣趋势的账户，然后我从服务中搜集了三个账户的所有数据: @natgeo 、 @kosh_dp 和 @jamesrodriguez10 。

国家地理

https://www.instagram.com/p/B5G_U_IgVKv/

2017 年，我正在进行一个项目，在那里我注意到国家地理 Instagram 账户中的一个异常。对于 2016 年 8 月这一个月，每张照片的点赞数突然莫名其妙地急剧增加，但随后一个月过去就又回到了基线。我想把这个峰值建模为这个月的一次营销活动，以增加喜欢，然后看看我是否可以预测未来营销活动的效果。

以下是 Natgeo 每个帖子的点赞数。这种趋势明显在增加，而且随着时间的推移，差异也在增加。有许多异常值具有极高的喜欢度，但在 2016 年 8 月出现了一个高峰，在那个月发布的所有照片的喜欢度都远远高于周围的帖子:

我不想猜测为什么会这样，但是为了这个模型，让我们假设 Natgeo 的营销部门进行了一些为期一个月的活动，旨在增加喜欢。首先，让我们建立一个忽略这一事实的模型，这样我们就有了一个可以比较的基线:

Natgeo likes per photo over time

预言家似乎被那道钉弄糊涂了。它正试图将其添加到年度季节性成分中，从每年八月的峰值中的蓝色实线可以看出这一点。先知希望这是一个重复发生的事件。为了告诉 Prophet 2016 年发生了一些其他年份不会重复的特殊事情，让我们为这个月创造一个节日:

promo = pd.DataFrame({'holiday': "Promo event",
                      'ds' : pd.to_datetime(['2016-08-01']),
                      'lower_window': 0,
                      'upper_window': 31})
future_promo = pd.DataFrame({'holiday': "Promo event",
                      'ds' : pd.to_datetime(['2020-08-01']),
                      'lower_window': 0,
                      'upper_window': 31})promos_hypothetical = pd.concat([promo, future_promo])

promo数据帧仅包含 2016 年 8 月的活动，而promos_hypothetical数据帧包含 Natgeo 假设考虑在 2020 年 8 月推出的额外宣传片。当添加假日时，Prophet 允许有一个较低的窗口和一个较高的窗口，基本上是假日事件包含的天数，例如，如果您想将黑色星期五与感恩节一起包含，或者将平安夜与圣诞节一起包含。我在“假期”后添加了 31 天，将整个月都包含在事件中。下面是代码和新的趋势图。注意，我只是在调用先知对象时发送holidays=promo:

prophet = Prophet(holidays=promo)
prophet.add_country_holidays(country_name='US')
prophet.fit(df)
future = prophet.make_future_dataframe(periods=365, freq='D')
forecast = prophet.predict(future)
fig = prophet.plot(forecast)
a = add_changepoints_to_plot(fig.gca(), prophet, forecast)
plt.show()
fig2 = prophet.plot_components(forecast)
plt.show()

Natgeo likes per photo over time, with a marketing campaign in August 2016

太棒了。现在 Prophet 并没有增加每年 8 月份的愚蠢涨幅，而是在 2016 年出现了不错的涨幅。现在让我们再次运行这个模型，但是使用那个promos_hypothetical数据框架，来估计如果 Natgeo 在 2020 年运行一个相同的活动会发生什么:

Natgeo likes per photo over time with a hypothetical marketing campaign upcoming in 2020

这演示了在添加非自然事件时如何预测行为。例如，计划的商品销售可以是今年的模型。现在让我们转到下一个账户。

阿纳斯塔西娅·科什

https://www.instagram.com/p/BfZG2QCgL37/

Anastasia Kosh 是一名俄罗斯摄影师，她在自己的 Instagram 上发布古怪的自拍照，并为 YouTube 制作音乐视频。几年前我住在莫斯科时，我们是同一条街上的邻居；当时她在 Instagram 上有大约 1 万名粉丝，但在 2017 年，她的 YouTube 账户在俄罗斯迅速走红，她已经成为莫斯科青少年中的名人。她的 Instagram 账户呈指数级增长，很快接近 100 万粉丝。这种指数式增长对 Prophet 来说似乎是一个很好的挑战。

这是我们要建模的数据:

这是乐观增长的经典曲棍球棒形状，除了在这种情况下它是真实的！用线性增长对其建模，就像我们对上述其他数据建模一样，会导致不切实际的预测:

Anastasia Kosh likes per photo over time, with linear growth

这条曲线会一直延伸到无穷大。显然，一张照片在 Instagram 上获得的赞数是有上限的。理论上，这将等于服务上的唯一帐户的数量。但实际上，并不是每个客户都会看到或喜欢这张照片。这就是分析师的一点领域知识派上用场的地方。我决定用逻辑增长来建模，这需要先知被告知一个上限(先知称之为cap)和一个下限:

cap = 200000
floor = 0
df['cap'] = cap
df['floor'] = floor

通过我自己对 Instagram 的了解和一点点试错，我决定了 20 万个赞的上限和 0 个赞的下限。值得注意的是，Prophet 确实允许这些值被定义为时间的函数，所以它们不必是常数。在这种情况下，常量值正是我所需要的:

prophet = Prophet(growth='logistic',
                  changepoint_range=0.95,
                  yearly_seasonality=False,
                  weekly_seasonality=False,
                  daily_seasonality=False,
                  seasonality_prior_scale=10,
                  changepoint_prior_scale=.01)
prophet.add_country_holidays(country_name='RU')
prophet.fit(df)
future = prophet.make_future_dataframe(periods=1460, freq='D')
future['cap'] = cap
future['floor'] = floor
forecast = prophet.predict(future)
fig = prophet.plot(forecast)
a = add_changepoints_to_plot(fig.gca(), prophet, forecast)
plt.show()
fig2 = prophet.plot_components(forecast)
plt.show()

我将增长定义为逻辑增长，关闭所有季节性(在我的图中似乎没有太多季节性)，并调整了一些调整参数。我还添加了俄罗斯的默认假日，因为大多数 Anastasia 的追随者都在那里。当在df上调用.fit方法时，Prophet 看到了cap和floor列，并知道将它们包含在模型中。非常重要的是，当您创建预测数据框架时，您要将这些列添加到其中(这就是上面代码块中的future数据框架)。我们将在下一节再次讨论这个问题。但是现在我们的趋势图看起来真实多了！

Anastasia Kosh likes per photo over time, with logistic growth

最后，让我们看看我们的最后一个例子。

哈梅斯·罗德里格斯

https://www.instagram.com/p/BySl8I7HOWa/

哈梅斯·罗德里格斯是一名哥伦比亚足球运动员，在 2014 年和 2018 年世界杯上都有出色表现。他的 Instagram 账号从一开始就有稳定的增长；但是在研究之前的分析时，我注意到在两届世界杯期间，他的账户出现了突然而持久的粉丝激增。与国家地理杂志(National Geographic)账户中可以被建模为假期的峰值相反，罗德里格斯的增长在两场比赛后没有回到基线，而是重新定义了一个新的基线。这是完全不同的行为，需要不同的建模方法来捕捉它。

这是整个帐户生命周期中哈梅斯·罗德里格斯的每张照片的喜欢数:

仅用我们在本教程中使用的技术很难清晰地建模。他在 2014 年夏天的第一届世界杯期间经历了趋势基线的上升，然后在 2018 年夏天的第二届世界杯期间经历了峰值，并可能改变了基线。用默认模型来模拟这种行为并不十分有效:

James Rodríguez likes per photo over time

不是可怕的车型；它只是没有很好地模拟这两场世界杯比赛的行为。如果我们像上面对 Anastasia Kosh 的数据所做的那样，将这些比赛建模为假期，我们确实会看到模型的改进:

wc_2014 = pd.DataFrame({'holiday': "World Cup 2014",
                      'ds' : pd.to_datetime(['2014-06-12']),
                      'lower_window': 0,
                      'upper_window': 40})
wc_2018 = pd.DataFrame({'holiday': "World Cup 2018",
                      'ds' : pd.to_datetime(['2018-06-14']),
                      'lower_window': 0,
                      'upper_window': 40})world_cup = pd.concat([wc_2014, wc_2018])prophet = Prophet(yearly_seasonality=False,
                  weekly_seasonality=False,
                  daily_seasonality=False,
                  holidays=world_cup,
                  changepoint_prior_scale=.1)
prophet.fit(df)
future = prophet.make_future_dataframe(periods=365, freq='D')
forecast = prophet.predict(future)
fig = prophet.plot(forecast)
a = add_changepoints_to_plot(fig.gca(), prophet, forecast)
plt.show()
fig2 = prophet.plot_components(forecast)
plt.show()

James Rodríguez likes per photo over time, with holidays added for the World Cups

我仍然不喜欢模型适应变化趋势线的速度太慢，尤其是在 2014 年世界杯前后。这只是一个过于平稳的过渡。不过，通过添加额外的回归变量，我们可以迫使 Prophet 考虑一个突然的变化。

在这种情况下，我为每场锦标赛定义了两个时间段，赛中和赛后。以这种方式建模假设在锦标赛之前，将会有一个特定的趋势线，在锦标赛期间，该趋势线将会有一个线性变化，在锦标赛之后，将会有另一个变化。我将这些周期定义为 0 或 1，开或关，并让 Prophet 根据数据训练自己来学习幅度:

df['during_world_cup_2014'] = 0
df.loc[(df['ds'] >= pd.to_datetime('2014-05-02')) & (df['ds'] <= pd.to_datetime('2014-08-25')), 'during_world_cup_2014'] = 1
df['after_world_cup_2014'] = 0
df.loc[(df['ds'] >= pd.to_datetime('2014-08-25')), 'after_world_cup_2014'] = 1df['during_world_cup_2018'] = 0
df.loc[(df['ds'] >= pd.to_datetime('2018-06-04')) & (df['ds'] <= pd.to_datetime('2018-07-03')), 'during_world_cup_2018'] = 1
df['after_world_cup_2018'] = 0
df.loc[(df['ds'] >= pd.to_datetime('2018-07-03')), 'after_world_cup_2018'] = 1

请注意，我正在更新future数据帧，以包含以下这些“假日”事件:

prophet = Prophet(yearly_seasonality=False,
                  weekly_seasonality=False,
                  daily_seasonality=False,
                  holidays=world_cup,
                  changepoint_prior_scale=.1)prophet.add_regressor('during_world_cup_2014', mode='additive')
prophet.add_regressor('after_world_cup_2014', mode='additive')
prophet.add_regressor('during_world_cup_2018', mode='additive')
prophet.add_regressor('after_world_cup_2018', mode='additive')prophet.fit(df)
future = prophet.make_future_dataframe(periods=365)future['during_world_cup_2014'] = 0
future.loc[(future['ds'] >= pd.to_datetime('2014-05-02')) & (future['ds'] <= pd.to_datetime('2014-08-25')), 'during_world_cup_2014'] = 1
future['after_world_cup_2014'] = 0
future.loc[(future['ds'] >= pd.to_datetime('2014-08-25')), 'after_world_cup_2014'] = 1future['during_world_cup_2018'] = 0
future.loc[(future['ds'] >= pd.to_datetime('2018-06-04')) & (future['ds'] <= pd.to_datetime('2018-07-03')), 'during_world_cup_2018'] = 1
future['after_world_cup_2018'] = 0
future.loc[(future['ds'] >= pd.to_datetime('2018-07-03')), 'after_world_cup_2018'] = 1forecast = prophet.predict(future)
fig = prophet.plot(forecast)
a = add_changepoints_to_plot(fig.gca(), prophet, forecast)
plt.show()
fig2 = prophet.plot_components(forecast)
plt.show()

James Rodríguez likes per photo over time, with additional regressors

在这里，蓝线是我们应该关注的。红线显示的只是趋势，减去了额外的回归因素和假期的影响。看看蓝色趋势线在世界杯期间是如何急剧跳跃的。这正是我们的领域知识告诉我们会发生的行为！在罗德里格斯打入他在世界杯上的第一个进球后，突然有成千上万的新粉丝来到了他的账户上。让我们看一下分量图，看看这些额外的回归量有什么具体的影响:

James Rodríguez component plot for the World Cup regressors

这告诉我们，在 2013 年和 2014 年初，世界杯对罗德里格斯每张照片的点赞没有影响。在 2014 年世界杯期间，他的平均每张照片点赞数大幅上升，并在比赛结束后继续上升(这可以解释为他在赛事期间获得了如此多的活跃粉丝)。2018 年世界杯期间也发生了类似的事件，但没有那么引人注目，大概是因为此时已经没有多少足球迷发现他的账户并关注他了。

感谢你一直关注这篇文章！我希望你现在明白如何使用假期，线性与逻辑增长率，以及额外的回归变量来丰富你的预言家预测。脸书用 Prophet 构建了一个非常有价值的工具，将曾经非常困难的概率预测变成了一组简单的参数，并且有很大的调整余地。祝你预测好运！

预测未来许多天的价格

原文：https://towardsdatascience.com/forecasting-prices-for-5-days-ahead-2460406c4ea2?source=collection_archive---------10-----------------------

在这篇文章中，我将提出一些想法来预测未来 5 天冰糖的价格。

Forecasted values (orange) against real values (green)

在公司内部，有几个风险因素需要考虑和管理。例如，价格变化时，经理必须决定卖出的最佳时机；产品需求预测改善库存管理。评估的另一个重要因素是对经常性费用的预测，因为表达性价值会影响公司的营运资本。

通常是对价格、需求、费用等的预测。考虑到经理的领域知识以及影响决策的特定主题的基本信息。

使用预测模型时，数据在训练集和测试集之间分离，系数在训练数据中生成，并且通常在测试集中逐步测试。因此，只预测一个时间段，然后用真实值更新，依此类推，直到测试集结束。

在这篇文章中，将会创建一个 ARIMA 模型，在这个模型中，从训练集中获得的系数将以稍微不同的方式进行测试。将对连续 5 天进行预测，并用实际数据更新该系列(连续 5 天)，同样，将为 5 个周期创建新的预测，直到测试数据结束。

这种方法对于那些希望创建具有多个提前期的预测模型的人来说很有意思，因为它可以很容易地用于将预测的外推与真实数据进行比较，从而验证模型的主要特征是什么以及哪些方面可以改进。

该模型的主要目的是预测 Esalq 通知的 50 公斤冰糖未来 5 天的每日价格，数据可在此处下载****(请从巴西网站下载数据，英文版没有巴西雷亚尔(BRL)、的价格，所有代码的 Jupyter 笔记本可在我的 GitHub 中找到)。

从现在开始，我们将使用从我的帖子中摘录的以下步骤: 创建时间序列预测的基本原则 。(如果您还不熟悉创建这些预测的基础，我强烈建议您阅读一下):

• 将系列分为训练集和测试集，
• 将系列转换成静态的，
• 搜索相关的滞后，
• 创建模型，
• 在训练集上绘制比较图，
• 评估模型和测试集中的错误，
• 对模型提出可能的改进建议。

创建培训和测试集

我们有大约 16 年的价格历史，我将保留最后 250 个工作日(大约 1 年)来测试模型。在预测值和当前值之间将有 50 个 5 天的比较。请注意，所有以前的数据都将用于训练模型:

平稳性检验和验证相关滞后

为了创建时间序列预测，序列必须是平稳的。所以，它必须具备以下相对恒定的条件:均值、方差和自相关。

下面我将使用函数来绘制序列，它的分布，自相关，偏自相关。此外， Dickey Fuller 统计测试将用于检查平稳性:

#creating a function to plot the graph and show the test result:
def check_stationarity(y, lags_plots=48, figsize=(22,8)):
 “Use Series as parameter”

 y = pd.Series(y)
 fig = plt.figure()ax1 = plt.subplot2grid((3, 3), (0, 0), colspan=2)
 ax2 = plt.subplot2grid((3, 3), (1, 0))
 ax3 = plt.subplot2grid((3, 3), (1, 1))
 ax4 = plt.subplot2grid((3, 3), (2, 0), colspan=2)y.plot(ax=ax1, figsize=figsize, color=’teal’)
 ax1.set_title(‘Crystal Sugar Price’)
 plot_acf(y, lags=lags_plots, zero=False, ax=ax2, color=’teal’);
 plot_pacf(y, lags=lags_plots, zero=False, ax=ax3, method=’ols’, color=’teal’);
 sns.distplot(y, bins=int(sqrt(len(y))), ax=ax4, color=’teal’)
 ax4.set_title(‘Price Distribution’)plt.tight_layout()

 print(‘Dickey-Fuller test results:’)
 adfinput = adfuller(y)
 adftest = pd.Series(adfinput[0:4], index=[‘Statistical Test’,’P-Value’,’Used Lags’,’Observations Number’])
 adftest = round(adftest,4)

 for key, value in adfinput[4].items():
 adftest[“Critical Values (%s)”%key] = value.round(4)

 print(adftest)

为了验证训练数据上的值是否是平稳的，我们将使用 5%的 P 值作为基础，如果该测试的 P 值小于 5%,我们可以认为序列是平稳的:

当我们分析图形时，注意到价格有轻微的上升趋势，分布图显示数据不符合高斯正态分布，并且自相关的特征显示所有滞后都有显著的相关性(这是具有某种趋势的系列的标志)，最终统计值显示 P 值为 8.7%。

因此，这一系列不符合被认为是稳定的标准。

我将做第一个改变，以消除趋势，只关注每天的变化:

Dickey Fuller 测试在 P 值上返回零，这意味着该序列变得稳定。除此之外，在价格图中，趋势已经消失，价格相对稳定。****

在自相关图中，值逐渐减小，与没有第一差异的图不同。

部分自相关下降得更突然，表明前三个值对当前值有更大的影响。显然，在预测中没有显著的季节相关性需要考虑。

沿着这些思路，我们可以考虑差分序列中具有前 3 个滞后的自回归模型，因此我们将对该模型使用以下术语: ARIMA (3，1，0)

创建模型

最初，将基于训练数据创建模型，因此将生成 3 个自相关滞后的系数，并将用于测试测试集中连续 5 天的外推。

有几种方法可以推断未来，主要的有:

• 为要预测的每个特定日期创建一个模型，然后在最终模型中添加所有内容，
• 递归，预测第一天并以此值为基础预测第二天。

我们将使用本例中的最后一个模型，在这种情况下，我们将预测未来 5 天，与真实数据进行比较，将后者添加到模型中，以便进行新的外推，然后我们将计算模型的误差，并分析预测的相关事实。

训练模型

# Training the model
model = ARIMA(train, order=(3,1,0)).fit()
train_pred = model.predict()

现在模型已经定型，让我们比较下图中的实际数据和预测数据:

在训练集中，模型设法捕捉市场的主要时刻。唯一一次预测不能捕捉到变化的时候是当有更大的变化时。分析误差，RMSE 为 BRL 4.967/袋。

获取系数

训练模型后，我们可以获得差值的乘法系数以及模型的常数，这些值将作为外推未来 5 天预测的基础。

# Taking the constant and the coefficients of the lags to use in the test base:
const, l1, l2, l3 = model.paramsprint(f'Constant value {round(const, 4)}')
print(f'Cofficients of Lag 1: {round(l1,4)}, Lag 2: {round(l2,4)} and Lag 3: {round(l3,4)}')**Constant value 0.0063
Cofficients of Lag 1: 0.2332, Lag 2: 0.2297 and Lag 3: 0.2309**

下一步是使用这些乘法系数来测试测试集中的外推。为了更好地可视化，我们使用了一个函数将每个五天的预测与真实数据进行了比较。虽然我们可以在下面的图表中看到，这种预测的主要特点是它遵循正在协商的最新趋势。仍然有改进的余地，因为模型本身不能预测趋势的可能变化。

在我们的测试集上分析 RMSE:

¬test_error = sqrt(mean_squared_error(test['Crystal'], test['Crystal Pred']))
print(f'The RMSE on the test set was BRL {round(test_error,4)}/bag')**The RMSE on the test set was BRL 1.1522/bag**

下一部分将详细分析错误:

分析错误

将分别分析每一步的平均误差。通常，由于不确定性，误差(在本例中为 RMSE)会随着外推周期而增加:

The RMSE for each step is:
[0.48524694 0.69167404 0.81373202 1.00020173 1.12683735]

正如预测的那样，在每一步，RMSE 都会增加，这是由于未来的不确定性，而且从第二步开始的预测是基于预测值而不是真实值。

以下是一些有助于提高模型准确性的想法:

• 创建各种模型的合奏(LSTM，Garch，MLP，TBATS，CNN 等)。)并制定一个平均价格，
• 有些模型可能对某些特定步骤有更好的预测，您可以使用每个步骤中误差最小的模型，这样可以得到每个步骤都有不同模型且误差较小的最终模型。
• 分析不同的误差指标，如平均误差、偏差等。为了验证每种度量的具体特征，
• 测试可能的季节性滞后以预测周期性运动，
• 每周对数据重新取样，
• 添加相关的外部数据，
• 将数据标准化，
• 对数变换或使用 Box-Cox 变换。

用简单的模型，甚至人工智能技术来预测未来几个时期，无疑是一种极其重要的工具。

面对竞争日益激烈的市场，创建一个具有良好准确性的模型对公司来说可能是一项优势，因为它可以捕捉可能的市场运动，这对管理员的决策可能是必不可少的。

像这样的模型可以用作创造现金流、管理生产和库存以及识别市场机会的重要工具。

我希望你喜欢这篇文章！如果您有任何问题或需要更多信息，请随时联系我。下面是在 LinkedIn 上与我联系的链接。

来源:

** [## 基于机器学习的用电量多步时间序列预测

鉴于智能电表的兴起和太阳能等发电技术的广泛采用…

machinelearningmastery.com](https://machinelearningmastery.com/multi-step-time-series-forecasting-with-machine-learning-models-for-household-electricity-consumption/) [## 多步时间序列预测的 4 种策略

时间序列预测通常在只需要一步预测的情况下进行讨论。当你需要的时候…

machinelearningmastery.com](https://machinelearningmastery.com/multi-step-time-series-forecasting/)**