TowardsDataScience 博客中文翻译 2020（七百二十九）

原文：TowardsDataScience Blog

协议：CC BY-NC-SA 4.0

PyCaret 2.2 已经发布了，有什么新消息吗？

原文：https://towardsdatascience.com/pycaret-2-2-is-here-whats-new-ad7612ca63b?source=collection_archive---------3-----------------------

PyCaret 2.2 现在可以使用 pip 下载。https://www.pycaret.org

我们很高兴宣布 py caret 2.2—2020 年 10 月的更新。

PyCaret 是一个用 Python 编写的开源、低代码机器学习库，可以自动化机器学习工作流。它是一个端到端的机器学习和模型管理工具，可以加快机器学习实验周期，让你更有效率。

与其他开源机器学习库相比，PyCaret 是一个替代的低代码库，可以用来替换数百行代码。这使得实验快速有效。

发行说明:【https://github.com/pycaret/pycaret/releases

文档:https://pycaret.readthedocs.io/en/latest/

正在安装 PyCaret

安装 PyCaret 非常容易，只需要几分钟。我们强烈建议使用虚拟环境来避免与其他库的潜在冲突。参见下面的示例代码，创建一个 conda 环境并在该 conda 环境中安装 pycaret:

**# create a conda environment** 
conda create --name yourenvname python=3.6 **# activate environment** 
conda activate yourenvname **# install pycaret** 
pip install pycaret**# create notebook kernel linked with the conda environment** python -mipykernel install --user --name yourenvname --display-name "display-name"

PyCaret 的默认安装是 pycaret 的精简版本，它只安装这里列出的硬依赖项。要安装 pycaret 的完整版本，请使用以下代码:

**# install the full version**
pip install pycaret[full]

当您安装 pycaret 的完整版本时，这里列出的所有可选依赖项也会被安装。

安装每夜构建

PyCaret 进化非常快。通常，您希望获得最新的特性，但是希望避免从源代码编译 PyCaret 或者等待下一个版本。幸运的是，您现在可以使用 pip 安装 pycaret-nightly。

**# install the nightly build** pip install pycaret-nightly

**# or install the full version of the nightly build** pip install pycaret-nightly[full]

PyCaret 2.2 功能概述

👉支持 GPU 的培训

PyCaret 2.2 提供了使用 GPU 进行选择模型训练和超参数调整的选项。API 的使用没有变化，但是，在某些情况下，必须安装额外的库，因为它们没有安装默认的精简版本或完整版本。以下模型现在可以在 GPU 上训练。

• 极端梯度增压(无需进一步安装)
• CatBoost(无需进一步安装)
• 光线渐变增强器(需要安装 GPU:https://lightgbm.readthedocs.io/en/latest/GPU-Tutorial.html
• 逻辑回归、岭分类器、随机森林、K 近邻分类器、K 近邻回归器、支持向量机、线性回归、岭回归、套索回归、K-均值聚类和基于密度的空间聚类(需要累积量> = 0.15【https://github.com/rapidsai/cuml】T2

要在 GPU 上启用光照梯度增强机器，您必须安装 LightGBM 的 GPU 启用版本。官方的分步指南是这里的。

如果你使用 Google Colab，你可以为 GPU 安装 Light Gradient Boosting Machine，但首先，你必须卸载 Light GBM-CPU 版本。在此之前，请确保在您的 Colab 会话中启用了 GPU。使用以下代码安装支持 GPU 的 LightGBM:

**# uninstall lightgbm CPU**
pip uninstall lightgbm -y**# install lightgbm GPU**
pip install lightgbm --install-option=--gpu --install-option="--opencl-include-dir=/usr/local/cuda/include/" --install-option="--opencl-library=/usr/local/cuda/lib64/libOpenCL.so"

截至今天，谷歌 Colab 上还不支持 cuML 0.15。这在未来可能会改变，但现在，你可以免费使用预装了 cuML 0.15 的 blazingSQL 笔记本。

https://blazingsql.com/

登录到您的帐户后，启动 Python 3 笔记本并使用以下代码安装 pycaret:

**# install pycaret on blazingSQL**
!/opt/conda-environments/rapids-stable/bin/python -m pip install --upgrade pycaret

或者，如果你在本地机器上有 GPU，或者你计划使用任何其他带有 GPU 的云服务，你可以遵循 cuML 的官方安装指南。

假设安装成功，在 GPU 上训练模型唯一需要做的就是在初始化设置函数时启用 GPU。

**# import dataset** from pycaret.datasets import get_data
data = get_data('poker')**# initialize the setup** from pycaret.classification import *clf = setup(data, target = 'CLASS', use_gpu = True)

就是这样。您现在可以像在 CPU 上一样使用 pycaret 了。它会自动使用 GPU 进行模型训练，如果可能的话，还会使用 CPU 等效算法。甚至在开始训练之前，您可以通过使用以下命令来检查 GPU 上启用了哪些模型:

models(internal=True)[['Name', 'GPU Enabled']]

模型输出(内部=真)

CPU 与 GPU 的性能指标评测比较(以秒为单位)

基准比较 CPU 与 GPU，数据集:扑克(100K x 88)

👉超参数调谐

现在有了新的超参数调谐方法。在 PyCaret 2.1 之前，在 PyCaret 中调整模型的超参数的唯一方法是使用 scikit-learn 中的随机网格搜索。2.2 中添加的新方法有:

• scikit-learn(网格)
• sci kit-优化(贝叶斯)
• tune-sklearn(随机、网格、贝叶斯、超点、bohb)
• optuna(随机，热塑性弹性体)

为了使用这些新方法，增加了两个新参数**‘搜索 _ 库’和‘搜索 _ 算法’**。

**# train dt using default hyperparameters**
dt = create_model('dt')**# tune hyperparameters with scikit-learn (default)**
tuned_dt_sklearn = tune_model(dt)**# tune hyperparameters with scikit-optimize**
tuned_dt_skopt = tune_model(dt, search_library = 'scikit-optimize')**# tune hyperparameters with optuna**
tuned_dt_optuna = tune_model(dt, search_library = 'optuna')**# tune hyperparameters with tune-sklearn**
tuned_dt_tuneskl = tune_model(dt, search_library = 'tune-sklearn')

tune_model 输出(dt，默认超参数 AUC = 0.7401

搜索 _ 算法依赖于搜索 _ 库。以下搜索算法可用于相应的搜索库:

• scikit-learn → random(默认)，网格
• sci kit-优化→贝叶斯(默认)
• tune-sklearn →随机(默认)，网格，贝叶斯，超点，bohb
• optuna →随机，tpe(默认)

具有 partial_fit 属性的估计器也支持提前停止。在发行说明中阅读更多关于它的。

不同调谐器的基准比较

PyCaret 数据集:糖尿病中可用搜索库/搜索算法的基准比较

👉内存和性能改进

PyCaret 2.2 关注的是性能和功能。大量代码被重构，以在不影响用户体验的情况下改善内存占用和优化性能。

一个例子是所有的数字数据从以前的 64 位动态转换为 32 位，大大减少了内存占用。性能提升的另一个例子是跨所有功能的交叉验证，与以前按折叠顺序训练相比，现在跨多个内核自动并行化。

我们在著名的纽约出租车数据集的 500 万个样本行上比较了 PyCaret 所有发布版本的性能。下图比较了完成设置初始化所需的时间:

来自纽约出租车数据集的 5M 行的基准性能比较

所有的比较都是在具有 8 个 CPU 内核的 AMD64 机器上进行的。

👉添加自定义指标

现在，您可以完全自定义(添加或删除)在交叉验证期间评估的指标。这意味着您不再局限于 PyCaret 的默认模型评估指标。增加了三个新功能获取指标、添加指标、和移除指标。用法超级简单。参见示例代码:

**# import dataset** from pycaret.datasets import get_data
data = get_data('juice')**# initialize the setup** from pycaret.classification import *clf = setup(data, target = 'Purchase')**# check all metrics used for model evaluation**get_metrics()**# add Log Loss metric in pycaret**
from sklearn.metrics import log_loss
add_metric('logloss', 'LogLoss', log_loss, greater_is_better=False)**# compare baseline models**best = compare_models()

添加自定义指标后的 compare_models 输出

请注意，在 TT 被添加到 compare_models 得分网格之前，一个新列“log loss”(所有新指标被添加到右侧)，因为我们使用 add_metric 函数添加了指标。您可以使用 scikit-learn 中可用的任何指标，或者您可以使用 make_scorer 函数创建自己的指标。您可以使用以下命令删除度量:

remove_metric('logloss')

👉迭代插补

迭代插补是一种使用回归和分类估计器来输入缺失数据的技术，将每个特征建模为其他特征的函数。每个特征以循环方式进行估算，以前的预测用于新的预测。为了提高估算的质量，这个过程要重复几次。与简单的插补相比，它可以创建更接近真实值的合成值，但需要额外的处理时间。

忠于 PyCaret 的精神，用法非常简单:

**# initialize setup**
from pycaret.classification import *clf = setup(data, target = 'Class', imputation_type="iterative")

默认情况下，它将使用光梯度增强机器作为分类特征(分类)和数值特征(回归)的估计器，这可以使用设置中的分类 _ 迭代 _ 估算和数值 _ 迭代 _ 估算参数进行更改。

迭代插补与简单插补的基准比较

为了比较简单均值插补和迭代插补的结果，我们使用了包含大量缺失值的马疝气数据集。下图比较了使用不同插补方法的逻辑回归的表现。

安东尼·鲍姆的一篇博文:https://www . LinkedIn . com/pulse/iterative-attumation-py caret-22-安东尼-鲍姆/

与简单平均插补相比，使用 KNN 迭代插补作为分类和数值特征的估计量，平均 AUC 得分提高了 0.014 (1.59%)。要了解更多关于这个功能的信息，你可以在这里阅读完整的博文。

👉折叠策略

PyCaret 2.2 提供了定义折叠策略的灵活性。在 PyCaret 2.1 之前，您不能定义交叉验证策略。它使用**‘StratifiedKFold’进行分类，使用‘k fold’**进行回归，这限制了 PyCaret 在某些用例中的使用，例如时间序列数据。

为了克服这个问题，在设置功能中增加了一个新参数**‘fold _ strategy’****。它可以采用以下值:**

• kfold为 KFold CV
• stratifiedkfold为分层 KFold CV
• groupkfold为集团 KFold CV
• timeseries为 TimeSeriesSplit CV 或者
• 与 scikit-learn 兼容的自定义 CV 生成器对象。

👉比较保留集上的模型

**如果你之前用过 PyCaret，你一定熟悉它最常用的函数 **compare_models。此功能使用交叉验证来训练和评估模型库中所有可用估计器的性能。然而，问题是如果你正在处理非常大的数据集，compare_models 可能要花很长时间才能完成。原因是它适合模型库中每个估计器的 10 倍。对于分类，这意味着总共有 15 x 10 = 150 个评估者。

在 PyCaret 2.2 中，我们在 compare_models 函数中引入了一个新参数 cross_validation ，当该参数设置为 False 时，将对维持集上的所有指标进行评估，而不是进行交叉验证。而单独依赖维持指标可能并不可取，尤其是当数据集太小时。在处理大型数据集时，这无疑可以节省大量时间。

为了量化影响，我们比较了 compare_models 在两种场景中的性能(交叉验证=真，交叉验证=假)。用于比较的数据集是这里是 (45K x 50)

交叉验证(耗时 7 分 13 秒):

compare_models 的输出(cross_validation = True)

无交叉验证(耗时 1 分 19 秒):

compare_models 的输出(cross_validation = False)

👉自定义转换

谈到灵活性，这是一个全垒打。一个新的参数 custom_pipeline 被添加到 setup 函数中，该函数可以接受任何转换器并附加到 PyCaret 的预处理管道中。所有自定义变换都在 train_test_split 之后分别应用于每个 CV 折叠，以避免目标泄漏的风险。用法非常简单:

**# import dataset** from pycaret.datasets import get_data
data = get_data('juice')**# build custom transformations** from imblearn.over_sampling import SMOTE
from sklearn.decomposition import PCA
custom_pp = [("PCA",PCA()),("smote", SMOTE())]**# initialize setup** from pycaret.classification import *clf = setup(data, target = 'Purchase', custom_pipeline = custom_pp)

👉独立的训练和测试设备

这是期待已久的，也是自第一次发布以来最受欢迎的特性之一。现在可以通过单独的测试集，而不是依赖 pycaret 内部的 train_test_split 。设置中增加了一个新参数‘测试数据’。当一个数据帧被传入 test_data 时，它被用作一个测试集，并且 train_size 参数被忽略。test_data 必须加标签。请参见下面的示例代码:

**# loading dataset**import pandas as pd
train_data = pd.read_csv('/path/train.csv')
test_data = pd.read_csv('/path/test.csv')**# initializing setup**from pycaret.classification import *
clf = setup(data = train_data, test_data = test_data)

👉禁用预处理

**如果您不想使用 PyCaret 的默认预处理管道，或者您已经有了转换后的数据集，而只想使用 PyCaret 的建模功能，这在以前是不可能的，但现在我们已经为您做好了准备。只需关闭设置中的’**预处理’参数。当预处理设置为 False 时，除了在 custom_pipeline 参数中传递的 train_test_split 和自定义转换之外，不应用任何转换。

**# initializing setup**from pycaret.classification import *
clf = setup(data = train_data, preprocess = False) 

然而，当在设置中关闭预处理时，您必须确保您的数据是建模就绪的，即没有缺失值、没有日期/时间戳、分类数据被编码等。)

👉其他变化

• 在 plot_model 中添加了新的图“提升”、“增益”和“树”。
• CatBoost 现在与 plot_model 功能兼容。它需要 catboost > = 0.23.2。
• 为了使使用和开发更容易，根据最佳实践，所有更新的 pycaret 函数都添加了类型提示。用户可以通过使用支持类型提示的 IDE 来利用它们。

要了解 PyCaret 2.2 中所有更新的更多信息，请参见发行说明。

使用 Python 中的轻量级工作流自动化库，您可以实现的目标是无限的。如果你觉得这有用，请不要忘记给我们 GitHub repo 上的⭐️。

想了解更多关于 PyCaret 的信息，请在 LinkedIn 和 Youtube 上关注我们。

重要链接

用户指南
文档
官方教程 示例笔记本
其他资源

想了解某个特定模块？

单击下面的链接查看文档和工作示例。

分类
回归
聚类
异常检测
自然语言处理 关联规则挖掘

PyCaret:用 Python 实现更好的机器学习

原文：https://towardsdatascience.com/pycaret-better-machine-learning-with-python-58b202806d1e?source=collection_archive---------19-----------------------

Scikit 的替代——使用 Python 学习机器

挑选任何数据科学书籍或在线课程。我敢打赌，唯一涵盖的机器学习库是 Scikit-Learn。这是一个很好的起点，不要误会我的意思，但在 2020 年，我们需要更加自动化和节省时间的东西。

安德烈·本兹在 Unsplash 上拍摄的照片

这就是 PyCaret 的用武之地。这是一个相对较新的库，所以没有太多的教程可用——尽管官方文档做了一项惊人的工作。下面是一句话的总结:

PyCaret 是 Python 中的一个开源、低代码机器学习库，允许您在几秒钟内从准备数据到在您选择的笔记本环境中部署模型。[1]

如果你更喜欢视频，或者只是想巩固你的知识，请随意观看我们关于这个主题的视频。源代码包含:

本文的目标读者是熟悉机器学习概念的人，他们也知道如何使用 Scikit-Learn 等库中的各种算法。理想的读者意识到自动化的需要，并且不想花费数天或数周时间寻找最佳算法及其超参数。

这将是一篇相当短的实践性文章，因为其目的只是展示该库能够做什么。未来的文章将更多地关注回归/分类/聚类任务。

好了，我们开始吧！

Python 演示

安装过程因操作系统而异。我在 macOS 上，所以在安装库之前需要一些步骤。

一旦完成，你就可以开始了。我们将从两个导入开始:

from pycaret.classification import * 
from pycaret.datasets import get_data

第一种方法是开始分类任务的推荐方法——尽管你们中的一些人可能不习惯使用 **import *** 语法。如果你不喜欢，就避开它。第二个导入允许我们使用内置数据集。

我们将使用内置的糖尿病数据集，以及 get_data() 函数来获取它。接下来，我们需要做一些设置，并告诉 PyCaret 目标变量是什么:

diabetes = get_data(‘diabetes’)
exp_clf = setup(diabetes, target=’Class variable’)

执行此代码将导致以下结果:

这个结果数据帧相当长，并且报告了大量关于数据的信息。我们现在可以继续机器学习部分。

寻找最佳算法

这一步很简单。键入以下内容:

compare_models()

这是结果:

是的， compare_models() 函数完成了所有这些工作。它还用黄色突出显示得分最高的单元格，这样更容易看到。逻辑回归似乎在准确性方面做得最好，但是 XGBoost 算法总体上执行得最好——所以这是我们将使用的算法。

*如何？*再一次，再简单不过了:

xgb = create_model(‘xgboost’)

在这个链接上，你会找到所有的模型缩写，所以你可以使用 XGBoost 以外的东西。接下来，让我们直观地看看我们的模型表现如何。

可视化模型性能

我们可以使用 plot_model() 函数来可视化模型性能:

plot_model(xgb)

它默认显示 ROC 曲线，但这很容易改变。这是混淆矩阵:

plot_model(xgb, ‘confusion_matrix’)

和分类报告:

plot_model(xgb, ‘class_report’)

*整齐，对吧？*这里是链接到你可以制作的所有可视化效果。让我们继续模型的可解释性。

SHAP 的模型解释

SHAP，或 *SHapley 附加解释，*是一种解释机器学习模型输出的方法。通过绘制每个样本的每个特征的 SHAP 值，我们可以使用它来查看哪些特征最重要。

PyCaret 又一次让事情变得非常简单:

interpret_model(xgb)

上图根据所有样本的 SHAP 值大小的总和对要素进行排序，并使用 SHAP 值来显示每个要素对模型输出的影响分布[2]。相当长的句子，所以通读几遍。颜色代表特征值，红色代表高，蓝色代表低。

简而言之，它指出，在口服葡萄糖耐量试验的 2 小时内血浆葡萄糖浓度水平较高(无论这意味着什么)，导致患糖尿病的几率较高(由于红色)。

我计划在其他文章中更深入地讨论 SHAP，所以我在这里就不赘述了。让我们继续在测试集上进行模型评估。

模型评估和预测

PyCaret 在加载时将数据分成训练和测试部分(70:30)，因此我们不必手动完成。

我们现在可以根据以前未见过的数据来评估模型:

predictions = predict_model(xgb)

上述代码会产生以下结果:

好了，现在当这个已经完成了，让我们看看如何保存和加载模型。

模型保存和加载

在保存模型之前，我们需要完成它:

finalize_model(xgb)

猜猜保存模型有多简单:

save_model(xgb, ‘diabetes_xgboost’)

模型以 pickle 格式保存。

尝试猜测加载已保存模型的函数的名称:

model = load_model(‘diabetes_xgboost’)

你准备再次使用它。

我们今天已经讲了很多东西，所以让我们在下一部分总结一下。更多类似的文章即将发布，敬请关注。

结论

今天我们简要介绍了 PyCaret 所提供的内容。我的想法是把这个变成一个系列，既有基于文章的，也有基于视频的，并且更深入地覆盖这个库。

当使用像这样的库时，像 Scikit-Learn 这样的库很快就会成为过去。不是因为他们不好，而是这些人在各方面都很优秀。

我希望你喜欢它，并设法从中获得一些有用的东西。感谢阅读。

喜欢这篇文章吗？成为 中等会员 继续无限制学习。如果你使用下面的链接，我会收到你的一部分会员费，不需要你额外付费。

[## 通过我的推荐链接加入 Medium-Dario rade ci

作为一个媒体会员，你的会员费的一部分会给你阅读的作家，你可以完全接触到每一个故事…

medium.com](https://medium.com/@radecicdario/membership)

参考

[1]https://pycaret.org

[2]https://github.com/slundberg/shap

PyCaret:机器学习的综合体

原文：https://towardsdatascience.com/pycaret-the-machine-learning-omnibus-dadf6e230f7b?source=collection_archive---------33-----------------------

…满足您所有机器学习需求的一站式商店

照片由阿什利·特雷瑟万在 Unsplash 拍摄

任何机器学习项目之旅都是从加载数据集开始，到结束(继续？！)最终确定用于预测未知数据和生产部署的最佳模型或模型集合。

作为机器学习的实践者，我们知道在达到最佳预测性能结果的过程中会有几个中途停留点。这些中间步骤包括探索性数据分析(EDA)、数据预处理——缺失值处理、异常值处理、改变数据类型、编码分类特征、数据转换、特征工程/选择、采样、训练-测试分离等。举几个例子，在我们开始建模、评估和预测之前。

我们最终导入了几十个 python 包来帮助我们做到这一点，这意味着要熟悉每个包中多个函数调用的语法和参数。

你是否希望有一个单一的包可以用一致的语法界面处理整个端到端的旅程？我当然有！

输入 PyCaret

这些愿望通过PyCaret包得到了回应，现在随着pycaret2.0的发布，它变得更加令人敬畏。

从这篇文章开始，我将发布一系列关于pycaret如何帮助我们快速完成一个 ML 项目的各个阶段的文章。

装置

安装轻而易举，几分钟内即可完成，所有依赖项也已安装完毕。建议使用虚拟环境安装，如 python3 virtualenv 或 conda environments 以避免与其他预安装的软件包冲突。

pip install pycaret==2.0

一旦安装完毕，我们就可以开始了！我们将包导入到我们的笔记本电脑环境中。这里我们将讨论一个分类问题。类似地，可以为涉及回归、聚类、异常检测、NLP 和关联规则挖掘的场景导入相应的 PyCaret 模块。

我们将使用来自 kaggle.com 的titanic数据集。你可以从这里下载数据集。

让我们使用head()函数检查数据集的前几行:

设置

pycaret的setup()函数只需一行代码就可以完成大部分——所有的修正工作，这些繁重的工作通常需要几十行代码才能完成！

我们只需要传递数据帧并指定目标特性的名称作为参数。setup 命令生成以下输出。

setup推断出了数据集中要素的数据类型，非常有用。如果我们同意，我们需要做的就是点击Enter。否则，如果您认为由setup推断的数据类型不正确，那么您可以在底部的字段中键入quit，然后返回到setup函数进行更改。我们将很快看到如何做到这一点。现在，让我们点击Enter，看看会发生什么。

输出续。,

产出续下

输出结束

咻！在一行看似无害的代码中，似乎发生了很多事情！让我们来盘点一下:

• 已检查缺失值
• 已识别的数字和分类特征
• 从原始数据集创建训练和测试数据集
• 具有平均值的连续要素中的估算缺失值
• 具有常数值的分类特征中的估算缺失值
• 完成标签编码
• …似乎还有许多其他选择，包括异常值处理、数据缩放、特征变换、降维、多重共线性处理、特征选择和处理不平衡数据等。！

但是嘿！第 11 和 12 行是什么？训练和测试数据集中的特征数量是 1745？似乎是标签编码变得疯狂的一个例子，很可能来自于分类特征，如name、ticket和cabin。在本文和下一篇文章中，我们将探讨如何根据我们的需求控制设置，以主动解决这种情况。

定制setup

首先，我们如何从模型构建中排除像上面三个特征这样的特征？我们在setup函数的ignore_features参数中传递我们想要排除的变量。需要注意的是，ID 和 DateTime 列在推断时会被自动设置为在建模时被忽略。

请注意下面的pycaret，在要求我们确认时，已经放弃了上面提到的 3 个特性。让我们点击Enter并继续。

在结果输出中(截断版本如下所示)，我们可以看到设置后，数据集形状现在更易于管理，只对剩余的更相关的分类特征进行标签编码:

在本系列的下一篇文章中，我们将详细讨论进一步的数据预处理任务，我们可以使用pycaret的这个setup函数通过传递额外的参数在数据集上实现这些任务。

但是在我们走之前，让我们用compare_model()函数来对pycaret惊人的模型比较能力做一个闪影。

在各种分类指标上比较模型性能。

嘣！所有需要做的只是compare_models()获得 15 个分类建模算法的结果，这些算法在交叉验证的各种分类指标之间进行比较。一目了然，我们可以看到CatBoost分类器在大多数指标上表现最佳，其中Naive-Bayes在召回率上表现良好，而Gradient Boosting在精确度上表现良好。每个指标中表现最佳的型号会由pycaret自动突出显示。

根据我们感兴趣的模型评估指标,pycaret有助于我们直接放大表现最佳的模型，我们可以使用超参数进一步调整这些模型。在接下来的文章中会有更多的介绍。

总之，我们已经简要了解了pycaret如何通过最少的代码结合关键数据预处理阶段的广泛而全面的定制，帮助我们快速完成 ML 项目生命周期。

您可能也会对我的其他文章感兴趣，这些文章介绍了在探索性数据分析(EDA)和可视化中使用最少的代码来交付最大的结果的优秀包。

[## 使用 SmartEDA 开发的 EDA

探索性数据分析——更智能、更快速的方式…

towardsdatascience.com](/eda-in-r-with-smarteda-eae12f2c6094) [## 拼凑——下一代 ggplots

进一步扩展 ggplot2 的多功能性…

towardsdatascience.com](/patchwork-the-next-generation-of-ggplots-1fcad5d2ba8a)

感谢您的阅读，并希望听到您的反馈。干杯！

PyCaret:机器学习综合教程—第 2 部分

原文：https://towardsdatascience.com/pycaret-the-machine-learning-omnibus-part-2-87c7d0756f2b?source=collection_archive---------31-----------------------

…满足您所有机器学习需求的一站式商店。

数据预处理和模型比较

阿什利·特雷瑟万在 Unsplash 上拍摄的照片

这是关于 python 中令人惊叹的pycaret包的系列文章的第二篇，它用少得可笑的代码行实现了对 ML 项目生命周期的几乎每个阶段的快速跟踪和自动化。如果您错过了第一部分，请单击下面的链接，我们在其中简要介绍了初始设置，它仅用一行代码就完成了数据预处理的所有方面，并将我们带到建模阶段。

[## PyCaret:机器学习的综合体

…满足您所有机器学习需求的一站式商店

towardsdatascience.com](/pycaret-the-machine-learning-omnibus-dadf6e230f7b)

在本文中，我们将查看几个可以传递给setup()函数的参数，以进一步控制由pycaret完成的预处理。默认情况下，setup函数只需要数据帧和我们想要预测其类别标签的目标特征。但是，由函数自动推断的要素数据类型可能并不总是正确的。在某些情况下，我们可能需要介入。例如，在我们正在使用的泰坦尼克号数据集中，setup函数正确地将Pclass(乘客等级)、SibSp(船上的兄弟姐妹&配偶)和Parch(船上的父母和子女)与sex一起推断为分类特征。pycaret将自动对分类特征进行一键编码，在这种情况下，将对Pclass、Sibsp、Parch和sex进行编码。然而，除了sex之外的这些特征都有其固有的等级顺序(普通性),对它们进行标记编码以捕捉等级顺序更为合适。

因此，在下面的输出中，在运行setup时，我们键入Quit并通过将下面列出的附加参数传递给setup函数来解决这个问题。注意，我们使用ignore_features参数从建模中排除了Name、Ticket和Cabin特征。

序数 _ 特征()

ordinal_features参数接受一个字典，其中包含我们想要标记为序数的特性以及相应的类别级别。

常化

我们还可以将normalize参数设置为 True(默认为 false)并使用normalize_method参数指定方法(比如说 minmax，默认为 zscore)。

宁滨数字特征

在这种情况下，让我们也对数字特征Age和Fares进行分类，因为它们值的每个单位变化不太可能影响我们对目标的预测。

设置

综上所述，我们修改后的setup命令如下所示:

当pycaret提示我们确认它的推断时，我们点击Enter。使用单个setup功能，pycaret已经完成:

• 缺失值插补
• 分类特征的编码(见下面的第 11 行和第 12 行，现在我们有 31 个特征)
• 规范化数字特征
• 收集数字特征
• 将数据集分为训练和测试

产出续下

产出续下

厉害！对吗？这还不是全部。通过向setup传递额外的参数，您可以根据数据集的性质探索许多设置定制。在上面的输出中，无论您在哪里看到 False/None，我们都保持默认的设置配置不变。您可以根据数据集的性质和想要进行的转换来启用相关的选项。配置这些附加参数的语法很简单。你可以在这里找到分类的完整列表和文档，在这里找到回归的。

因此，用一个函数就完成了预处理！

使用 compare_model()进行模型比较

compare_model()函数使我们能够生成 15 个分类模型，并一次性比较多个分类指标的交叉验证性能！模型还会根据精度指标(默认)自动排序，可以使用sort参数将精度指标更改为您喜欢的其他指标。

模型比较

此外，您可以单击任何其他首选指标列标题，根据新指标对模型排名进行交互式排序。这非常有助于快速确定我们首选指标的最佳性能模型，以便我们可以继续对所选模型进行超参数调整。

compare_models()返回指定度量上的最佳执行模型，该模型可以存储为模型对象。

通过将附加参数传递给以下对象，可以进一步定制compare_models()的输出:

• 省略型号(blacklist)
• 仅考虑特定型号(whitelist)
• 指定交叉验证折叠(fold，默认为 10)
• 四舍五入小数(round，默认为 4)
• 指定模型排序的度量标准(sort，默认为精确度)
• 要返回的顶级模型的数量(n_select，默认为 1。如果大于 1，则返回顶级模型的列表)

最后，我们看到，使用pycaret，只需两行代码，我们就能够在 ML 项目中快速完成预处理和初步模型比较阶段。

在接下来的文章中，我们将了解pycaret的强大功能，它可以使用遵循****_model()语法的函数来创建模型、调整模型超参数、集合模型、堆栈模型等等。

如果您喜欢这篇文章，您可能还想看看下面关于探索性数据分析(EDA)和用最少的代码和最大的输出进行可视化的文章。

[## 使用 SmartEDA 开发的 EDA

探索性数据分析——更智能、更快速的方式…

towardsdatascience.com](/eda-in-r-with-smarteda-eae12f2c6094) [## 拼凑——下一代 ggplots

进一步扩展 ggplot2 的多功能性…

towardsdatascience.com](/patchwork-the-next-generation-of-ggplots-1fcad5d2ba8a)

感谢您的宝贵时间！注意安全！！

PyCharm:在 Ubuntu 中管理多个 Python 版本

原文：https://towardsdatascience.com/pycharm-managing-multiple-python-versions-in-ubuntu-11224789845e?source=collection_archive---------22-----------------------

使用 PyCharm 运行不同的 Python 版本

来源:来自 Pixabay 的 MustangJoe 的照片

运行 Python 项目时的一个潜在问题是版本冲突。例如，您可能安装了 Python 3.6 作为 Python 的系统级版本，但是需要构建一个需要 Python 3.8 的特定程序。

然而，将 Python 3.8 安装为系统范围的版本会有使您的 Python 3.6 代码与新版本不兼容的风险。

在这种情况下，需要使用虚拟环境，即一个独立的环境，它运行一个不同于系统范围内运行的 Python 版本。

PyCharm 作为多个 Python 版本的 IDE(集成开发环境)非常有效。使用 Ubuntu 18.04，让我们看看如何通过 PyCharm 在虚拟环境中设置 Python 3.8。

安装 PyCharm

从 JetBrains 网站上，你会看到 PyCharm(在撰写本文时)既提供了付费的专业版py charm，也提供了名为社区的免费开源版。

通过单击下载，该文件将被下载为 tar.gz 文件。在撰写本文时，文件名是 pycharm-community-2020.2.3。

然后，可以将这个 tar.gz 文件解压缩到您选择的文件夹中:

资料来源:Ubuntu 18.04

如另一个教程中所述，导航到 bin 文件夹并打开一个终端。PyCharm 可以通过以下方式安装:

sh pycharm.sh

然后会打开一个相关的窗口，在这里可以定制 PyCharm 的默认设置。

在 PyCharm 中创建项目

如前所述，Python 3.6 是本例中 Python 的默认版本。在 PyCharm IDE 中，一旦选择了创建项目，将显示以下内容:

来源:PyCharm IDE

为了简单起见，这个项目被称为 pythonProject ，但是您可以为这个项目指定任何名称。从上面可以看出，python3.6 被选为基础解释器。在这种情况下，Python 3.6 和 Python 2.7 是系统上安装的仅有的两个版本。

来源:PyCharm IDE

一旦 Python 控制台启动，我们可以看到 Python 3.6.9 被加载。

也就是说，由于我们在自己的独立环境中使用这个版本的 Python，请注意相关的库必须单独安装到这个环境中。

例如，尝试运行 numpy 库第一次失败，因此需要在 PyCharm 项目中手动安装它，如下所示。

来源:PyCharm IDE

然而，当在标准的 Linux 终端(使用系统范围的版本)中运行时，库运行起来没有问题。

来源:PyCharm IDE

配置新的 Python 版本

如前所述，在许多情况下，用户需要使用不同于默认系统范围版本的 Python 版本。

对于这个例子，让我们假设我们希望安装 Python 3.8.0 。

正如来自的 Serhat Teker 的这个教程中所描述的，这个 Python 版本可以如下下载:

ssudo wget https://www.python.org/ftp/python/3.8.0/Python-3.8.0.tgz

然后，tar 文件可以按如下方式解压缩:

sudo tar xzf Python-3.8.0.tgz

在这种情况下，我们不希望 Python 3.8 成为默认版本。假设系统上的许多其他 Python 代码以不同的版本命名(在这种情况下是 Python 3.6.9 ),那么直接升级 Python 可能会导致严重的依赖性问题，并导致库被弃用。

因此，必须使用 sudo make altinstall 而不是标准的sudo make install——后者将覆盖默认的 Python 版本。在这方面，必须导航到新下载的 Python 文件夹，配置并启用优化，然后使用 altinstall。

cd Python-3.8.0
sudo ./configure --enable-optimizations
sudo make altinstall

在终端中检查 Python 版本时，您会注意到 Python 3.6.9 的默认版本保持不变:

来源:Ubuntu 终端

但是，Python 3.8 也可以作为一个单独的版本来调用:

来源:Ubuntu 终端

我们可以通过终端中的 which 命令找到每个 Python 版本的安装位置:

来源:Ubuntu 终端

相应地，可以通过 PyCharm 设置新目录，基本解释器现在设置为 Python 3.8，如下所示:

来源:PyCharm IDE

相应地，当打开新项目的 Python 控制台时，我们可以看到 Python 3.8.0 现在已经加载。

来源:PyCharm IDE

如上例所示，需要为此环境单独安装库，并且需要配置依赖关系。

然而，PyCharm 符合我们试图实现的目标——在虚拟环境中安装一个单独的 Python 版本，而不影响系统范围的版本。

结论

在本例中，您看到了:

• 如何在 Ubuntu 中安装多个版本的 Python
• 创建不干扰系统范围 Python 版本的虚拟环境
• 配置 PyCharm IDE 以使用不同的 Python 版本

在这方面，PyCharm 为使用不同的 Python 版本提供了一个非常有用的平台。安装的虚拟环境在技术上独立于 PyCharm 本身，可以简单地在终端上运行。然而，我们不能否认 PyCharm 确实使这个过程更加直观。

非常感谢您的阅读！非常感谢任何问题或反馈。你可以在michael-grogan.com找到更多我的数据科学内容。

免责声明:本文是在“原样”的基础上编写的，没有担保。它只是为了教育目的而写的，不应该被解释为任何形式的专业建议。我们不保证上述代码的执行适合您自己的系统，使用风险由您自行承担。

pydantic

原文：https://towardsdatascience.com/pydantic-688e897cfd3a?source=collection_archive---------5-----------------------

如果你篡改数据，你必须知道这个库

作者图片

拥有复杂的嵌套数据结构很难。在 Python 中存储这种数据的传统方法是嵌套字典。虽然 Python 字典很神奇，但通常会出现两个问题:(1)作为开发人员，我如何知道传递的字典中应该有哪种数据;( 2)我如何防止输入错误？

如果您使用类，您的编辑器或 mypy 可以防止属性名中的拼写错误。对于字典，任何有效的字符串都可以是一个键。

解决这两个问题的方法是使用一个库: pydantic 。它是一个验证和解析库，将数据映射到 Python 类。

先决条件

安装 pydantic via

pip install pydantic

对于本文，我假设您的数据是[people.json](https://gist.github.com/MartinThoma/517d20998501afc4fff72be032782d41)中的一个人的网络。他们有一个 ID，一个名字，一个由他们的 ID 给出的朋友列表，一个生日，以及他们银行账户上的金额。

我们想通过这个例子来计算每个人比他们的中位朋友多拥有多少财富。

我们想要解决的问题是:

• 丑陋无:把None作为朋友的值是相当丑陋的。在某些情况下，None不同于空列表。在这种情况下，我们就假设它的意思是一样的。用其他东西替换缺失数据被称为插补，有几种技术。
• 丑 Any :类型注释有Dict[str, Any]，因为它被认为太复杂甚至不可能知道代表一个人的字典也能映射哪些值。
• 难看的索引:只是语法，但是["id"]比.id长 50%。这就是帮 / 蒙克存在的原因。
• 错别字:如果你在字典的任何一个字符串索引中出现了错别字，没有静态代码分析工具可以帮你检测出来。没有自动完成将正常工作。
• 后期错误:好的软件将职责模块化。一个模块负责业务逻辑，一个模块负责输入/输出。当我写“模块”时，我指的是一个代码单元——它可以是一个 Python 模块、一个类、一个函数、一个方法。简单地用一个函数来处理这些不同类型的复杂性是不好的。获得正确的业务逻辑本质上不同于确保您的输入/输出看起来像预期的那样。一个是关于定义和验证适当的接口，另一个是关于理解领域。
• 文档:新开发人员必须定期阅读更大应用程序的代码。通常，没有时间让对应用程序了如指掌的人详细解释每一个部分。大多数时候，甚至没有一个人什么都知道——即使项目是由一个开发人员完成的。我总是忘记部分代码。因此记录是至关重要的。在 Python 中，记录预期的参数类型及其返回值至关重要。有Dict[str, Any]总比没有好，但比Person差得多。

创建一个 pydantic 模型

我们为一个人的 ID 创建了一个新的类型，仅仅是因为 PersonId 比 T2 更有意义。然后我们子类化[pydantic.BaseModel](https://pydantic-docs.helpmanual.io/usage/models/#basic-model-usage)

将其用于输入解析

接下来，我们使用[parse_file_as](https://pydantic-docs.helpmanual.io/usage/models/#parsing-data-into-a-specified-type)来读取 JSON 文件:

请注意， datetime 和 Decimal 是自动解析的——您仍然应该经常查看它是如何完成的。尽早进行数据验证是有好处的，这样错误也会在已知的地方尽早发生。这意味着 pydantic 轻推你到一个好的设计。我喜欢它😍

约束类型

约束类型是某个取值范围内的整数/浮点数或与正则表达式匹配的字符串😃

缺少数据:使用默认值

如果您的 JSON 可能缺少一些您想要的属性，那么您需要使用默认值。典型的默认值是None，这意味着您需要将类型更改为Optional[what it was before]。这通常很难看，因为您需要在代码的后面检查None。

特别是对于列表，您可能要考虑使用空列表。你这样做:

对于不可变的数据类型，比如字符串、整数、浮点数、元组，可以简单地赋值。对于可变列表，您需要使用带有default_factory的 Field，每次都会生成一个新列表。如果你还不知道的话，学习一下为什么可变缺省是邪恶的。

附加数据:忽略、允许或禁止

有时在开发时不可能知道 JSON 对象具有哪些属性。尽管如此，你还是需要把它们传来传去。这是非常不幸的，应该受到挑战，但这是可能发生的。

Pydantic 称那些为临时演员。如果忽略它们，读取的 pydantic 模型将不会知道它们。被忽略的额外参数被丢弃。允许他们意味着接受这个不幸的设计是必要的。允许的额外内容将是解析对象的一部分。最后，禁止额外参数意味着如果出现额外参数，将抛出pydantic.ValidationError 异常。

这是通过向 pydantic 模型添加一个名为Config的子类来配置的:

重命名属性

名字很重要。可读性计数。在 Python 中，变量应该遵循一个snake_case命名方案，而在 JavaScript 中，变量应该遵循一个camelCase命名方案。为了实现这两个目标，pydantic 提供了[allow_population_by_field_name](https://pydantic-docs.helpmanual.io/usage/model_config/)作为配置参数。

验证器

有时候，简单的类型是不够的。你想检查更复杂的东西。

文件已经给出了一个很好的例子:

你可以检查几乎任何东西，只要你只需要类本身。请不要对数据库运行查询来进行一致性检查；例如，如果您希望用户名是唯一的或类似的。尽管您可能会运行此操作，但通过创建“数据容器”来触发数据库查询将是意想不到的。

在我们的例子中，我们可能想要阻止人们成为自己的朋友:

除了抛出异常，我们还可以简单地修复它:

使用 Pydantic 进行基于属性的测试

在 Unsplash 上科学高清摄影

基于属性的测试自动生成被测函数的输入，并确保满足某个属性。在最简单的情况下，这个属性是测试中的函数不会崩溃。如果你想了解更多关于这种类型的测试，请阅读我的关于基于属性的测试和假设的文章。

顺便说一下，这个测试实际上指出了一个潜在的问题:

Falsifying example: test_two_people(
    person_a=Person(id=0, name='', bank_account=Decimal('NaN'), birthdate=datetime.date(2000, 1, 1), friends=[]),
    person_b=Person(id=0, name='', bank_account=Decimal('NaN'), birthdate=datetime.date(2000, 1, 1), friends=[]),
)

更整洁的东西

Pydantic 非常棒:

• 您可以从它的模型中生成一个模式( source )，
• mypy 插件提供了更好的类型检查
• 对字典的序列化可以用.dict()方法完成，对 JSON 字符串的序列化可以用.json()方法完成。

操作安全

照片由本工程图在 Unsplash 上拍摄

我通常担心的一部分是普遍支持。以下是 pydantic 拥有健康社区的一些指标:

• ✔️ GitHub : 4.5k 星，404 叉，172 贡献者
• ✔️用法: FastAPI 使用它。微软的 Onefuzz 使用它。 AWS Lambda Powertools 也是如此。许多机器学习项目
• ✘ PyPI 项目:很遗憾，这个只有一个维护者。如果有两个，我会更喜欢一点。只是为了防止有人丢失密码。
• ✔️Self-decided 项目状态:维护者认为 pydantic 可以投入生产/稳定。

摘要

pydantic 是一个非常棒的数据解析和验证库。它可以很好地支持你在 Python 中获得更好的类型注释。用它！

非线性模型参数不确定性量化的 PyMC3 和贝叶斯推理(二)

原文：https://towardsdatascience.com/pymc3-and-bayesian-inference-for-parameter-uncertainty-quantification-towards-non-linear-models-a03c3303e6fa?source=collection_archive---------19-----------------------

参数和模型不确定性量化

导航的不确定性(作者图片)

在本系列的第 2 部分中，我们将深入 PyMC3，并学习如何将其应用于非线性模型的不确定性量化(UQ)。为了验证，我们还将看到贝叶斯方法如何与 UQ 参数的频率主义方法进行比较。我假设你已经为你的模型优化了参数。我将使用的模型来自本系列的第 1 部分。

我们将在本系列中使用的库如下:

将要使用的库

对于那些从第 1 部分跳到这里的人，或者那些需要重温第 1 部分的人，非线性模型是集总电容对象(多么拗口)的一阶瞬态热响应的解决方案。

集总电容物体一阶瞬态热响应的解

Python 中的函数

我们假设 T_0 = 100C，T_env = 30C。我们试图估计冷却参数τ和模型的不确定性。这个问题的日常例子是冷却保温良好的热饮！

热饮…也许我们把它放在一个非常好的保温瓶里(图片由作者提供)

PyMC3 下的贝叶斯推理

为了开始对模型参数的 UQ 进行贝叶斯推断，我们首先使用 PyMC3 创建一个模型，这里我在第 1 行将其称为“model_g”。

引用 Therese Donovan 和 Ruth M. Micket 的话[1]，

“如果不使用先验分布，就无法进行贝叶斯分析”。

贝叶斯框架。p 代表概率。

这里我们的先验是 P(τ，ϵ)，我们将把这个先验分解成两个先验，每个参数一个。而我们对τ的分布知之甚少或没有多少先验知识。然而，我们可以提供一个弱信息先验[2]。通过简单地画出τ的不同值，我们可以看出τ一定是有界的。我们将假设τ为均匀分布，下界为 0，上界为 100。我们可以看到，这很容易涵盖这个问题的可能值，然后一些。从第 1 部分中，我们发现τ的最佳参数是 14.93，从图中可以看出这是有意义的。

调查 tau 的前科

我们需要的下一个先验是误差的标准差，ϵ.我们知道这个参数不能小于 0，因为标准差总是大于或等于 0。半正态分布似乎适合于ϵ.的先验分布我们再次通过选择标度值为 5 的半正态分布来选择弱信息先验。我使用ϵ而不是第一部分中使用的σ，以避免与也用于衡量半正态分布的σ术语混淆

调查ϵ的前科

将所有内容结合在一起，第 5 行描述了预测温度“T_pred”来自使用我们的非线性牛顿冷却函数的平均值(μ)的正态分布。这个模型是我们的可能性，代表 p(数据| τ，ϵ)。我们还让 PyMC3 知道我们模型的标准偏差是ϵ.的先验分布最后，我们传递的最后一个参数是观察到的“实验”温度数据。

使用第 8 行，我们让 PyMC3 通过抽取 2000 个后验样本，使用 NUTS 采样来完成繁重的工作。我非常非常推荐去看看 Richard 的 McElreath 关于采样算法的网页。我们通常无法显式求解上述贝叶斯框架方程中的分子，因此我们必须使用采样算法来帮助我们探索它。

一旦执行了采样，我们就可以检索摘要

az.summary(trace_g)

跟踪摘要

在第 1 部分中，我们发现τ_opt = 14.93，σ的估计值=1.42。查看上表，我们的采样方法的平均值也提供了非常相似的结果，这是有道理的。我们之前做的优化也允许我们检查我们的贝叶斯推理。

我们参数的联合图

joint_plt = az.plot_joint(data_spp, var_names=['τ', 'ϵ'], kind='kde', fill_last=False);

使用上面的代码，我们可以通过对后验概率 P(τ，ϵ |数据)进行抽样来生成一个联合图。

我们还可以获得一个跟踪图，以验证我们的采样没有停滞，并且我们的采样充分探索了参数空间。我们想要像“毛毛虫”一样的痕迹。

trace_fig = az.plot_trace(trace_g,
 var_names=[ 'τ', 'ϵ'],
 figsize=(12, 8));

追踪图

我们还可以通过使用以下方法获得参数的最高后验密度(HPD ),这里我们选择 99%。

posterior_params = pm.plot_posterior(data_spp,
                                     point_estimate='mean',
                                     credible_interval=0.99,
                                     textsize=20,
                                     round_to=4);

如果没有噪声，我们的平均值将完全等于“未知”参数。由于存在噪声，就像真实数据一样，我们的参数值存在不确定性，我们需要量化其中的不确定性。

从第 1 部分可以看出，实际值τ=15，σ=1.5。不仅平均后验概率接近这些值，而且 HPD 很容易包含真实值，这正是我们想要的。

贝叶斯推理的好处之一是，除了计算参数的平均值，我们还可以得到参数的分布，如上所示。我们的非线性模型的 UQ 就像从我们的分布中取样一样简单，这是在上面的 sample_posterior.py 中完成的，为了清楚起见，下面使用方法

pm.sample_posterior_predictive()

计算模型中的不确定性，不仅从参数中，而且从数据中，是通过从我们的后验样本完成的。使用 NumPy 中的百分位数方法，从我们的后验预测中生成分位数是非常容易的。

将所有这些放在一起，让我们不仅绘制平均后验反应，而且绘制 95%后验预测。

参数和模型 UQ 的贝叶斯推断

频繁主义方法

如果我们想验证我们的贝叶斯推断，我们还应该将我们的结果与频率主义方法进行比较。由于模型是非线性的，这将需要更多的手工计算，但是对于我们将要使用的库来说，这应该不会太难。

首先，我们需要为模型“grad_newton_values”构建一个灵敏度矩阵。由于模型中只有 1 个参数，这个矩阵只是一个列向量。我们将使用自动签名库来自动区分关于τ的模型。

模型误差δk 可以从我们的灵敏度阵列中计算出来。利用δk、t-逆数值(1.99)和我们对第 1 部分中σ的估计(1.42)，我们可以计算参数τ [3]的置信区间。

我们发现，在 95%的置信区间下，τ在 14.571 和 15.292 之间。非常接近贝叶斯方法！

接下来，我们想比较贝叶斯和频率主义方法之间的 UQ 传播在我们的模型。这包括模型参数和数据的不确定性。下面我为我们的非线性模型生成了置信区间和预测区间。

我用星形影线绘制了 95%的贝叶斯后验预测，用纯紫色绘制了 95%的频率预测区间。这两种方法几乎完全重叠，因为在紫色阴影区域下面可以看到星形影线。我们看到在边界上，最小的绿色。看起来我们的贝叶斯方法和频率主义方法是一致的。

根据这些预测，我们可以提供作为时间函数的饮料温度的不确定性。例如，在 25 分钟时，如果我们多次运行该实验，我们预计大多数热饮的温度在 40℃到 45℃之间，有 95%的把握。

两种方法的预测区间几乎相同

结论

我们对参数的估计以及模型本身进行了 UQ。我们已经发现 PyMC3 库对于处理非线性模型是多么强大，以及它与 Python 的集成是多么好。贝叶斯方法为我们如何开发和设计模型提供了极大的自由度。如果需要，我们可以将我们的可能性或温度模型从正态分布改为另一个参数模型。PyMC3 的过程不会有太大变化。

参考文献

[1] T. Donovan 和 R. M. Mickey，贝叶斯统计初学者:循序渐进的方法。牛津大学出版社。

[2] O. Martin，用 Python 进行贝叶斯分析:使用 PyMC3 和 ArviZ 进行统计建模和概率编程介绍，第 2 版。帕克特出版有限公司，2018 年。

[3] R. C. Smith，不确定性量化:理论、实施和应用。暹罗，2013 年。

非线性模型的 PyMC3 和贝叶斯推理:第一部分

原文：https://towardsdatascience.com/pymc3-and-bayesian-inference-for-parameter-uncertainty-quantification-towards-non-linear-models-d1d50250f16a?source=collection_archive---------24-----------------------

参数和模型不确定性量化

不确定的世界(作者图片)

像 Stan 这样的概率编程语言已经为科学家和工程师提供了进入贝叶斯统计世界的便利手段。作为一个熟悉 Python 编程的人，我一直在努力寻找完美的贝叶斯库来补充 Python 美妙而简洁的语法。

学习新材料和新图书馆的最大挑战之一是找到足够的资源和例子。我找不到一篇关于使用 Python 和贝叶斯推理对非线性模型进行参数和模型不确定性量化(UQ)的文章。因此，我希望这篇文章能帮助其他人将这些方法应用到他们自己的研究或工作中。

这将是一篇由两部分组成的文章，目的是让更多的读者能够阅读它，同时也让浏览和阅读文章变得更加容易。

第 1 部分将快速讨论贝叶斯推理的两个公共库:PyStan 和 PyMC3。我们还将建立一个非线性函数用于贝叶斯推理。这将包括通过 SciPy 库进行参数优化。我们将使用最佳参数作为我们贝叶斯推理的起点，并检查我们的贝叶斯方法。如果你已经熟悉了参数优化，可以直接跳到第 2 部分中的 PyMC3 部分。

第 2 部分是我们开始使用 PyMC3 并深入研究语法和结果的地方。为了验证，我们还将看到贝叶斯方法如何与频率主义方法进行参数不确定性量化的比较。

我不打算对 Baye 定理和整个贝叶斯框架进行过多的推测，有大量关于 Medium 的文章可以满足您的好奇心。本文的主要成果是让读者了解如何利用 PyMC3 来量化线性模型(尤其是非线性模型)中使用的参数和模型的不确定性。Jupyter 笔记本可以在这里找到。

PyStan

Stan 以 PyStan 的形式提供给 Python，但是语法相当令人不快。让我们看一个在 PyStan 中创建模型的例子，这个模型被分配给一个 string 类型的变量。然后这个模型会在 PyStan 的帮助下编译成 C++。

上面是一个简单的二次模型示例，它有两个参数，α和β，响应为 y。当开始使用更复杂的模型时，该模型非常冗长，很难调试。肯定有更 Pythonic 化的东西…

PyMC3

我相信 PyMC3 对于刚进入 Python 概率编程领域的人来说是一个完美的库。PyMC3 使用原生 Python 语法，使得调试更容易，也更直观。让我们深入一个例子，看看这个库的威力。

我们将使用的库

在这个例子中，我将使用一个简单的非线性方程的日常例子，希望这些概念更容易理解和记忆。

让我们从利用牛顿冷却定律来冷却一个热物体开始。我们假设物体中不存在热梯度(集总电容)。

让我们假设热物体是一种饮料(图片由作者提供)

对于那些对常微分方程(ODE)生疏的人来说，不用担心。解决方法如下。

T(t)表示在给定时间 t 饮料的温度。τ是我们感兴趣的参数，它决定了这个问题的冷却速率。T_env 是我们的环境温度，T_0 是 t=0 时饮料的温度。我们假设 T_env 不随时间变化。

在 Python 中:

让我们也假设我们提前知道，先验，环境温度是 30 摄氏度，物体开始于 100 摄氏度。最后，我们假设未知参数τ=15。

我们将首先绘制出热饮冷却的理论解，如蓝色实线所示。这里我们可以对方程和参数有一个直观的认识。

我们还将生成一些带有人工噪声的数据来模拟现实生活中的数据，这里显示为红色散射点。我们将假设数据是 IID 。饮料温度的响应是高斯型的，即 T(t)~ N(μ，σ)。平均值μ将遵循上面的解。我们将很快解决从数据本身不知道τ的问题。响应中的标准偏差为σ=1.5。这就是我们如何产生假噪声的，我们也能够从数据中估计出来。

根据这些模拟数据，我们想要估计冷却参数τ，并估计其值的不确定性。理想情况下，我们的估计参数和不确定性将封装真实值。同样，我们也希望模型标准误差中的不确定性包含σ。

为了阐明术语“封装”的含义，假设我们估计一个参数γ= 100 ^ 10。我们希望真实参数在 90–110 的范围内，与 100±10 相同。

优化:估计τ

在估算τ的不确定性之前，我们必须首先优化最符合数据的τ值。用更专业的术语来说，我们希望找到τ，使得残差平方和最小。

n 是数据点的数量

其中括号左侧的函数是我们的数据，括号右侧的函数是我们的 ODE 的解。

我们发现τ_opt = 14.93，非常接近τ=15 的真实参数的估计值。回到我们最初的问题，记住响应，饮料的温度，是高斯的，即 T(t)~ N(μ，σ)。我们终于有了τ的估计，我们可以检验正态性。

请记住，如果我们的响应来自高斯分布，那么我们的残差(数据减去预测)也必须是高斯分布。

residuals = T_data -  newton_cooling(tau_opt, T_0, T_env, time_data)

想一想，线性变换，或者更基本的，Z 值的变化。

残差，o =数据和我们的模型之间的差异，看起来很大，平均值大约为 0。它们看起来是随机的，同质的，大部分数据在两条红色虚线之间。

让我们看看我们是如何估算σ的

sigma_error_square = (1/(len(time_data)-1))*ss_min
sigma_error = np.sqrt(sigma_error_square)
print(sigma_error)

我们发现σ的估计值是 1.42，非常接近 1.5

剩余剧情。两条红色虚线代表 0 1.96*σ或 95%的数据。

不确定性量化

如果你还和我在一起，太好了！现在我们已经有了τ和σ的估计值，我们终于可以开始研究贝叶斯统计和 PyMC3 了。点击查看该系列的第二部分。

PyQGIS 连接表-厄瓜多尔自然保护区的国际游客

原文：https://towardsdatascience.com/pyqgis-join-table-international-visitors-in-natural-protected-areas-in-ecuador-7deb7ed8d868?source=collection_archive---------41-----------------------

如何在 PyQGIS 3.10 中使用连接表创建 choropleth 地图

PyQGIS 为管理 GIS 过程中的自动化提供了一个很好的解决方案。厄瓜多尔国家保护区系统(SNAP)以电子表格形式获取国际游客的信息，但不一定是地图。因此，如果需要地图，可以使用一个 python/PyQGIS 脚本来创建地图，该脚本可以输出适当的图层，包括受保护区域的游客数量。此处可在厄瓜多环境信息地理门户下载快照图层。

在这篇文章/教程中，我将向你展示如何在 Qgis 3.10 中将连接表自动化为 Shapefile。这个脚本最好的一点是，它通过 QInputDialog() 与用户(在本例中是您)进行交互，因此您必须使用的输入只是您的数据的 URL。

所以，我们来看看这个剧本的产物。如上所述，该脚本生成了一个新图层，其中包含添加到 自然保护区 的 shapefile 中的 游客数量 的附加信息。如果你不熟悉加入，看看这里的。

图片由作者提供，choropleth 地图包含国际游客数量的相关信息

首先，您必须在 Qgis 中打开 python 控制台并打开编辑器。

图片由作者提供，Python 控制台，编辑器窗口“无标题-0”位于右侧

我们的第一个任务是打开 shapefile，并在 Qgis 中打开表。在本教程中，您唯一需要做的就是复制和编辑数据的 URL。我的例子的 URL 如下:C:\ Desktop \ Pyqgis \ data \ natural _ areas . shp

首先，将下面的代码复制并粘贴到 python 编辑器中，然后运行。然后，您将看到一个新的对话窗口打开。

#Dialog 
qid1 = QInputDialog()
title = ‘Write the URL’
label = ‘Shapefile url: ‘
mode = QLineEdit.Normal
default = ‘<shp url here>’
snapUrl, ok = QInputDialog.getText(qid1, title, label,mode, default)#Adding Layer 
snap_layer = iface.addVectorLayer(snapUrl,’Natural Areas’,’ogr’)

图片作者，输入框添加您的 shapefile 的 URL

然后，您将看到 shapefile 如何像图层一样添加到内容列表中。添加 shapefile 时，您只需注意正确书写您计算机的 URL，包括扩展名 shp 。
( 例如\…\…\areas.shp )

作者的图片，屏幕上的 shapefile。

如您所见，您正在与 python 控制台进行交互。基本上就是在问你开什么，你在回复一个网址。下一步是执行相同的过程，但是包括显示的表格，在本例中是 。csv ，要与***Natural areas . shp .***进行连接记住，要正确处理连接表，两边需要有相似的属性。在这种情况下，shapefile 和表都有一个包含自然区域代码的列。
按照下面的代码打开表格。

#Dialog CSV june
qid2 = QInputDialog()
title = ‘url of table’
label = ‘csv url: ‘
mode = QLineEdit.Normal
default = ‘<csv url here>’
juncsvUrl, ok = QInputDialog.getText(qid2, title, label,mode, default)#Adding csv juneurljun = r’file:\\’+juncsvUrl+’?delimiter={}’.format(“,”)
juncsv = QgsVectorLayer(urljun, “June 2019”, “delimitedtext”)
iface.addVectorLayer(juncsvUrl, “June visitors”, “ogr”)

图片作者，输入收件箱添加您的 CSV 表的网址

然后，您将看到该表被添加到屏幕上。

图片由作者提供，屏幕上的 CSV 表。

目前，我们已经获得了连接表所需的数据。现在我们运行下面的代码来处理连接。你会看到一个新的消息框打开，要求一个网址，你必须写下你要保存新层的文件夹的网址。不要忘记包括扩展名。例如:①…\…\join.shp )。

如果你用的是自己的数据。请注意，在行字段中必须包含 shapefile 的列名，在字段 2 中必须包含 CSV 表的列名。

#Dialog JOIN 
qid4 = QInputDialog()
title = ‘Output layer of Join’
label = ‘url + \name.shp: ‘
mode = QLineEdit.Normal
default = ‘<new layer url here>’
out, ok = QInputDialog.getText(qid4, title, label,mode, default)#Join Table June
processing.run(“native:joinattributestable”, \
{‘INPUT’:snapUrl,\
‘FIELD’:’csnap’,\
‘INPUT_2’:juncsvUrl,\
‘FIELD_2’:’codt_snap’,\
‘FIELDS_TO_COPY’:[],\
‘METHOD’:1,\
‘DISCARD_NONMATCHING’:True,\
‘PREFIX’:’’,\
‘OUTPUT’: out,})jun_visitors = iface.addVectorLayer(out,’Visitors in June 2019',’ogr’)

作者的图像，新图层 URL 的消息框。

最后，您会看到新层被添加。包含访问者数量的那个。

图像由作者创建，屏幕上的新图层包含表中的数据。

这是您应该知道的关于将表连接到 shapefiles 的全部内容。您可以轻松地重写代码，并使其与您自己的数据一起工作。本教程是我关于 PyQGIS 的一个简短课程的一部分。如果您对更大的项目有任何疑问，请访问我的个人资料。

Pyro 自上而下预测|应用案例

原文：https://towardsdatascience.com/pyro-top-down-forecasting-application-case-4781eb2c8485?source=collection_archive---------48-----------------------

将一段时间内的点连接起来，满怀信心地进行预测(-区间)。

跨越时间的连接点— 照片由伊斯雷尔·帕拉西奥在 Unsplash 上拍摄

你有没有想过如何解释时间序列预测中的不确定性？你有没有想过应该有一种方法从以前看到的数据中生成数据点，并对确定性做出判断？我知道我有。
如果你想建立捕捉概率并保持信心的模型，我们建议使用概率编程框架，如 Pyro 。
在之前的一篇文章中，我们研究了 NGBoosting，并将其应用于 Kaggle 上的 M5 预测挑战。简单回顾一下——M5 预测挑战要求我们预测沃尔玛商品的销售将如何发展。它提供了来自不同州的不同商店的不同类别的商品的大约 4-5 年的数据，并要求我们预测我们没有信息的 28 天。作为对挑战和数据集的概述，我们仍然推荐这款令人惊叹的笔记本。上次我们总结说，我们用 Pyro 得到了更好的结果，这里是一个简单的演示我们是如何做到这一点的。

随着时间的推移，有不同的方法来建模和预测数据。有自顶向下模型、状态空间模型和层次模型——仅举几个例子。在这篇文章中，我们看到了如何通过一个非常粗糙和基本的自顶向下的模型进行预测。我们充分利用了 Pyro 中现有的预测器对象。我们只是告诉模型如何做出明智的预测。之后，我们将所有数据转储到其中。就是这样——很简单！没有花哨的先验，没有花哨的潜在分布假设，只有数据和概率框架。我们还需要做一些假设，我们会带你经历这些。

Sidenote: One more elegant way to do forecasting is with hierarchical models. These models allow accounting for different distributions over different categories in your data. You can assign individual priors and make use of information that you have about your model. The elegance of such an approach will be covered in another post.

最后，我们将所有的时间序列信息转储到模型中，这是一段时间内的销售，同时将销售作为独立的随机事件保存。我们不用别的。
Pyro 编程框架，特别是 ForecastingModel ，是我们模型的引擎，我们还需要一个驱动程序。
推动我们的模型朝着正确方向发展的因素是:

prediction = regressor + trend + seasonal + motion + bias

这一行简单的代码结合了我们随时间采样的 Pyro 对象。现在，让我们看看所有这些部分的含义:

构建模型

实际模型零件由回归器捕获。这一部分寻找合适的模型权重，并将它们与直接来自数据输入的特征进行比较。每个权重都位于经典的高斯曲线上。这意味着我们从以零为中心、标准差为 1 (μ=0，σ =1)的正态分布中采样。

weight = pyro.sample("weight", dist.Normal(0, 1).expand(
        [feature.size(-1)]).to_event(1))
 regressor = (weight * feature).sum(-1)

第二个参数是趋势。我们参考我们的概率窥镜，然后决定它应该来自对数正态分布，例如 l∞LN，μ=-2，σ =1。一个关键的因素是我们考虑了我们提供给模型的时间特征。

trend_coef = pyro.sample("trend", dist.LogNormal(-2, 1))
trend = trend_coef * time

现在没有季的趋势是什么？潮流来来去去，但季节是我们可以依靠的。我们知道基础数据描述了一段时间内的购买情况。我们认为销售额的变化是通过一周中的某一天来捕捉的。因此我们说

with pyro.plate("day_of_week", 7, dim=-1):
    seasonal = pyro.sample("seasonal", dist.Normal(0, 5))
seasonal = periodic_repeat(seasonal, duration, dim=-1)

这里我们使用了另一个简洁的烟火物体:盘子。这种结构用于模拟随机事件的独立性。这意味着我们正在观察一个实例的发生。这个观察来自我们数据的每个时间步。彻底的解释可以在烟火 SVI 教程中找到。出于简单的目的，我们将其分解，并说我们认为数据在七天的窗口内上升和下降。我们将其转化为独立采样的事件。最后但并非最不重要的是，不仅仅有一个七天的间隔，但是我们在我们所有给定的时间输入 aka 上周期性地重复这个。我们的持续时间是时间张量的大小。

最后，每一个专业的 ML 模型都需要一些可训练的参数，以说明有时事情并不像你期望的那样好——我们称之为偏差。我们决定我们的偏差(b)来自一个正态分布，这样 b∞N 其中μ=0 和σ =10 或 loc=0 和 scale=10，如果你想谈论 Pyro。

bias = pyro.sample("bias", dist.Normal(0, 10))

在我们继续之前，让我们看看参数的分布。

图 1——根据颜色绘制的所有采样参数的密度图。y 轴被截断，趋势对数正态分布的最大值为> 2.5。

你，细心的读者，可能会在这一点上得出结论。你有你需要的一切，我们同意。我们现在所拥有的足以给我们一个估计，并对未来几天做出预测。不过我们可以做得更好。到目前为止，我们预测的是一个非常僵硬的表达。模型被告知要符合输入，但我们不允许留有余地。为了便于讨论，让我们假设加州的所有人都决定在周四下午而不是周六去购物；或者危机发生，每个人都去为世界末日做准备。我们的模型可能已经正确预测了所有以前的星期六，但是那个星期四真的把它打乱了。为了解决这个问题，我们引入了一个叫做重新参数化的东西。有关该主题的详细概述，请参见[2]。
为了使它更加健壮，我们通过另一个重新参数化来进行重新参数化——我们这样做的原因在这里解释。
将所有东西打包在一起，我们得到一个漂移，它定义了我们模型的运动。这就翻译成了烟火:

drift_stability = pyro.sample("drift_stability", dist.Uniform(1, 2))
drift_scale = pyro.sample("drift_scale", dist.LogNormal(-17, 5))
with self.time_plate:
   with poutine.reparam(config={"drift": LocScaleReparam()}):
      with poutine.reparam(config={"drift": SymmetricStableReparam()}):
         drift = pyro.sample("drift", dist.Stable(drift_stability, 0, drift_scale))
motion = drift.cumsum(dim=-1)

现在我们有了模型运行所需的一切。

你的工作模型

对于实际的模型拟合部分，我们使用 Pyro 的随机变量推理机或 SVI。这个算法允许我们在合理的时间内有效地计算后验分布。在我们的例子中，这样做是为了最大化证据下限或 ELBO 。考虑 ELBO 的一个简单方法是观察两个不同的分布，一个当前分布和一个先前分布。如果分布相差很大，我们会得到一个低的 ELBO 值。另一方面，如果我已经拟合了一个好的分布，我的下一个可能也一样好，与已经观察到的一个很好地对齐，并且我得到一个高的 ELBO 值。有些人甚至会说两种分布的差异被最小化。要详细阅读它，我们可以参考 SVI 的原始论文(见【3】)。现在让我们假设这是可行的，因为聪明的研究人员已经充分实现了它。

We will go into the interplay of Stochastic Variational Inference and ELBO (or minimizing KL) in another article.

现在，我们已经建立了一个工作模型，选择了我们的算法来计算我们的后验分布，我们必须设置超参数，加载数据并开始训练。对于数据处理，我们广泛关注 github 上的 M5 初学者工具包实现，这是 Pyro PPL 团队慷慨提供的。这是一种简单的方法，可以对销售数据等进行汇总。请自行实现。这是了解张量的可靠方法。对于我们的初学者模型，我们只关心一段时间内的销售数据。
我们现在设置参数并实例化预测器，如下所示:

forecaster_opt = {
    "learning_rate": 0.1,
    "learning_rate_decay": 0.1,
    "clip_norm": 10,
    "num_steps": 3501,
    "log_every": 100,
}
forecaster = Forecaster(TopDownModel(), data, covariates[:-28], **forecaster_opt)
samples = forecaster(data, covariates, num_samples=1000).exp().squeeze(-1).cpu()
pred = samples.mean(0)

在上面的代码中，我们安装了预测器，它使用了一个 DCTAdam 优化器。我们提供了培训的协变量，并保留了 28 天，我们在下面的行中将其用于实际预测。实际的预测是样本的平均值，但是我们不要在这里妄下结论。当我们运行上述代码并正确实现我们的模型时，我们应该得到以下输出—实际训练:

预测者训练过程用指定的模型结束了 3501 步的所有训练数据。随着时间的推移，损耗被最小化。

我们看到我们的损失减少了。这很好，它告诉我们，随着时间的推移，我们预测的计算分布会变得更好。这里的损失是相对于数据和部分预测者文档的负 ELBO 的商(见【4】)。

真正不确定的事情——自信和一切

我们现在在记忆中保存了一个最终模型。我们可以用它来计算损失，用它来对抗新数据，戳戳它，看看它会产生什么。
例如，我们可以对一段时间内的销售额进行采样，如下图所示:

图 2 —真实数据和随时间推移从拟合模型中采样的数据。红线表示最近可用(真实)的 14 天，而模型的抽样预测为浅绿色。报告的预测值是所有样本的平均值(绿线)。我们可以预测所需的 28 天。每个时间步长 10%-90%的置信区间用红色阴影显示。请参考下面的代码。

为了计算各个步骤的置信度，我们要求我们的模型有 100 个样本，并计算这些样本的分位数(0.1 和 0.9)。这就是你在图 2 中看到的红色区域。既然这个模型由你支配，如果你想使它更精确，你可以要求更多的样品。为了更好地拟合模型，图中的销售额是累计的记录值。要获得更真实范围内的值，请使用 numpy 的 exp 转换这些值。

总结

我们已经展示了如何构建 Pyro 模型并使用预测器实现。您可以使用您所知道的关于问题的事实来设置模型参数并找到合适的分布。Pyro 附带了很多功能。其中之一是 SVI，一种计算后验概率的算法，这使得模型训练首先成为可能。在我们拟合我们的模型后，我们可以从中取样，预测数据点并计算我们预测的置信区间。总而言之，概率编程和 Pyro 之类的框架提供了一个强大的工具集，在处理不确定性和对现实世界数据建模时可以使用。

完整的型号代码

class TopDownModel(ForecastingModel):
  """
  Top-Down Hierarchical Forecasting Model
  """
  def model(self, zero_data, covariates):
    # check univariate data
    assert zero_data.size(-1) == 1 
    duration = zero_data.size(-2) time, feature = covariates[..., 0], covariates[..., 1:] bias = pyro.sample("bias", dist.Normal(0, 10))
    trend_coef = pyro.sample("trend", dist.LogNormal(-2, 1))
    trend = trend_coef * time weight = pyro.sample("weight", dist.Normal(0, 1).expand(
        [feature.size(-1)]).to_event(1))
    regressor = (weight * feature).sum(-1) # weekly seasonality as independent events
    with pyro.plate("day_of_week", 7, dim=-1):
      seasonal = pyro.sample("seasonal", dist.Normal(0, 5))
    seasonal = periodic_repeat(seasonal, duration, dim=-1) drift_stability = pyro.sample("drift_stability", dist.Uniform(1, 2))
    drift_scale = pyro.sample("drift_scale", dist.LogNormal(-17, 5)) # introduce drift
    with self.time_plate:
          # We combine two different reparameterizers: the inner SymmetricStableReparam
          # is needed for the Stable site, and the outer LocScaleReparam improves inference.
          with poutine.reparam(config={"drift": LocScaleReparam()}):
              with poutine.reparam(config={"drift": SymmetricStableReparam()}):
                  drift = pyro.sample("drift",
                                      dist.Stable(drift_stability, 0, drift_scale))
    motion = drift.cumsum(dim=-1) # predict
    prediction = regressor + trend + seasonal + motion + bias
    # Pyro Forecast is multivariate - univariate timeseries is needed
    prediction = prediction.unsqueeze(-1) # heavy tail nose to account for outliers
    stability = pyro.sample("noise_stability", dist.Uniform(1, 2).expand([1]).to_event(1))
    skew = pyro.sample("noise_skew", dist.Uniform(-1, 1).expand([1]).to_event(1))
    scale = pyro.sample("noise_scale", dist.LogNormal(-5, 5).expand([1]).to_event(1))
    noise_dist = dist.Stable(stability, skew, scale)
    with poutine.reparam(config={"residual": StableReparam()}):
      self.predict(noise_dist, prediction)

抽样和绘制置信区间

# sample from the fitted model to forecast data
samples = forecaster(data[T0:T1], covariates[T0:T2], num_samples=100)
samples.clamp_(min=0)p10, p50, p90 = quantile(samples[:, :], (0.1, 0.5, 0.9)).squeeze(-1)fig, ax = plt.subplots()ax.set_title("aggregated (log1p) sales over 42 days")
ax.plot(data[T1:T2].squeeze(-1).cpu().numpy(), "r-",lw=4, label="truth")
ax.fill_between(np.arange(0, (T2-T1)), p10.cpu().numpy(), p90.cpu().numpy(), alpha=0.3, color="red", label="confidence")
ax.plot(samples.squeeze(-1).T.cpu().numpy()[:, 10:13], color="green", alpha=0.35, label="sample")
ax.plot(np.mean(samples.squeeze(-1).T.cpu().numpy()[:, :], axis=1),color="green", lw=2, alpha=0.75, label="mean")
ax.set_xlabel("t in days")
ax.set_ylabel("aggregated sales")
labels = np.arange(len(covariates)-42, len(covariates), step=8)
ax.set_xticklabels(labels)
plt.ylim((9.5, 11.2))
plt.legend(loc='lower right')
plt.show()

承认

非常感谢 Pyro 文档和开发团队。GitHub 上关于主题和实现的全面写作确实使 Pyro 用户能够立即投入使用。我们的模型是不同模型的混合实验，我们测试了哪些部分可以很好地协同工作。

参考

1. Pyro 文档—预测教程 I 、 II 和 III
2. 米（meter 的缩写））戈里诺瓦，d .摩尔，m .霍夫曼。概率程序的自动重新参数化。2019 发表于 ArXiv
3. 米（meter 的缩写））霍夫曼等人在随机变分推理。 2013 年。在 JMLR
4. 优步技术公司。2018.烟火预报员见文件

PySnpTools

原文：https://towardsdatascience.com/pysnptools-9038f994fdd7?source=collection_archive---------56-----------------------

在 Python 中读取和操作基因组数据

国立癌症研究所在 Unsplash 拍摄的照片

PySnpTools 是一个 Python 库，用于读取和操作 Python 中的基因组数据。它允许用户有效地选择和重新排序个体(行)和 SNP 位置(列)。然后，它只读取选定的数据。PySnpTools 最初是为了支持快速 LMM——一种全基因组关联研究(GWAS)工具——而开发的，现在支持大内存和集群规模的工作。

2022 年 3 月 22 日更新:更新示例以匹配本笔记本。

介绍

当我们开发快速 LMM 包(1；2;3;4;5;6)，我们厌倦了重新编写代码来支持越来越多的文件格式(例如，PLINK 的 Bed，Dat，Pheno 等。) (7;8).此外，我们注意到自己重复地进行类似的操作，例如，

• 读取 SNP 子集(列)的数据，
• 对我们的表型数据中的个体(行)进行重新排序，以匹配它们在 SNP 数据中的顺序，以及
• 填补缺失数据和规范化数据。

受 NumPy 和 Pandas 的启发，我们创建了 PySnpTools，这是一个开源库，使这些操作变得简单。

我们在西雅图的 PyData 2015 大会上展示了 PySnpTools 的第一个公开版本。

在会议上，Travis oli phant(NumPy 的主要创建者)推荐了 PySnpTools 读取和操作基因组数据的方法(9)。

会议结束后不久，希拉里·菲纽肯(麻省理工学院和哈佛大学布罗德研究所菲纽肯实验室的负责人)写信给我们:

“我一直喜欢 PySnpTools，并向其他使用 Python (10)的统计遗传学家推荐它！”

自会议以来，随着 FaST-LMM 发展到支持多达 100 万个样本的数据集(6)，我们扩展了 PySnpTools 以简化:

• 大于内存的数据集
• 在多个处理器或任何集群上运行循环，以及
• 在本地或从/向任何远程存储器读写文件。

最近，我们已经使用多线程来加快读床速度，并增加了对 BGEN 格式的支持。

本文讲述了如何安装 PySnpTools，并描述了 PySnpTools 关注的基因组数据。然后给出了核心用法的例子。最后，它列出了 PySnpTools 的其他特性，并提供了与其他数据集方法的简要比较。

安装 PySnpTools

要使用 PySnpTools:

pip install pysnptools

你可以在https://github.com/fastlmm/PySnpTools找到 Jupyter 笔记本、完整的 API 文档和示例以及源代码。

基因组数据

我们感兴趣的基因组数据通常由 500，000 到 150 万列组成——每个 SNP 或变体一列(这是已知人类不同的基因组位置)。数据包括每个人的一行。一千到一百万行是典型的。初始数据中的值可能是 0、1、2 或缺失(代表在一个基因组位置为一个个体测量的次要等位基因的数量)。标准化后，值是 64 位或 32 位浮点数，缺少的值用 NaN(“非数字”)表示。

图 1 显示了 PySnpTools 对基因组数据的内存表示。被称为iid的两个字符串标识每个个体。一个名为**sid**的字符串识别每个 SNP。一种叫做**val**的浮动，告诉一个 SNP 个体的等位基因数。最后，名为**pos**的三重浮点表示每个 SNP 的位置(染色体数目、遗传距离和碱基对位置)。对于**val**中的一个给定位置，PySnpTools 可以轻松找到对应的**iid**、**sid**、**pos.** 而且，对于任何感兴趣的**iids** 或**sids** ，PySnpTools 都可以轻松高效地在**val**中找到对应的位置。

图一。PySnpTools 的基因组数据的内存表示

对于表型和协变量数据(例如，年龄、性别、身高、体重、疾病的存在)，PySnpTools 使用相同的表示法。对于这样的数据，**sid** 给出特征的名称(如，“高度”)，而**pos**被忽略。**val** 数组作为浮点数，也可以通过 0.0 和 1.0 表示二进制数据。我们通过 hot-one 编码处理分类数据。

读取基因组文件:一个例子

PySnpTools 将从文件中读取基因组数据的 Python 对象称为 SnpReader。下面是如何为流行的 Bed 文件格式创建 SnpReader:我们告诉它要读取什么文件以及如何读取。(这个 SnpReader 还没有从它的文件中读取任何内容。)

输入:

*# Download BED files and return local file name*
**from** **pysnptools.util** **import** example_file
all_bed_file = example_file("doc/ipynb/all.*",endswith="*.bed"**from pysnptools.snpreader import Bed** snpreader = Bed(all_bed_file, count_A1=True)
print(snpreader)

输出:

Bed("[...]/all.bed",count_A1=True)

我们可以向 SNP 阅读器询问个体的数量和 SNP 的数量。我们也可以，例如，要求它为第一个体的**iid** 。我们可以另外向它索要最后一个 SNP 的**sid** 。(Bed SnpReader 只读取回答这些问题所需的小文件。)

输入:

print(snpreader.iid_count,snpreader.sid_count)
print(snpreader.iid[0])
print(snpreader.sid[-1])

输出:

500 5000
['cid0P0' 'cid0P0']
snp124_m0_.23m1_.08

接下来，我们可以将所有基因组数据读入内存，创建一个新的 SnpReader，称为 SnpData。因为 SnpData 是 SnpReader，我们可以再次询问个体和 SNP 的数量。

输入:

snpdata = snpreader.read()
print(snpdata)
print(snpdata.iid_count, snpdata.sid_count)

输出:

SnpData(Bed('all.bed',count_A1=True))
500 5000

SnpData 是一个 SnpReader，它包含一个其他 SNP reader 没有的额外的**val** 属性。**val** 属性是基因组数据的(内存中)NumPy 数组。

例如，我们可以显示前 7 个个体和前 7 个 SNPs 的基因组数据。我们还可以找到所有基因组数据的平均值。

输入:

**import** **numpy** **as** **np** print(snpdata.val[:7,:7])
print(np.mean(snpdata.val))

输出:

[[0\. 0\. 1\. 2\. 0\. 1\. 2.]
 [0\. 0\. 1\. 1\. 0\. 0\. 2.]
 [0\. 0\. 1\. 2\. 1\. 0\. 0.]
 [0\. 0\. 0\. 2\. 0\. 0\. 0.]
 [0\. 0\. 0\. 0\. 0\. 0\. 2.]
 [0\. 0\. 1\. 0\. 0\. 0\. 2.]
 [0\. 0\. 2\. 1\. 0\. 1\. 2.]]
0.521412

或者，我们可以在一行中做所有的事情。

输入:

print(np.mean(Bed(all_bed_file,count_A1=**True**).read().val))

输出:

0.521412

SnpData

SnpData，特殊的内存 SnpReader，由任何 SnpReader 的**read()** 方法创建。我们还可以从头开始创建 SnpData。

这里我们为三个个体和两个 SNP 创建**snpdata1** 。我们用 NaN 标记一个缺失值，然后求平均值(忽略缺失值)。

输入:

**from pysnptools.snpreader import SnpData**
snpdata1 = SnpData(iid=[['f1','c1'],['f1','c2'],
                       ['f2','c1']],
                   sid=['snp1','snp2'],
                   val=[[0,1],[2,.5],[.5,np.nan]])
print(np.nanmean(snpdata1.val))

输出:

0.8

在读取前选择并重新排序数据

假设我们只关心前 7 个个体和前 7 个 SNP 的基因组数据。PySnpTools 使从磁盘中读取所需的数据变得很容易。我们在 read 方法之前使用类似 NumPy 的索引。

输入:

*snpreader = Bed(all_bed_file,count_A1=**True**)
snpdata77 = snpreader[:7,:7].read()
print(snpdata77.val)*

输出:

*[[0\. 0\. 1\. 2\. 0\. 1\. 2.]
 [0\. 0\. 1\. 1\. 0\. 0\. 2.]
 [0\. 0\. 1\. 2\. 1\. 0\. 0.]
 [0\. 0\. 0\. 2\. 0\. 0\. 0.]
 [0\. 0\. 0\. 0\. 0\. 0\. 2.]
 [0\. 0\. 1\. 0\. 0\. 0\. 2.]
 [0\. 0\. 2\. 1\. 0\. 1\. 2.]]*

支持所有类似 NumPy 的花哨索引:切片、布尔、整数列表、负整数(11)。此外，PySnpTools 允许同时对行和列进行有趣的索引。

将索引附加到任何 SnpReader，创建一个新的 SnpReader。这里有一个极端的例子。它说，“创造一个阅读器

• 来自 bed 格式的文件“all.bed”
• 个体以相反的顺序并且对于每第二个 SNP
• 对于前 5 个这样的个体和前 5 个这样的 SNP
• 对于第一个和最后一个个体以及第一个、第二个和第五个 SNP，
• 然后读。"

和以前一样，只从磁盘中读取最终所需的数据。

输入:

*print(Bed(all_bed_file,count_A1=**True**)[::-1,::2][:5,:5]
       [[0,-1],[True,True,False,False,True]]
       .read().val)*

输出:

*[[0\. 0\. 1.]
 [0\. 1\. 1.]]*

通过单个(行)和 SNP(列)标识符进行索引

我们已经看到了如何使用基于位置的索引来操作 SnpReaders。如果我们转而希望通过**iid** (个体的标识符)或**sid**(SNPs 的标识符)来操纵 SnpReaders 呢？这是一个更像熊猫的场景。在这种情况下，我们可以使用**iid_to_index()** 或**sid_to_index()** 的方法。

这个例子显示了如何读取三个感兴趣的 SNP，每个 SNP 由**sid**而不是位置识别。

输入:

***desired_sid_list = ['snp1750_m0_.02m1_.04',
              'snp0_m0_.37m1_.24','snp122_m0_.26m1_.34']
snpreader = Bed(all_bed_file,count_A1=**True**)
desired_snpreader = \
    snpreader[:,snpreader.sid_to_index(desired_sid_list)]
print(desired_snpreader.iid_count,
      desired_snpreader.sid_count)***

输出:

***500 3***

读者和伪读者

PySnpTools 包括这些常见文件格式的 SNP reader:Bed、Pheno、Dat、Ped、Dense。它还定义了这些新的二进制格式:SnpNpz、SnpHdf5 和 SnpMemMap。最后一点很有趣，因为它使用 NumPy 的内存映射数组来提供对大于内存容量的数据的内存访问。PySnpTools 还定义了一种称为 DistributedBed 的格式，这种格式提供了对密集存储在几十个(或几百或几千个)Bed 文件中的数据的随机访问。

此外，PySnpTools 还包括 SnpGen，可以称之为“伪阅读器”。对于任何使用 SnpGen 的程序来说，它看起来像一个文件阅读器，但它不是从磁盘读取数据，而是动态生成随机的基因组数据。基于用户提供的种子的这一代是确定性的。

除了阅读 SNPS，PySnpTools 还包括大量内核数据的读者和伪读者。(核数据表示个体之间成对的相似性。)

如果 PySnpTools 不支持感兴趣的数据格式，您有两种选择。如果您的数据相对较小，您可以使用其他 Python 工具读取它，然后创建内存中的 SnpData(或 SnpMemMap)对象。如果需要更高的性能，可以为该格式编写一个新的 SnpReader 模块。在这两种情况下，为 SnpReader 编写的任何程序都将自动使用新的数据格式。

超越阅读

最近，当我们努力将快速 LMM 扩展到 100 万个样本时，我们向 PySnpTools 添加了更多功能。PySnpTools 现在包括以下工具:

• 标准化 SNP 数据和内核，
• 交叉(和排序)来自任何数量 SnpReaders 的个体(在一行中)，
• 有效地处理大整数范围，
• 编写循环，以便它们可以在多个处理器或任何集群(由一个模块定义)上运行，以及
• 向任何本地或任何分布式文件系统(由模块定义)读写文件。

与其他数据集格式相比

在 PyData 2015 上，Joshua Bloom 比较了用于数据科学的数据集工具(12)。表 1 用 PySnpTools 更新了他的表。基本上，与 xarray 非常相似，PySnpTools 用不可变的索引行和列标识符来修饰 2d NumPy 数组。这给了我们方便和高效率，但比其他方法更特殊的目的。PySnpTools 还延迟了所有的复制(和读取)，直到**read()** 方法被显式调用，这给了我们效率和可预测性。

表 1。Bloom 的数据科学数据集，加上 PySnpTools

摘要

PySnpTools 是一个 Python 库，用于读取和操作基因组数据。受 NumPy 和 Pandas 的启发，PySnpTools 可以被认为是一种向基因组数据添加熊猫式索引的方法。要安装 PySnpTools:

pip install pysnptools

你可以在https://github.com/fastlmm/PySnpTools找到 Jupyter 笔记本、完整的 API 文档和示例以及源代码。

感谢

感谢 FaST-LMM 团队中的每个人，他们的使用和反馈为 PySnpTools 的设计提供了信息。特别感谢 Christoph Lippert，他编写了我们最初的 BED reader，以及 David Heckerman，他鼓励了 PySnpTools 的开发。

参考

1.用于全基因组关联研究的快速线性混合模型。【Lippert，c .】等* 2011，《自然方法》，8 833–835 页。***

2.对扩展的威康信托数据进行详尽的上位 SNP 关联分析。 Lippert，c .等人 2013 年，科学报告 3，1099。

3.基于核的遗传变异集关联测试的能力和计算效率更高。 李柏特，C，等 2014，生物信息学 30，22。

4.对全基因组关联研究线性混合模型的进一步改进。 Widmer，c .，等 2015，科学报告，4 6874。

5.遗传力估计的线性混合模型。 赫克曼，D. 等人。 2016 年美国国家科学院院刊 113 期(27)。

6.用于全基因组关联研究的 Ludicrous Speed 线性混合模型。 卡迪，C. &海克尔马，D. 2019，bioRxiv 154682。

7.珀塞尔，s .普林克。【http://zzz.bwh.harvard.edu/plink/.在线】

8. PLINK:全基因组关联和基于人群的工具集。 Purcell，S. 等 2007，美国人类遗传学杂志，81。

9.奥列芬特，t。2015.

10.个人通信。2015.

11.分度。* NumPy 用户指南。【在线】https://docs . scipy . org/doc/numpy/reference/arrays . indexing . html***

12.机器学习的系统观点。 布鲁姆，J. 西雅图:s.n .，2015。PyData。

PySpark 和 SparkSQL 基础知识

原文：https://towardsdatascience.com/pyspark-and-sparksql-basics-6cb4bf967e53?source=collection_archive---------0-----------------------

使用 PySpark 实施大数据

如何用 Python 编程实现 Spark

【来源】

Python 通过 Spark Python API 揭示了 Spark 编程模型来处理结构化数据，该 API 被称为 PySpark 。

这篇文章的目的是演示如何用 PySpark 运行 Spark 并执行常见的函数。

Python 编程语言需要安装 IDE。使用 Python 和 Anaconda 的最简单的方法，因为它安装了足够的 IDE 和重要的包。

1.下载 Anaconda 并安装 PySpark

借助这个 链接 ，可以下载 Anaconda。下载合适的 Anaconda 版本后，点击它继续安装程序，安装程序在 Anaconda 文档 中有详细说明。

安装完成后，将会打开 Anaconda Navigator 主页。要使用 Python，只需点击“笔记本”模块的“启动”按钮。

Anaconda Navigator 主页(图片由作者提供)

为了能够通过 Anaconda 使用 Spark，应该遵循下面的包安装步骤。

Anaconda 提示终端

康达安装 pyspark

康达安装 pyarrow

PySpark 和 PyArrow 包安装完成后，只需关闭终端，返回 Jupyter Notebook，在代码顶部导入所需的包。

import pandas as pd
from pyspark.sql import SparkSession
from pyspark.context import SparkContext
from pyspark.sql.functions 
import *from pyspark.sql.types 
import *from datetime import date, timedelta, datetime
import time

2.正在初始化 SparkSession

首先，需要初始化一个 Spark 会话。在 SparkSession 的帮助下，可以创建 DataFrame 并将其注册为表。而且执行 SQL 表，可以缓存表，可以读取 parquet/JSON/CSV/Avro 数据格式的文件。

sc = SparkSession.builder.appName("PysparkExample")\    
.config ("spark.sql.shuffle.partitions", "50")\    .config("spark.driver.maxResultSize","5g")\    
.config ("spark.sql.execution.arrow.enabled", "true")\    .getOrCreate()

关于 SparkSession 各个参数的详细解释，敬请访问py spark . SQL . spark session。

3.创建数据框

数据帧可以被接受为标题列的分布式列表集合，类似于关系数据库中的表。在本文中，我们将在 PySpark API 上使用 DataFrame 操作来处理数据集。

你可以从这个链接下载 Kaggle 数据集。

3.1。来自星火数据来源**

可以通过读取文本、CSV、JSON 和拼花文件格式来创建数据帧。在我们的例子中，我们将使用一个. json 格式的文件。您还可以通过使用如下所示的相关读取函数来查找和读取文本、CSV 和拼花文件格式。

**#Creates a spark data frame called as raw_data.****#JSON**
dataframe = sc.read.json('dataset/nyt2.json')**#TXT FILES#** 
dataframe_txt = sc.read.text('text_data.txt')**#CSV FILES#** 
dataframe_csv = sc.read.csv('csv_data.csv')**#PARQUET FILES#** 
dataframe_parquet = sc.read.load('parquet_data.parquet')

4.重复值

使用 dropDuplicates()函数可以消除表中的重复值。

dataframe = sc.read.json('dataset/nyt2.json') 
dataframe.show(10)

(图片由作者提供)

在应用 dropDuplicates()函数之后，我们可以观察到重复项被从数据集中删除了。

dataframe_dropdup = dataframe.dropDuplicates() dataframe_dropdup.show(10)

(图片由作者提供)

5.问题

查询操作可以用于各种目的，例如用 “select” 子集化列，用 “when” 添加条件，用 *"like "过滤列内容。*下面举例说明一些最常用的操作。有关查询操作的完整列表，请参见 Apache Spark 文档 。

5.1。【选择】操作**

可以通过属性(“作者”)或索引(dataframe[‘author’])来获取列。

**#Show all entries in title column**
dataframe.select("author").show(10)**#Show all entries in title, author, rank, price columns**
dataframe.select("author", "title", "rank", "price").show(10)

第一个结果表只显示“作者”选择，第二个结果表显示多列(作者的图像)

5.2。【当】操作

在第一个示例中，选择了“title”列，并添加了一个带有“when”条件的条件。

**# Show title and assign 0 or 1 depending on title**dataframe.select("title",when(dataframe.title != 'ODD HOURS', 
1).otherwise(0)).show(10)

显示 10 行指定条件(图片由作者提供)

在第二个示例中，应用了“isin”操作，而不是“when ”,后者也可用于为行定义一些条件。

**# Show rows with specified authors if in the given options**dataframe [dataframe.author.isin("John Sandford", 
"Emily Giffin")].show(5)

结果集显示 5 行指定的条件(图片由作者提供)

5.3。【喜欢】操作**

在“Like”函数的括号中,%字符用于过滤掉所有包含“the”单词的标题。如果我们正在寻找的条件是完全匹配的，那么不应该使用%字符。

**# Show author and title is TRUE if title has " THE " word in titles**dataframe.select("author", "title",
dataframe.title.like("% THE %")).show(15)

包含“THE”单词的标题的结果集。(图片由作者提供)

5.4。**

从单词/内容的开头开始扫描，并在括号中指定标准。并行地，EndsWith 从末尾开始处理单词/内容。这两个函数都区分大小写。

dataframe.select("author", "title", dataframe.title.startswith("THE")).show(5)dataframe.select("author", "title", dataframe.title.endswith("NT")).show(5)

结果集有 5 行 startsWith 和 endsWith 操作。(图片由作者提供)

5.5。【子串】操作

Substring 函数提取指定索引之间的文本。在以下示例中，从索引号(1，3)、(3，6)和(1，6)中提取文本。

dataframe.select(dataframe.author.substr(1, 3).alias("title")).show(5)dataframe.select(dataframe.author.substr(3, 6).alias("title")).show(5)dataframe.select(dataframe.author.substr(1, 6).alias("title")).show(5)

分别显示子字符串(1，3)、(3，6)、(1，6)的结果。(图片由作者提供)

6.添加、更新和删除列

DataFrame API 中也提供了数据操作函数。下面，您可以找到添加/更新/删除列操作的示例。

6.1。添加列

**# Lit() is required while we are creating columns with exact values.**dataframe = dataframe.withColumn('new_column', 
F.lit('This is a new column'))display(dataframe)

新列被添加到数据集的末尾(图片由作者提供)

6.2。更新列**

对于 DataFrame API 的更新操作，withColumnRenamed()函数使用两个参数。

**# Update column 'amazon_product_url' with 'URL'**dataframe = dataframe.withColumnRenamed('amazon_product_url', 'URL')dataframe.show(5)

“Amazon_Product_URL”列名更新为“URL”(图片由作者提供)

6.3。移除列**

可以通过两种方式删除列:在 drop()函数中添加列名列表，或者通过在 drop 函数中指向来指定列。两个示例如下所示。

dataframe_remove = dataframe.drop("publisher", "published_date").show(5)dataframe_remove2 = dataframe \ .drop(dataframe.publisher).drop(dataframe.published_date).show(5)

“publisher”和“published_date”列以两种不同的方法删除。(图片由作者提供)

7.检查数据

存在几种类型的函数来检查数据。下面，你可以找到一些常用的。要深入了解，请访问 Apache Spark 文档 。

**# Returns dataframe column names and data types**
dataframe.dtypes**# Displays the content of dataframe**
dataframe.show()**# Return first n rows**
dataframe.head()**# Returns first row**
dataframe.first()**# Return first n rows**
dataframe.take(5)**# Computes summary statistics**
dataframe.describe().show()**# Returns columns of dataframe**
dataframe.columns**# Counts the number of rows in dataframe**
dataframe.count()**# Counts the number of distinct rows in dataframe**
dataframe.distinct().count()**# Prints plans including physical and logical**
dataframe.explain(4)

8.“分组”操作

通过在函数中添加列名，使用 GroupBy()函数应用分组过程。

**# Group by author, count the books of the authors in the groups**dataframe.groupBy("author").count().show(10)

作者按出版的书籍数量分组(图片由作者提供)

9.“过滤”操作

通过使用 filter()函数并在其中添加一个条件参数来应用过滤。该函数区分大小写。

**# Filtering entries of title
# Only keeps records having value 'THE HOST'**dataframe.filter(dataframe["title"] == 'THE HOST').show(5)

标题列被过滤，内容只有“主机”并显示 5 个结果。(图片由作者提供)

10.缺少和替换值

对于每个数据集，在数据预处理阶段总是需要替换现有值、删除不必要的列和填充缺失值。py spark . SQL . dataframenafunction库帮助我们在这方面操作数据。下面补充一些例子。

**# Replacing null values**
dataframe.na.fill()
dataFrame.fillna()
dataFrameNaFunctions.fill()**# Returning new dataframe restricting rows with null** valuesdataframe.na.drop()
dataFrame.dropna()
dataFrameNaFunctions.drop()**# Return new dataframe replacing one value with another**
dataframe.na.replace(5, 15)
dataFrame.replace()
dataFrameNaFunctions.replace()

11.分配

可以增加或减少 RDD 中现有的分区级别。增加可以通过使用 repartition(self，numPartitions) 函数来实现，该函数会产生一个获得更多分区数量的新 RDD。减少可以用 coalesce(self，numPartitions，shuffle=False) 函数处理，该函数产生一个新的 RDD，其分区数量减少到指定的数量。更多信息请访问**Apache Spark docs。**

**# Dataframe with 10 partitions**
dataframe.repartition(10).rdd.getNumPartitions()**# Dataframe with 1 partition**
dataframe.coalesce(1).rdd.getNumPartitions()

12.以编程方式运行 SQL 查询

还可以通过在 SparkSession 上启用“sql”操作来使用原始 SQL 查询，从而以编程方式运行 SQL 查询，并将结果集作为 DataFrame 结构返回。更多详细信息，敬请访问Apache Spark docs。

**# Registering a table**
dataframe.registerTempTable("df")sc.sql("select * from df").show(3)sc.sql("select \               
CASE WHEN description LIKE '%love%' THEN 'Love_Theme' \               WHEN description LIKE '%hate%' THEN 'Hate_Theme' \               WHEN description LIKE '%happy%' THEN 'Happiness_Theme' \               WHEN description LIKE '%anger%' THEN 'Anger_Theme' \               WHEN description LIKE '%horror%' THEN 'Horror_Theme' \               WHEN description LIKE '%death%' THEN 'Criminal_Theme' \               WHEN description LIKE '%detective%' THEN 'Mystery_Theme' \               ELSE 'Other_Themes' \               END Themes \       
from df").groupBy('Themes').count().show()

13.输出

13.1。数据结构**

DataFrame API 使用 RDD 作为基础，它将 SQL 查询转换为低级 RDD 函数。通过使用。rdd 操作，一个数据帧可以转换成 RDD。也可以将 Spark Dataframe 转换成一串 RDD 和熊猫格式。

**# Converting dataframe into an RDD**
rdd_convert = dataframe.rdd**# Converting dataframe into a RDD of string** dataframe.toJSON().first()**# Obtaining contents of df as Pandas** 
dataFramedataframe.toPandas()

不同数据结构的结果(图片由作者提供)

13.2。写&保存到文件**

作为数据帧加载到我们代码中的任何数据源类型都可以很容易地转换并保存为其他类型，包括。parquet 和. json .更多保存、加载、写函数细节请访问Apache Spark doc。

**# Write & Save File in .parquet format**
dataframe.select("author", "title", "rank", "description") \
.write \
.save("Rankings_Descriptions.parquet")

拼花文件是在。写吧。save()函数被处理。(图片由作者提供)

**# Write & Save File in .json format**
dataframe.select("author", "title") \
.write \
.save("Authors_Titles.json",format="json")

JSON 文件是在。写吧。save()函数被处理。(图片由作者提供)

13.3。停止火花会话**

可以通过运行 stop() 功能来停止 Spark 会话，如下所示。

**# End Spark Session**
sc.stop()

代码和 Jupyter 笔记本都在我的 GitHub 上。

非常感谢您的提问和评论！

参考文献:

1. http://spark.apache.org/docs/latest/
2. https://docs.anaconda.com/anaconda/

Pyspark 数据操作教程

原文：https://towardsdatascience.com/pyspark-data-manipulation-tutorial-8c62652f35fa?source=collection_archive---------13-----------------------

入门

绝对火花初学者入门

图片由免费提供-照片来自 Pixabay

这个教程适合你吗？

本教程是为有一些 Python 经验的数据人员编写的，他们绝对是 Spark 初学者。它将帮助您安装 Pyspark 并启动您的第一个脚本。您将了解弹性分布式数据集(rdd)和数据帧，它们是 Pyspark 中的主要数据结构。我们讨论一些基本的概念，因为我认为它们会让你以后不会感到困惑和调试。您将学习数据转换以及从文件或数据库中读取数据。

为什么火花

学习 Spark 的主要原因是，您将编写可以在大型集群中运行并处理大数据的代码。本教程只讨论 Python API py Spark，但是你应该知道 Spark APIs 支持 4 种语言:Java、Scala 和 R 以及 Python。由于 Spark core 是用 Java 和 Scala 编程的，所以那些 API 是最完整和最有原生感觉的。

Pyspark 的优势在于 Python 已经有了许多数据科学的库，您可以将它们插入到管道中。再加上 Python 摇滚的事实！！！能让 Pyspark 真正有生产力。例如，如果你喜欢熊猫，你知道你可以用一个方法调用将 Pyspark 数据帧转换成熊猫数据帧。然后你可以对数据做一些事情，并用 matplotlib 绘制出来。

随着时间的推移，你可能会发现 Pyspark 几乎和 pandas 或 sklearn 一样强大和直观，并在你的大部分工作中使用它。试试 Pyspark 吧，它可能会成为你职业生涯中的下一件大事。

Pyspark 快速启动

有很多关于如何创建 Spark 集群、配置 Pyspark 向它们提交脚本等等的文章。所有这些都是在 Spark 上进行高性能计算所需要的。然而，在大多数公司中，他们会有数据或基础架构工程师来维护集群，而您只需编写和运行脚本。因此，您可能永远不会安装 Spark 集群，但是如果您需要，您可以在问题出现时处理它。OTOH，如果我告诉你，只需一个命令，你就可以在几分钟内开始工作，那会怎么样？

pip install pyspark

这将在您的机器上下载并安装 Pyspark，这样您就可以学习其中的诀窍，然后，当您得到一个集群时，使用它将会相对容易。

为了支持本教程，我创建了这个 GitHub 库:

[## ArmandoRiveroPi/py spark _ 教程

pyspark GitHub 中的数据操作和机器学习介绍是超过 5000 万开发人员的家园，他们致力于…

github.com](https://github.com/ArmandoRiveroPi/pyspark_tutorial)

示例中的大部分代码最好放在tutorial _ part _ 1 _ data _ wrangling . py文件中。

在加快速度之前，先给你点提示。我发现如果我不在主 Python 脚本的开头设置一些环境变量，Pyspark 会抛出错误。

import sys, os
environment = ['PYSPARK_PYTHON', 'PYSPARK_DRIVER_PYTHON']
*for* var *in* environment:
    os.environ[var] = sys.executable

之后，您可以创建 spark 会话

*from* pyspark.sql *import* SparkSession
session = SparkSession.builder.getOrCreate()

相信我，这是你开始工作所需要的。

基本概念

火花建筑

好了，我们可以进入正题了。首先，你必须知道 Spark 是复杂的。它不仅包含大量用于执行许多数据处理任务的代码，而且在执行分布在数百台机器上的代码时，它有望高效无误。在开始的时候，你真的不需要这种复杂性，你需要的是把事情抽象出来，忘记那些具体的细节。不过，重要的是司机和工人之间的区别。驱动程序是 Python (single)流程，您在其中为工人发布订单。工人是听从司机命令的火花过程。您可能会认为这些工人“在云中”(Spark 集群)。Spark 假设大数据将分布在工人中，这些工人有足够的内存和处理能力来处理这些数据。预计驱动程序没有足够的资源来保存这么多的数据。这就是为什么你需要明确地说你什么时候想要移动数据到驱动程序。我们一会儿会谈到这一点，但要注意，在实际应用中，你应该非常小心你带给驱动程序的数据量。

RDDs

弹性分布式数据集是 Spark 中最重要的数据结构之一，是数据帧的基础。您可以将它们视为“分布式”阵列。在许多方面，它们的行为类似于列表，一些细节我们将在下面讨论。

那么，如何创建 RDD 呢？最直接的方法是“并行化”Python 数组。

rdd = session.sparkContext.parallelize([1,2,3])

要开始与您的 RDD 互动，请尝试以下操作:

rdd.take(num=2)

这将为驾驶员带来 RDD 的前两个值。计数方法将返回 RDD 的长度

rdd.count()

如果你想把所有的 RDD 数据作为一个数组发送给驱动程序，你可以使用 collect

rdd.collect()

但是要小心，正如我们之前说过的，在实际应用中，这可能会使驱动程序崩溃，因为 RDD 的大小可能是千兆字节。一般来说，更喜欢使用 take 传递您想要的行数，只有当您确定 RDD 或数据帧不太大时才调用 collect 。

数据帧

如果你了解熊猫或 R 数据帧，你会非常清楚 Spark 数据帧代表什么。它们用命名列表示表格(矩阵)数据。在内部，它们是通过行对象的 RDD 实现的，这有点类似于名为 tuple 的 Python。数据框架的最大优势在于，它们使您能够将 SQL 思维付诸实践。我们稍后将讨论数据帧操作，但是让我们开始创建一个数据帧，以便您可以使用它。

df = session.createDataFrame(
  [[1,2,3], [4,5,6]], [‘column1’, ‘column2’, ‘column3’]
)

首先是数据矩阵，然后是列名。

可以像 RDD 案一样，尝试拿、数、收的方法；take and collect 将为您提供一个行对象列表。但是对我来说，最友好的显示方式应该是 show :

df.show(n=3)

它将用前 n 行打印数据帧的表格表示。

不变

Python 列表的一个关键区别是 rdd(以及数据帧)是不可变的。在并发应用程序和函数式语言中，经常需要不可变数据。

让我们讨论一下 Python 中不变性的含义。假设你这样做:

a = list(range(10))
a.append(11)

在这里，解释器首先创建一个名称“a ”,它指向一个列表对象，在第二行中，它修改了同一个对象。然而，当你这样做的时候

st = “my string”
st += “ is pretty”

解释器首先创建指向字符串对象的名称“st”，然后创建一个带有“my string is pretty”的全新字符串对象(在内存中的不同位置)，并将名称“**ST”**指向该新对象。发生这种情况是因为字符串是不可变的，你不能就地修改它们。

同样，rdd 和数据帧也不能就地修改，所以当您这样做时

my_rdd.map(*lambda* *x*: x*100)

我的 _rdd 不变。你需要做什么

my_rdd = my_rdd.map(*lambda* *x*: x*100)

以获得指向转换后的 rdd 对象的名称“my_rdd”。

转换和动作

与普通 Python 相比，一个更令人震惊的区别可能会让您在开始时感到困惑。有时您会注意到一个非常繁重的操作会立即发生。但是后来你做了一些小事情(比如打印 RDD 的第一个值),这似乎要花很长时间。

让我试着把这个问题简化很多来解释。在 Spark 中，变换和动作是有区别的。当你改变一个数据帧时，那是一个转换，然而，当你实际使用数据时(例如 df.show(1))，那是一个动作。转换是延迟加载的，当你调用它们时它们不会运行。当您通过操作使用它们的结果时，它们就会被执行。然后，所有需要的转换将被一起计划、优化和运行。所以当你看到看起来瞬时的操作时，即使它们很重，那也是因为它们是转换，它们将在以后运行。

数据帧操作

用户定义函数

用户定义的函数允许您使用 Python 代码对数据帧像元进行操作。你创建一个常规的 Python 函数，把它包装在一个 UDF 对象中，并把它传递给 Spark，它会让你的函数在所有的工作器中可用，并安排它的执行来转换数据。

*import* pyspark.sql.functions *as* funcs
*import* pyspark.sql.types *as* typesdef multiply_by_ten(number):
    return number*10.0multiply_udf = funcs.udf(multiply_by_ten, types.DoubleType())transformed_df = df.withColumn(
    'multiplied', multiply_udf('column1')
)
transformed_df.show()

首先你创建一个 Python 函数，它可以是一个对象中的方法，也是一个函数。然后创建一个 UDF 对象。恼人的部分是您需要定义输出类型，这是我们在 Python 中不太习惯的。要真正有效地使用 UDF，你需要学习这些类型，特别是复合映射类型(如字典)和数组类型(如列表)。这样做的好处是，您可以将这个 UDF 传递给 dataframe，告诉它将对哪一列进行操作，这样您就可以在不离开旧 Python 的舒适环境的情况下完成奇妙的事情。

然而，UDF 的一个主要限制是，尽管它们可以接受几列作为输入，但是它们不能整体改变行。如果您想要处理整行，您将需要 RDD 地图。

RDD 制图

这很像使用 UDF，你也传递一个常规的 Python 函数。但是在这种情况下，该函数将接收一个完整的行对象，而不是列值。预计它也将返回一整行。这将赋予您对行的最终权力，但有几个注意事项。首先:Row 对象是不可变，所以您需要创建一个全新的行并返回它。第二:你需要将数据帧转换成 RDD，然后再转换回来。幸运的是，这些问题都不难克服。

让我向您展示一个函数，它将对数转换您的数据帧中的所有列。作为一个不错的尝试，它还将每个列名转换为“log(column_name)”。我发现最简单的方法是用 row.asDict()将行放入字典。这个函数有点 Python 的魔力，比如字典理解和关键字参数打包(双星形)。希望你能全部理解。

*import* pyspark.sql.types *as* types*def* take_log_in_all_columns(row: types.Row):
     old_row = row.asDict()
     new_row = {f'log({column_name})': math.log(value) 
                *for* column_name, value *in* old_row.items()}
     *return* types.Row(**new_row)

这本身不会做任何事情。你需要执行地图。

logarithmic_dataframe = df.rdd.map(take_log_in_all_columns).toDF()

您会注意到这是一个链式方法调用。首先调用 rdd，它将为您提供存储数据帧行的底层 RDD。然后在这个 RDD 上应用 map，在这里传递函数。要关闭，您可以调用 toDF()，它将行的 RDD 转换为数据帧。进一步讨论见这个栈溢出问题。

SQL 操作

由于数据帧表示表，自然地，它们被赋予了类似 SQL 的操作。为了让你兴奋，我将只提到其中的几个，但是你可以期待找到几乎同构的功能。调用 select 将返回一个只有一些原始列的数据帧。

df.select('column1', 'column2')

对的调用，其中将返回一个数据帧，其中只有列 1 的值为 3 的行。

df.where('column1 = 3')

这个对 join 的调用将返回一个 dataframe，也就是，嗯…，一个 df 和 df1 通过 column1 的连接，就像 SQL 中的内部连接一样。

df.join(df1, [‘column1’], how=’inner’)

如果您需要执行右或左连接，在我们的例子中，df 就像左边的表，而 df1 就是右边的表。外部连接也是可能的。

但是除此之外，在 Spark 中，您可以更直接地执行 SQL。您可以从数据帧中创建时态视图

df.createOrReplaceTempView(“table1”)

然后在视图上执行查询

df2 = session.sql("SELECT column1 AS f1, column2 as f2 from table1")

这些查询将返回一个具有相应列名和值的新数据帧。

数据帧列操作

列抽象支持像算术这样的直接操作，例如，假设您想将列 1 相加，再加上列 2 和列 3 的乘积。那你可以

df3 = df.withColumn(
    'derived_column', df['column1'] + df['column2'] * df['column3']
)

聚合和快速统计

对于这个例子，我们需要从一个 CSV 文件中读取。我从https://archive.ics.uci.edu/ml/datasets/adult下载了成人数据文件

我将文件移动到我的存储库中的数据文件夹，并将其重命名为成人.数据. csv。

ADULT_COLUMN_NAMES = [
     "age",
     "workclass",
     "fnlwgt",
     "education",
     "education_num",
     "marital_status",
     "occupation",
     "relationship",
     "race",
     "sex",
     "capital_gain",
     "capital_loss",
     "hours_per_week",
     "native_country",
     "income"
 ]

读取文件后，我使用这个列表设置列名，因为文件没有标题。

csv_df = session.read.csv(
     'data/adult.data.csv', header=*False*, inferSchema=*True*
 )*for* new_col, old_col *in zip*(ADULT_COLUMN_NAMES, csv_df.columns):
     csv_df = csv_df.withColumnRenamed(old_col, new_col)

之后，您可以尝试下面的一些快速描述性统计

csv_df.describe().show()

若要获取聚合，请将 groupBy 与 agg 方法一起使用。例如，下面将为您提供一个数据框架，其中包含按年龄组划分的平均工作时间和标准偏差。我们也按年龄对数据帧进行排序。

work_hours_df = csv_df.groupBy(
    'age'
).agg(
    funcs.avg('hours_per_week'),
    funcs.stddev_samp('hours_per_week')
).sort('age')

连接到数据库

您很可能将数据保存在关系数据库中，由 MySQL 或 PostgreSQL 之类的 RDBMS 处理。如果是这样的话，不要担心，Spark 有办法与多种数据存储进行交互。我敢打赌，你可以通过谷歌搜索到你可能有的大多数 IO 需求。

要连接到 RDBMS，您需要一个 JDBC 驱动程序，这是一个 Spark 可以用来与数据库对话的 jar 文件。我将向您展示一个 PostgreSQL 示例。

首先你需要下载驱动。

在我的 GitHub repo 中，我把 JDBC 驱动程序放在 bin 文件夹中，但是你可以在任何路径中找到它。此外，您应该知道本例中的代码位于 repo 中的另一个文件中。我必须这样做，因为与其他文件不同，这个文件会抛出错误，因为数据库凭证当然是假的。您需要用您的数据库来替换它们以使其运行。

现在，您将通过一个额外的步骤启动 Pyspark:

session = SparkSession.builder.config(
    'spark.jars', 'bin/postgresql-42.2.16.jar'
).config(
    'spark.driver.extraClassPath', 'bin/postgresql-42.2.16.jar'
).getOrCreate()

这将创建一个可识别驱动程序的会话。然后，您可以像这样从数据库中读取(用您的真实配置替换假配置):

url = f"jdbc:postgresql://your_host_ip:5432/your_database"
properties = {'user': 'your_user', 'password': 'your_password'}
*# read from a table into a dataframe* df = session.read.jdbc(
    url=url, table='your_table_name', properties=properties
)

然后，您可以以任何方式创建转换后的数据帧，并将数据写回数据库(可能在不同的表中)。

transformed_df.write.jdbc(
    url=url, table='new_table', mode='append', properties=properties
)

根据文档，写入模式有:

• append:将此数据帧的内容追加到现有数据中。
• 覆盖:覆盖现有数据。
• 忽略:如果数据已经存在，则忽略此操作。
• error 或“errorifexists”(默认情况):如果数据已经存在，则抛出异常。

结论

我希望你受到鼓励，马上开始学习 Pyspark。我的意思是，这将需要几周的时间来提高效率，这取决于你每天能花多少时间，但之后你就可以用 Pyspark 来执行大量的数据处理了。万一你不能选择你在工作中使用的，至少你可以在你的空闲时间学习它。我的介绍几乎没有触及表面，主要是希望给你一个诱人的咬。比熊猫或者 sklearn 难不了多少。

熟练掌握 Spark，它可以为您打开大数据的大门。

第 2 部分将介绍基本的分类和回归。

进一步阅读

PySpark 食谱作者 Raju Kumar Mishra。Apress，2018。

官方文件

【https://spark.apache.org/docs/latest/api/python/