TowardsDataScience 2023 博客中文翻译（二百二十三）

原文：TowardsDataScience

协议：CC BY-NC-SA 4.0

精通逻辑回归

原文：towardsdatascience.com/mastering-logistic-regression-3e502686f0ae

从理论到 Python 实现

Roi Yehoshua 博士

·发表于 Towards Data Science ·17 分钟阅读·2023 年 5 月 20 日

–

图片来自 Gerd Altmann 的 Pixabay

逻辑回归是最常见的机器学习算法之一。它可以用来预测事件发生的概率，例如预测来邮件是否为垃圾邮件，或者肿瘤是否为恶性肿瘤，基于给定的标记数据集。

由于其简单性，逻辑回归通常被用作评估其他更复杂模型的基准。

该模型的名称中包含“逻辑”一词，因为它使用逻辑函数（Sigmoid）将输入特征的线性组合转换为概率。

它的名称中也包含“回归”一词，因为它的输出是一个介于 0 和 1 之间的连续值，尽管它通常作为二分类器使用，通过选择一个阈值（通常为 0.5），并将概率大于阈值的输入分类为正类，而将低于阈值的输入分类为负类。

在本文中，我们将深入讨论逻辑回归模型，从头开始在 Python 中实现它，然后展示其在 Scikit-Learn 中的实现。

背景：二分类问题

回顾一下在 监督机器学习 问题中，我们会得到一个包含 n 个标记样本的训练集：D = {(x₁, y₁), (x₂, y₂), … , (xₙ, yₙ)}，其中 xᵢ 是一个 m 维向量，包含样本 i 的特征，yᵢ 代表该样本的标签。我们的目标是构建一个预测尽可能接近真实标签的模型。

在分类问题中，标签 yᵢ 可以取 k 个值之一，表示样本所属的 k 个类别。更具体地说，在二分类问题中，标签 yᵢ 只能取两个值：0（表示负类）和 1（表示正类）。

此外，我们区分两种类型的分类器：

1. 确定性分类器为每个样本输出一个硬标签，而不提供类别的概率估计。这类分类器的例子包括感知机、K-近邻和支持向量机（SVM）。

2. 概率分类器输出类别的概率估计，然后根据这些概率为给定样本分配标签（通常是具有最高概率的类别标签）。这类分类器的例子包括逻辑回归、朴素贝叶斯分类器和使用 sigmoid 或 softmax 作为输出层的神经网络。

逻辑回归模型

逻辑回归是一种处理二分类问题的概率分类器。给定一个样本 (x, y)，它输出样本属于正类的概率 p：

如果这个概率高于某个阈值（通常选择为 0.5），则样本被分类为 1，否则被分类为 0。

模型如何估计概率p？

逻辑回归的基本假设是样本属于正类的事件的对数赔率是其特征的线性组合。

对数赔率（也称为logit）是赔率比的对数，赔率比是样本属于正类的概率与样本属于负类的概率之间的比率：

对数赔率（logit）函数

我们在这里假设对数的底数是 e（即自然对数），尽管也可以使用其他底数。

logit 函数的图像如下所示：

logit 函数

如图所示，logit 函数将 (0, 1) 区间的概率值映射到 (-∞, +∞) 区间的实数值。

在逻辑回归中，我们假设对数赔率是特征的线性组合，即：

其中 w = (w₀, …, wₘ) 是模型的参数（或权重）。参数 w₀ 通常被称为截距（或偏置）。

对于 p = 0.5（即对数赔率等于 0）的点定义了两个类别之间的分隔超平面，其方程为：

分隔超平面的方程

权重向量 w 与此超平面正交。超平面上方的每个样本 (wᵗx > 0**)** 被分类为正样本，而超平面下方的每个样本 (wᵗx < 0**)** 被分类为负样本：

这使得逻辑回归成为一种线性分类器，因为它假设类别之间的边界是一个线性表面。其他线性分类器包括感知器和支持向量机（SVM）。

我们可以通过对对数几率方程两边取指数，找到 p 和参数 w 之间的直接关联：

其中 σ 是Sigmoid 函数（也称为逻辑函数）：

Sigmoid 函数用于将对数几率 (wᵗx) 转换为概率。它具有一个特征“ S”形曲线：

Sigmoid 函数

如可以看出，该函数将实数范围 (-∞, +∞) 映射为概率值范围 (0, 1)。

Sigmoid 函数具有一些优良的数学性质，这些性质在后面会很有用：

下图总结了从输入到最终预测的逻辑回归计算过程：

逻辑回归模型

对数损失

我们的目标是找到参数w，使模型的预测 p = σ(wᵗx) 尽可能接近真实标签 y。为此，我们需要定义一个损失函数，用来衡量模型预测与真实标签之间的差距。这个函数需要是可微分的，以便可以使用如梯度下降等技术进行优化（有关机器学习中的损失函数的更多信息，请参见这篇文章）。

逻辑回归使用的损失函数称为对数损失（或逻辑损失）。其定义如下：

对数损失函数

我们是如何得到这个函数的？这个函数的推导基于最大似然原理。更具体地说，我们可以证明对数损失是在标签具有伯努利分布（即一种二元随机变量的概率分布，其中 1 的概率为 p，0 的概率为 1 − p）假设下的负对数似然。

从数学上讲，我们将证明：

其中 P(y|p) 是在模型预测 p 给定的情况下获得标签 y 的概率（即，数据在我们的模型下的似然）。

证明：

给定一个具有参数 p 的贝努利分布模型（标签），样本属于正类的概率就是 p，即，

类似地，样本属于负类的概率是：

我们可以将这两个方程更紧凑地写成如下形式：

解释：当 y = 1 时，pʸ = p 和 (1 − p)¹⁻ʸ = 1，因此 P(y|p) = p。类似地，当 y = 0 时，pʸ = 1 和 (1 − p)¹⁻ʸ = 1 − p，因此 P(y|p) = 1 − p。

因此，给定模型的数据的对数似然是：

对数损失恰好是该函数的负值。因此，最大化对数似然等同于最小化对数损失。

以下图展示了当 y = 1 时的对数损失：

对数损失仅在预测完全准确时（p = 1 且 y = 1，或 p = 0 且 y = 0）为 0，并且当预测变差时（即，当 y = 1 且 p → 0 或 y = 0 且 p → 1）接近无穷大。

成本函数 计算整个数据集上的平均损失：

这个函数可以用向量化的形式表示如下：

其中 y = (y₁, …, yₙ) 是一个包含所有训练样本标签的向量，而 p = (p₁, …, pₙ) 是一个包含模型对所有训练样本的预测概率的向量。

这个成本函数是凸的，即，它有一个全局最小值。然而，由于对数函数引入的非线性，没有封闭形式的解决方案来找到最佳 w*。因此，我们需要使用迭代优化方法如梯度下降来找到最小值。

梯度下降

梯度下降是一种迭代方法，用于寻找函数的最小值，其中我们沿着梯度的相反方向采取小步，以接近最小值：

梯度下降

为了使用梯度下降法找到成本 J(w)* 的最小值，我们需要计算其相对于每一个权重的偏导数。J(w)* 对于给定权重 wⱼ 的偏导数为：

证明：

因此，梯度向量可以用向量化的形式表示如下：

梯度下降更新规则为：

其中α是一个学习率，控制步长（0 < α < 1）。

请注意，每当你使用梯度下降时，你必须确保数据集已经标准化（否则梯度下降可能会在不同方向上采取不同大小的步伐，这会导致不稳定）。

Python 实现

现在我们将从头实现逻辑回归模型，包括成本函数和梯度计算，使用梯度下降优化模型，模型评估以及绘制最终的决策边界。

为了演示，我们将使用Iris 数据集（BSD 许可证）。原始数据集包含 150 个属于三种花卉之一的鸢尾花样本（山鸢尾、变色鸢尾和维吉尼亚鸢尾）。我们将其转化为一个二分类问题，只使用前两种花卉（山鸢尾和变色鸢尾）。此外，我们只使用每朵花的前两个特征（花萼宽度和花萼长度）。

加载数据集

首先，我们导入所需的库并固定随机种子，以获得可重复的结果：

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

np.random.seed(0)

接下来，我们加载数据集：

from sklearn.datasets import load_iris

iris = load_iris()
X = iris.data[:, :2]  # Take only the first two features
y = iris.target

# Take only the setosa and versicolor flowers
X = X[(y == 0) | (y == 1)]
y = y[(y == 0) | (y == 1)]

让我们绘制数据：

def plot_data(X, y):
    sns.scatterplot(x=X[:, 0], y=X[:, 1], hue=iris.target_names[y], style=iris.target_names[y], 
                    palette=['r','b'], markers=('s','o'), edgecolor='k')
    plt.xlabel(iris.feature_names[0])
    plt.ylabel(iris.feature_names[1])
    plt.legend() 

plot_data(X, y)

鸢尾花数据集

如可以看到，数据集是线性可分的，因此逻辑回归应该能够找到两个类别之间的边界。

接下来，我们需要向特征矩阵X中添加一列 1，以表示偏置项（w₀）：

# Add a column for the bias
n = X.shape[0] 
X_with_bias = np.hstack((np.ones((n, 1)), X))

现在我们将数据集分成训练集和测试集：

from sklearn.model_selection import train_test_split

X_train, X_test, y_train, y_test = train_test_split(X_with_bias, y, random_state=0)

模型实现

我们现在准备实现逻辑回归模型。我们从定义一个辅助函数来计算 sigmoid 函数开始：

def sigmoid(z):
    """ Compute the sigmoid of z (z can be a scalar or a vector). """
    z = np.array(z)
    return 1 / (1 + np.exp(-z))

接下来，我们实现成本函数，该函数返回给定数据集（X，y）上具有参数w的逻辑回归模型的成本，以及相对于w的梯度。

def cost_function(X, y, w):
    """ J, grad = cost_function(X, y, w) computes the cost of a logistic regression model 
        with parameters w and the gradient of the cost w.r.t. to the parameters. """
    # Compute the cost
    p = sigmoid(X @ w)    
    J = -(1/n) * (y @ np.log(p) + (1-y) @ np.log(1-p)) 

    # Compute the gradient    
    grad = (1/n) * X.T @ (p - y)  
    return J, grad

请注意，我们正在使用之前展示的成本函数和梯度函数的向量化形式。

为了对这个函数进行合理性检查，我们来计算模型在某个随机权重向量上的成本和梯度：

w = np.random.rand(X_train.shape[1])
cost, grad = cost_function(X_train, y_train, w)

print('w:', w)
print('Cost at w:', cost)
print('Gradient at w:', grad)

我们得到的输出是：

w: [0.5488135  0.71518937 0.60276338]
Cost at w: 2.314505839067951
Gradient at w: [0.36855061 1.86634895 1.27264487]

梯度下降实现

我们现在将实现梯度下降，以找到最优的w，使成本函数在给定训练集上最小化。该算法最多会对训练集进行max_iter次迭代（默认为 5000），除非成本在上一次迭代后没有至少减少tol（默认为 0.0001），在这种情况下训练将立即停止。

def optimize_model(X, y, alpha=0.01, max_iter=5000, tol=0.0001):
    """ Optimize the model using gradient descent.
        X, y: The training set        
        alpha: The learning rate
        max_iter: The maximum number of passes over the training set (epochs)
        tol: The stopping criterion. Training will stop when (new_cost > cost - tol)
    """
    w = np.random.rand(X.shape[1])
    cost, grad = cost_function(X, y, w)

    for i in range(max_iter):
        w = w - alpha * grad
        new_cost, grad = cost_function(X, y, w)        
        if new_cost > cost - tol:
            print(f'Converged after {i} iterations')
            return w, new_cost
        cost = new_cost

    print('Maximum number of iterations reached')
    return w, cost

通常在这一点上你需要对数据集进行标准化，因为梯度下降对于具有不同尺度的特征效果不好。在我们的特定数据集中，由于两个特征的范围相似，因此标准化不是必需的。

现在让我们调用这个函数来优化我们的模型：

opt_w, cost = optimize_model(X_train, y_train)

print('opt_w:', opt_w)
print('Cost at opt_w:', cost)

算法在 1,413 次迭代后收敛，我们得到的w*是：

Converged after 1413 iterations
opt_w: [ 0.28014029  0.80541854 -1.48367938]
Cost at opt_w: 0.28389717767222555

还有其他优化器可以使用，这些优化器通常比梯度下降更快，例如共轭梯度（CG）和截断牛顿（TNC）。有关如何使用这些优化器的更多细节，请参见scipy.optimize.minimize。

使用模型进行预测

现在我们已经找到了模型的最佳参数，可以使用它进行预测。

首先，我们编写一个函数，它接受新样本的矩阵X并返回它们属于正类的概率：

def predict_prob(X, w):
    """ Return the probability that samples in X belong to the positive class
        X: the feature matrix (every row in X represents one sample)
        w: the learned logistic regression parameters  
    """
    p = sigmoid(X @ w)
    return p

该函数通过简单地计算Xᵗw的 sigmoid 值来计算模型的预测（即对矩阵中每一行x计算σ(wᵗx)）。

例如，让我们找出位于(6, 2)的样本属于 versicolor 类的概率：

predict_prob([[1, 6, 2]], opt_w)

array([0.89522808])

这个样本有 89.52%的概率属于 versicolor 花。这是合理的，因为这个样本位于 versicolor 花的区域内，远离类别之间的边界。

另一方面，位于(5.5, 3)的样本属于 versicolor 类的概率是：

predict_prob([[1, 5.5, 3]], opt_w)

array([0.56436688])

这次概率要低得多（仅 56.44%），因为这个样本接近类别之间的边界。

让我们编写另一个函数，它返回预测的类别标签而不是概率：

def predict(X, w):
    """ Predict whether the label is 0 or 1 for the samples in X using a threshold of 0.5
        (i.e., if sigmoid(X @ theta) >= 0.5, predict 1)
    """
    p = sigmoid(X @ w)
    y_pred = (p >= 0.5).astype(int)
    return y_pred

该函数简单地在正类的概率至少为 0.5 时预测 1，否则预测 0。

让我们用上面的样本测试这个函数：

predict([[1, 6, 2], [1, 5.5, 3]], opt_w)

array([1, 1])

如预期的那样，这两个样本都被分类为 1。

评估模型

接下来，让我们编写一个函数来计算模型在给定数据集上的准确性：

def evaluate_model(X, y, w):
    y_pred = predict(X, w)
    accuracy = np.mean(y == y_pred)
    return accuracy

该函数首先找到模型在给定数据集X上的预测标签，并将其与真实标签y进行比较。然后计算准确性，作为正确分类的平均数量：

让我们使用这个函数来找出模型在训练集和测试集上的准确性：

train_accuracy = evaluate_model(X_train, y_train, opt_w)
print(f'Train accuracy: {train_accuracy * 100:.3f}%')

Train accuracy: 98.667%

test_accuracy = evaluate_model(X_test, y_test, opt_w)
print(f'Test accuracy: {test_accuracy * 100:.3f}%')

Test accuracy: 100.000%

如预期的那样，由于数据集是线性可分的，得分非常高。

除了准确性，还有其他重要指标用于评估分类模型，如精确度、召回率和 F1 分数。这些指标将在未来的文章中讨论。

绘制决策边界

最后，由于我们的数据集是二维的，我们可以绘制模型找到的类别之间的边界线。为此，我们首先需要找到这条线的方程。

边界线由模型预测值恰好为 0.5 的点定义，即：

当 sigmoid 函数的输入等于 0 时，其值为 0.5，因此我们可以写成：

重新排列项后我们得到：

即，边界线的斜率为 -w₁/w₂，截距为 -w₀/w₂。我们现在可以编写一个绘制这条线的函数：

def plot_decision_boundary(X, y, w):
    """ Plot the decision boundary between the classes """
    plot_data(X, y)

    line_x = np.array(plt.gca().get_xlim())
    line_y = -1 / w[2] * (w[1] * line_x + w[0])
    plt.plot(line_x, line_y, c='k', ls='--')

plot_decision_boundary(X, y, opt_w)

类别之间的决策边界

我们可以看到，只有一个样本被模型误分类。训练模型更多的迭代（约 200,000 次）会找到一个完美分离两类的分界线。使用固定步长时，梯度下降的最优收敛速度非常慢。可以通过使用自适应学习率（例如，使用更激进的步长来补偿快速消失的梯度）来改进这一点。

Scikit-Learn 中的 LogisticRegression 类

尽管从头实现逻辑回归有其自身的教育意义，但更实际的选择是使用 Scikit-Learn 提供的现成的 LogisticRegression 类。该类使用比普通梯度下降更高效的解算器，并且还提供了额外的选项，如正则化和提前停止。

该类的重要超参数有：

• penalty — 指定要应用的正则化类型。可以是以下选项之一：None, ‘l2’（默认值）, ‘l1’ 和 ‘elasticnet’。有关正则化的更多信息，请参见 这篇文章。

• tol — 停止准则的容忍度（默认为 0.0001）。

• C — 正则化系数的倒数（默认为 1.0）。较小的值表示更强的正则化。

• solver — 用于优化的算法。可以选择以下选项之一：‘lbfgs’（默认值），‘liblinear’，‘newton-cg’，‘newton-cholesky’，‘sag’，‘saga’。有关这些优化器的更多信息，请阅读文档。

• max_iter — 解算器收敛的最大迭代次数（默认为 100）

• multi_class — 处理多分类问题的方法。可以选择以下选项之一：‘ovr’（一对其余，即为每个类别构建一个二分类器与其他类别对抗）、‘multinomial’（使用多项式逻辑回归）或 ‘auto’（默认）。

使用 LogisticRegression 类时，您无需手动将一列全是 1 的列添加到设计矩阵 X 中，因为这会由 Scikit-Learn 自动完成。因此，在构建模型之前，我们将原始数据（没有额外的全是 1 的列）分割成训练集和测试集：

from sklearn.model_selection import train_test_split

X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)

我们将使用默认设置创建一个 LogisticRegression 实例，并将其拟合到训练集上：

from sklearn.linear_model import LogisticRegression

clf = LogisticRegression()
clf.fit(X_train, y_train)

接下来，我们将在训练集和测试集上评估模型：

train_accuracy = clf.score(X_train, y_train)
print(f'Train accuracy: {train_accuracy * 100:.3f}%')

test_accuracy = clf.score(X_test, y_test)
print(f'Test accuracy: {test_accuracy * 100:.3f}%')

Train accuracy: 100.000%
Test accuracy: 100.000%

这次我们在训练集和测试集上都获得了完美的分数。我们还可以通过查询 n_iter_ 属性来检查收敛所需的迭代次数：

print(clf.n_iter_)

[15]

仅需 15 次迭代即可收敛！显然，LogisticRegression 使用的求解器（默认使用L-BFGS）比我们实现的梯度下降更高效。

我们可以像之前一样绘制模型找到的决策边界。然而，这次最佳系数存储在模型的两个不同属性中：

• coef_ 是一个数组，包含所有权重，除了截距项

• intercept_ 是截距项 (w₀)

因此，我们需要在调用plot_decision_boundary()函数之前，将这两个属性连接成一个数组：

opt_w = np.insert(clf.coef_, 0, [clf.intercept_])
plot_decision_boundary(X, y, opt_w)

LogisticRegression 找到的决策边界

正如预期的那样，LogisticRegression 找到的线完美地分隔了两个类别。

总结

让我们总结一下逻辑回归与其他分类模型的优缺点。

优点：

• 当数据是线性可分的时，算法保证能找到一个类间的分离超平面。

• 提供类概率估计

• 不容易过拟合（但通常对数据存在欠拟合）

• 高度可解释（与每个特征相关的权重表示其重要性）

• 高度可扩展（需要的参数数量与特征数量线性相关）

• 可以处理冗余特征（通过赋予它们接近 0 的权重）

• 超参数数量较少

缺点：

• 只能找到类别之间的线性决策边界

• 通常被更复杂的模型超越

• 仅支持二分类，但可以扩展到多分类。将逻辑回归扩展到多分类问题（称为多项式逻辑回归或softmax 回归）在这篇文章中有介绍。

• 无法处理缺失值

最终备注

除非另有说明，所有图像均由作者提供。

本文的代码示例可以在我的 github 上找到：github.com/roiyeho/medium/tree/main/logistic_regression

感谢阅读！

精通 Python 中的长短期记忆：释放 LSTM 在 NLP 中的力量

原文：towardsdatascience.com/mastering-long-short-term-memory-with-python-unleashing-the-power-of-lstm-in-nlp-381ec3430f50

一本关于理解和实现 LSTM 层用于 Python 自然语言处理的全面指南

Eligijus Bujokas

·发布于Towards Data Science ·17 分钟阅读·2023 年 11 月 28 日

–

Sven Brandsma拍摄的照片，来自Unsplash。

这项工作是我关于 RNN 和 Python NLP 的文章的延续。一个简单的递归层的深度学习网络的自然发展是带有长短期记忆（LSTM）层的深度学习网络。

与 RNN 和 NLP 一样，我将尝试详细解释 LSTM 层，并从头编写该层的前向传播代码。

所有代码可以在这里查看: github.com/Eligijus112/NLP-python

我们将使用与上一篇文章相同的¹数据集：

# Data wrangling
import pandas as pd

# Reading the data 
d = pd.read_csv('input/Tweets.csv', header=None)

# Adding the columns 
d.columns = ['INDEX', 'GAME', "SENTIMENT", 'TEXT']

# Leaving only the positive and the negative sentiments 
d = d[d['SENTIMENT'].isin(['Positive', 'Negative'])]

# Encoding the sentiments that the negative will be 1 and the positive 0
d['SENTIMENT'] = d['SENTIMENT'].apply(lambda x: 0 if x == 'Positive' else 1)

# Dropping missing values
d = d.dropna()

数据集中的随机行；作者拍摄的图片

请记住，SENTIMENT=1 表示负面情感，SENTIMENT=0 表示正面情感。

我们需要将文本数据转换为整数序列。不过，与上一篇文章不同的是，我们现在将创建一个字符序列，而不是单词序列。

例如，文本“Nice Game”可以转换为以下示例向量：

[1, 2, 3, 4, 5, 6, 7, 8, 3]

每个字符，包括空格和标点符号，都会有一个索引。

def create_word_index(
    x: str, 
    shift_for_padding: bool = False, 
    char_level: bool = False) -> Tuple[dict, dict]: 
    """
    Function that scans a given text and creates two dictionaries:
    - word2idx: dictionary mapping words to integers
    - idx2word: dictionary mapping integers to words

    Args:
        x (str): text to scan
        shift_for_padding (bool, optional): If True, the function will add 1 to all the indexes.
            This is done to reserve the 0 index for padding. Defaults to False.
        char_level (bool, optional): If True, the function will create a character level dictionary.

    Returns:
        Tuple[dict, dict]: word2idx and idx2word dictionaries
    """
    # Ensuring that the text is a string
    if not isinstance(x, str):
        try: 
            x = str(x)
        except:
            raise Exception('The text must be a string or a string convertible object')

    # Spliting the text into words
    words = []
    if char_level:
        # The list() function of a string will return a list of characters
        words = list(x)
    else:
        # Spliting the text into words by spaces
        words = x.split(' ')

    # Creating the word2idx dictionary 
    word2idx = {}
    for word in words: 
        if word not in word2idx: 
            # The len(word2idx) will always ensure that the 
            # new index is 1 + the length of the dictionary so far
            word2idx[word] = len(word2idx)

    # Adding the <UNK> token to the dictionary; This token will be used 
    # on new texts that were not seen during training.
    # It will have the last index. 
    word2idx['<UNK>'] = len(word2idx)

    if shift_for_padding:
        # Adding 1 to all the indexes; 
        # The 0 index will be reserved for padding
        word2idx = {k: v + 1 for k, v in word2idx.items()}

    # Reversing the above dictionary and creating the idx2word dictionary
    idx2word = {idx: word for word, idx in word2idx.items()}

    # Returns the dictionaries
    return word2idx, idx2word

让我们将数据拆分为训练集和测试集，并应用我们创建的函数：

# Spliting to train test 
train, test = train_test_split(d, test_size=0.2, random_state=42)

# Reseting the indexes 
train.reset_index(drop=True, inplace=True)
test.reset_index(drop=True, inplace=True)

print(f'Train shape: {train.shape}')
print(f'Test shape: {test.shape}')
Train shape: (34410, 4)
Test shape: (8603, 4)

# Joining all the texts into one string
text = ' '.join(train['TEXT'].values)

# Creating the word2idx and idx2word dictionaries
word2idx, idx2word = create_word_index(text, shift_for_padding=True, char_level=True)

# Printing the size of the vocabulary
print(f'The size of the vocabulary is: {len(word2idx)}')
The size of the vocabulary is: 274

我们的数据中有 274 个独特的字符。让我们打印word2idx字典中的前 10 项：

{'I': 1,
 ' ': 2,
 'd': 3,
 'o': 4,
 'w': 5,
 'n': 6,
 'l': 7,
 'a': 8,
 'e': 9,
 'G': 10
}

让我们将文本转换为序列：

# For each row in the train and test set, we will create a list of integers
# that will represent the words in the text
train['text_int'] = train['TEXT'].apply(lambda x: [word2idx.get(word, word2idx['<UNK>']) for word in list(x)])
test['text_int'] = test['TEXT'].apply(lambda x: [word2idx.get(word, word2idx['<UNK>']) for word in list(x)])

# Calculating the length of sequences in the train set 
train['seq_len'] = train['text_int'].apply(lambda x: len(x))

# Describing the length of the sequences
train['seq_len'].describe()
count    34410.000000
mean       103.600262
std         79.972798
min          1.000000
25%         41.000000
50%         83.000000
75%        148.000000
max        727.000000

回顾一下，按词级别拆分文本导致序列的平均长度为**~22** 个令牌。现在，我们有长度为**~103** 个令牌的序列。标准差非常高，因此我们将使用最大序列长度为200 进行填充。

def pad_sequences(x: list, pad_length: int) -> list:
    """
    Function that pads a given list of integers to a given length

    Args:
        x (list): list of integers to pad
        pad_length (int): length to pad

    Returns:
        list: padded list of integers
    """
    # Getting the length of the list
    len_x = len(x)

    # Checking if the length of the list is less than the pad_length
    if len_x < pad_length: 
        # Padding the list with 0s
        x = x + [0] * (pad_length - len_x)
    else: 
        # Truncating the list to the desired length
        x = x[:pad_length]

    # Returning the padded list
    return x

# Padding the train and test sequences 
train['text_int'] = train['text_int'].apply(lambda x: pad_sequences(x, 200))
test['text_int'] = test['text_int'].apply(lambda x: pad_sequences(x, 200))

到目前为止，训练集和验证集的数据集如下：

一段数据；作者拍摄的照片

为什么我们应该从普通的 RNN 切换到 LSTM 网络？问题有两个方面：

• 一个简单的 RNN 有所谓的消失梯度问题²或爆炸梯度问题，与网络中for 循环使用的权重相关。

• 网络往往会**“忘记”**长序列数据的初始步骤输入。

为了说明遗忘，请考虑以下示例：

在我们的数据中，平均而言，有 103 个时间步（文本中从左到右的令牌数量）。回顾 RNN 文章中的图表：

展开 n 步 RNN；作者拍摄的照片

我们有相同的权重W，用来乘以 ReLU 层的输出。然后，我们将该信号添加到下一个时间步，依此类推。如果我们为W选择一个相对较小的值（例如 0.5），并且我们有 103 步时间序列数据，从第一次时间步输入到最终输出的影响大致为0.5¹⁰³ * input1，这大致等于零。

第二次输入的信号将是0.5¹⁰² * input2，以此类推。

可以看出，随着时间步的增加，初始时间步的信息在最终输出中的占比越来越少。

为了应对遗忘过去的问题，伟大的思想家们提出了用于时间序列问题的 LSTM 层³。

从内部来看，LSTM 层使用两个激活函数：

• Sigmoid 函数

• Tanh 函数

关于这些函数要记住的关键事实是：

• sigmoid 激活函数接受实数平面上的任何值，并输出一个在 0 和 1 之间的值。

• tanh 函数接受实数平面上的任何值，并输出一个在 -1 和 1 之间的值。

def sigmoid(x: float) -> float: 
    """
    Function that calculates the sigmoid of a given value

    Args:
        x (float): value to calculate the sigmoid

    Returns:
        float: sigmoid of the given value in (0, 1)
    """
    return 1 / (1 + np.exp(-x))

def tanh(x: float) -> float: 
    """
    Function that calculates the tanh of a given value

    Args:
        x (float): value to calculate the tanh

    Returns:
        float: tanh of the given value in (-1, 1)
    """
    return (np.exp(x) - np.exp(-x)) / (np.exp(x) + np.exp(-x))

既然我们已经了解了 sigmoid 和 tanh 激活函数，让我们回到 LSTM 层。

LSTM 层由 2 部分组成（因此得名）：

• 长期记忆块

• 短期记忆块

在每个时间步（或令牌步），LSTM 层输出两个预测：长期预测和短期预测。LSTM 单元的高层次图示可以这样可视化：

展开的简单 LSTM 网络；作者拍摄的图表

在每个时间步骤，LSTM 层输出一个数字，这就是我们所称的短期记忆输出。它通常只是一个标量。此外，长期记忆标量也在 LSTM 层中计算，但它不被输出并传递到序列的第二步。值得注意的是，在每个时间步骤中，短期和长期记忆都会被更新。

现在让我们深入探讨 LSTM 层。LSTM 层的第一部分是所谓的忘记门操作：

忘记门；图表由作者提供

忘记门得名于我们计算希望保留的长期记忆百分比。这是因为 sigmoid 激活函数会输出一个介于 0 和 1 之间的数字，我们将这个数字乘以长期记忆并传递到网络中。

我们可以开始看到在训练时会更新的权重：w1, w2 和 b1。这些权重直接影响保持的长期记忆量。

请注意，在这个步骤中，短期记忆没有调整，而是传递到网络的第二步。

class ForgetGate: 
    """
    Class that implements the forget gate of an LSTM cell
    """
    def __init__(
            self, 
            w1: float = np.random.normal(), 
            w2: float = np.random.normal(),
            b1: float = np.random.normal(),
            long_term_memory: float = np.random.normal(), 
            short_term_memory: float = np.random.normal(), 
            ):
        """
        Constructor of the class

        Args:
            long_term_memory (float): long term memory
            short_term_memory (float): short term memory
            w1 (float): weight 1
            w2 (float): weight 2
            b1 (float): bias term 1
        """
        # Saving the input
        self.long_term_memory = long_term_memory
        self.short_term_memory = short_term_memory
        self.w1 = w1
        self.w2 = w2
        self.b1 = b1

    def forward(self, x: float) -> float: 
        """
        Function that calculates the output of the forget gate

        Args:
            x (float): input to the forget gate

        Returns:
            float: output of the forget gate
        """
        # Calculates the percentage of the long term memory that will be kept
        percentage_to_keep = sigmoid((self.w1 * x  + self.w2 * self.short_term_memory) + self.b1)

        # Updating the long term memory
        self.long_term_memory = self.long_term_memory * percentage_to_keep

        # The output of the forget gate is the new long term memory and the short term memory
        return self.long_term_memory, self.short_term_memory

# Initiating 
forget_gate = ForgetGate()

print(f'Initial long term memory: {forget_gate.long_term_memory}')
print(f'Initial short term memory: {forget_gate.short_term_memory}')

# Calculating the output of the forget gate
lt, st = forget_gate.forward(0.5)

print(f'Long term memory: {lt}')
print(f'Short term memory: {st}')

Initial long term memory: -0.8221542907288696
Initial short term memory: -0.5617438418718841
Long term memory: -0.37335827895028
Short term memory: -0.5617438418718841

接下来在 LSTM 层的是输入门：

输入门；图表由作者提供

输入门仅调整 LSTM 网络的长期记忆部分，但为此，它使用当前输入和当前短期记忆值。

从图表来看，在乘法步骤之前，我们有两个输出：一个来自 sigmoid 激活函数，另一个来自 tanh 激活层。宽泛地说，sigmoid 层输出要记住的记忆百分比（0，1），而 tanh 输出记住的潜在记忆（-1，1）。

然后我们将当前长期记忆（在忘记门中稍作调整）与输入门输出相加。

class InputGate:
    def __init__(
            self, 
            w3: float = np.random.normal(), 
            w4: float = np.random.normal(),
            w5: float = np.random.normal(),
            w6: float = np.random.normal(),
            b2: float = np.random.normal(),
            b3: float = np.random.normal(),
            long_term_memory: float = np.random.normal(), 
            short_term_memory: float = np.random.normal(), 
            ):
        """
        Constructor of the class

        Args:
            long_term_memory (float): long term memory
            short_term_memory (float): short term memory
            w3 (float): weight 3
            w4 (float): weight 4
            w5 (float): weight 5
            w6 (float): weight 6
            b2 (float): bias 2
            b3 (float): bias 3
        """
        # Saving the input
        self.long_term_memory = long_term_memory
        self.short_term_memory = short_term_memory
        self.w3 = w3
        self.w4 = w4
        self.w5 = w5
        self.w6 = w6
        self.b2 = b2
        self.b3 = b3

    def forward(self, x: float) -> float:
        """
        Function that calculates the output of the input gate

        Args:
            x (float): input to the input gate

        Returns:
            float: output of the input gate
        """
        # Calculating the memory signal 
        memory_signal = tanh((self.w3 * x + self.w4 * self.short_term_memory) + self.b2)

        # Calculating the percentage of memory to keep 
        percentage_to_keep = sigmoid((self.w5 * x + self.w6 * self.short_term_memory) + self.b3) 

        # Multiplying the memory signal by the percentage to keep
        memory_signal = memory_signal * percentage_to_keep

        # Updating the long term memory
        self.long_term_memory = self.long_term_memory + memory_signal

        # The output of the input gate is the new long term memory and the short term memory
        return self.long_term_memory, self.short_term_memory

# Creating the input gate object with the forget gates' output 
input_gate = InputGate(long_term_memory=lt, short_term_memory=st)

# Forward propagating 
lt, st = input_gate.forward(0.5)

print(f'Long term memory: {lt}')
print(f'Short term memory: {st}')

Long term memory: -1.028998511766425
Short term memory: -0.5617438418718841

从上面的代码片段可以看出，唯一变化的是长期记忆。

LSTM 层的最后一部分是输出门。输出门是我们将调整短期记忆的步骤：

输出门；图表由作者提供

逻辑与之前门的逻辑非常相似：sigmoid 激活函数计算要保留的记忆百分比，而 tanh 函数计算总体信号。

class OutputGate:
    def __init__(
            self, 
            w7: float = np.random.normal(), 
            w8: float = np.random.normal(),
            b4: float = np.random.normal(),
            long_term_memory: float = np.random.normal(), 
            short_term_memory: float = np.random.normal(), 
            ):
        """
        Constructor of the class

        Args:
            long_term_memory (float): long term memory
            short_term_memory (float): short term memory
            w7 (float): weight 7
            w8 (float): weight 8
            w9 (float): weight 9
            w10 (float): weight 10
            b4 (float): bias 4
            b5 (float): bias 5
        """
        # Saving the input
        self.long_term_memory = long_term_memory
        self.short_term_memory = short_term_memory
        self.w7 = w7
        self.w8 = w8
        self.b4 = b4

    def forward(self, x: float) -> float:
        """
        Function that calculates the output of the output gate

        Args:
            x (float): input to the output gate

        Returns:
            float: output of the output gate
        """
        # Calculating the short term memory signal 
        short_term_memory_signal = tanh(self.long_term_memory)

        # Calculating the percentage of short term memory to keep 
        percentage_to_keep = sigmoid((self.w7 * x + self.w8 * self.short_term_memory) + self.b4) 

        # Multiplying the short term memory signal by the percentage to keep
        short_term_memory_signal = short_term_memory_signal * percentage_to_keep

        # Updating the short term memory
        self.short_term_memory = short_term_memory_signal

        # The output of the output gate is the new long term memory and the short term memory
        return self.long_term_memory, self.short_term_memory

# Creating the output gate object
output_gate = OutputGate(long_term_memory=lt, short_term_memory=st)

# Forward propagating
lt, st = output_gate.forward(0.5)

print(f'Long term memory: {lt}')
print(f'Short term memory: {st}')

Long term memory: -1.028998511766425
Short term memory: -0.7233077589896045

如我们所见，输出门仅调整了短期记忆标量。

LSTM 层；图表由作者提供

上图显示了忘记门、输入门和输出门的合成图⁴。

当我们有一个 x 变量的输入序列时，使用 LSTM 层的内部循环是这样的：

1. 随机初始化短期和长期记忆。

2. 对于每个x1到xn：

2.1 通过 LSTM 层向前传播。

2.2 输出短期记忆

将长期和短期记忆保存到层中。

让我们将每个门封装到一个类中，并创建一个 Python 示例。

# Redefining the forget, input and output gates as functions 
def forget_gate(x: float, w1: float, w2: float, b1: float, long_term_memory: float, short_term_memory: float) -> Tuple[float, float]:
    """
    Function that calculates the output of the forget gate

    Args:
        x (float): input to the forget gate
        w1 (float): weight 1
        w2 (float): weight 2
        b1 (float): bias 1
        long_term_memory (float): long term memory
        short_term_memory (float): short term memory

    Returns:
        Tuple[float, float]: output of the forget gate
    """
    # Calculates the percentage of the long term memory that will be kept
    percentage_to_keep = sigmoid((w1 * x  + w2 * short_term_memory) + b1)

    # Updating the long term memory
    long_term_memory = long_term_memory * percentage_to_keep

    # The output of the forget gate is the new long term memory and the short term memory
    return long_term_memory, short_term_memory

def input_gate(x: float, w3: float, w4: float, w5: float, w6: float, b2: float, b3: float, long_term_memory: float, short_term_memory: float) -> Tuple[float, float]:
    """
    Function that calculates the output of the input gate

    Args:
        x (float): input to the input gate
        w3 (float): weight 3
        w4 (float): weight 4
        w5 (float): weight 5
        w6 (float): weight 6
        b2 (float): bias 2
        b3 (float): bias 3
        long_term_memory (float): long term memory
        short_term_memory (float): short term memory

    Returns:
        Tuple[float, float]: output of the input gate
    """
    # Calculating the memory signal 
    memory_signal = tanh((w3 * x + w4 * short_term_memory) + b2)

    # Calculating the percentage of memory to keep 
    percentage_to_keep = sigmoid((w5 * x + w6 * short_term_memory) + b3) 

    # Multiplying the memory signal by the percentage to keep
    memory_signal = memory_signal * percentage_to_keep

    # Updating the long term memory
    long_term_memory = long_term_memory + memory_signal

    # The output of the input gate is the new long term memory and the short term memory
    return long_term_memory, short_term_memory

def output_gate(x: float, w7: float, w8: float, b4: float, long_term_memory: float, short_term_memory: float) -> Tuple[float, float]:
    """
    Function that calculates the output of the output gate

    Args:
        x (float): input to the output gate
        w7 (float): weight 7
        w8 (float): weight 8
        b4 (float): bias 4
        long_term_memory (float): long term memory
        short_term_memory (float): short term memory

    Returns:
        Tuple[float, float]: output of the output gate
    """
    # Calculating the short term memory signal 
    short_term_memory_signal = tanh(long_term_memory)

    # Calculating the percentage of short term memory to keep 
    percentage_to_keep = sigmoid((w7 * x + w8 * short_term_memory) + b4) 

    # Multiplying the short term memory signal by the percentage to keep
    short_term_memory_signal = short_term_memory_signal * percentage_to_keep

    # Updating the short term memory
    short_term_memory = short_term_memory_signal

    # The output of the output gate is the new long term memory and the short term memory
    return long_term_memory, short_term_memory 

class simpleLSTM: 
    def __init__(
            self, 
            w1: float = np.random.normal(),
            w2: float = np.random.normal(),
            w3: float = np.random.normal(),
            w4: float = np.random.normal(),
            w5: float = np.random.normal(),
            w6: float = np.random.normal(),
            w7: float = np.random.normal(),
            w8: float = np.random.normal(),
            b1: float = np.random.normal(),
            b2: float = np.random.normal(),
            b3: float = np.random.normal(),
            b4: float = np.random.normal(),
            long_term_memory: float = np.random.normal(),
            short_term_memory: float = np.random.normal(),
            ):
        """
        Constructor of the class

        Args:
            long_term_memory (float): long term memory
            short_term_memory (float): short term memory
            w1 (float): weight 1
            w2 (float): weight 2
            w3 (float): weight 3
            w4 (float): weight 4
            w5 (float): weight 5
            w6 (float): weight 6
            w7 (float): weight 7
            w8 (float): weight 8
            b1 (float): bias 1
            b2 (float): bias 2
            b3 (float): bias 3
            b4 (float): bias 4
        """

        # Saving the input
        self.long_term_memory = long_term_memory
        self.short_term_memory = short_term_memory
        self.w1 = w1
        self.w2 = w2
        self.w3 = w3
        self.w4 = w4
        self.w5 = w5
        self.w6 = w6
        self.w7 = w7
        self.w8 = w8
        self.b1 = b1
        self.b2 = b2
        self.b3 = b3
        self.b4 = b4

    def forward(self, x: float) -> float:
        """
        Function that calculates the output of the simple LSTM cell

        Args:
            x (float): input to the simple LSTM cell

        Returns:
            float: output of the simple LSTM cell
        """
        # Calculating the output of the forget gate
        lt, st = forget_gate(x, self.w1, self.w2, self.b1, self.long_term_memory, self.short_term_memory)

        # Updating the long term memory
        self.long_term_memory = lt

        # Calculating the output of the input gate
        lt, st = input_gate(x, self.w3, self.w4, self.w5, self.w6, self.b2, self.b3, self.long_term_memory, self.short_term_memory)

        # Updating the long term memory
        self.long_term_memory = lt

        # Calculating the output of the output gate
        lt, st = output_gate(x, self.w7, self.w8, self.b4, self.long_term_memory, self.short_term_memory)

        # Updating the short term memory
        self.short_term_memory = st

        # The output of the simple LSTM cell is the new long term memory and the short term memory
        return self.long_term_memory, self.short_term_memory

    def forward_sequence(self, x: list) -> list:
        """
        Function that forward propagates a sequence of inputs through the simple LSTM cell

        Args:
            x (list): sequence of inputs to the simple LSTM cell

        Returns:
            list: sequence of outputs of the simple LSTM cell
        """
        # Creating a list to store the outputs
        outputs = []

        # Forward propagating each input
        for input in x: 
            # Forward propagating the input
            _, st = self.forward(input)

            # Appending the output to the list
            outputs.append(st)

        # Returning the list of outputs
        return outputs

# Creating the simple LSTM cell object
simple_lstm = simpleLSTM()

# Creating a sequence of x
x = [0.5, 0.6, 0.7, 0.8, 0.9]

# Forward propagating the sequence
outputs = simple_lstm.forward_sequence(x)

# Rounding 
outputs = [round(output, 2) for output in outputs]

# Printing the outputs
print(f'The outputs of the simple LSTM cell are: {outputs}')
The outputs of the simple LSTM cell are: [0.63, 0.41, 0.33, 0.28, 0.25]

现在我们将所有内容包裹在一个漂亮的 pytorch 示例中，使用 LSTM 层。语法与基本的 RNN 模型非常相似：

# Defining the torch model for sentiment classification 
class SentimentClassifier(torch.nn.Module):
    """
    Class that defines the sentiment classifier model
    """
    def __init__(self, vocab_size, embedding_dim):
        super(SentimentClassifier, self).__init__()

        self.embedding = nn.Embedding(vocab_size + 1, embedding_dim)
        self.lstm = nn.LSTM(input_size=embedding_dim, hidden_size=1, batch_first=True)
        self.fc = nn.Linear(1, 1)  # Output with a single neuron for binary classification
        self.sigmoid = nn.Sigmoid()  # Sigmoid activation

    def forward(self, x):
        x = self.embedding(x)  # Embedding layer
        output, _ = self.lstm(x)  # RNN layer

        # Use the short term memory from the last time step as the representation of the sequence
        x = output[:, -1, :]

        # Fully connected layer with a single neuron
        x = self.fc(x) 

        # Converting to probabilities
        x = self.sigmoid(x)

        # Flattening the output
        x = x.squeeze()

        return x

# Initiating the model 
model = SentimentClassifier(vocab_size=len(word2idx), embedding_dim=16)

# Initiating the criterion and the optimizer
criterion = nn.BCELoss() # Binary cross entropy loss
optimizer = torch.optim.Adam(model.parameters(), lr=0.001)

# Defining the data loader 
from torch.utils.data import Dataset, DataLoader

class TextClassificationDataset(Dataset):
    def __init__(self, data):
        self.data = data

    def __len__(self):
        return len(self.data)

    def __getitem__(self, idx):
        # The x is named as text_int and the y as airline_sentiment
        x = self.data.iloc[idx]['text_int']
        y = self.data.iloc[idx]['SENTIMENT']

        # Converting the x and y to torch tensors
        x = torch.tensor(x)
        y = torch.tensor(y)

        # Converting the y variable to float 
        y = y.float()

        # Returning the x and y
        return x, y

# Creating the train and test loaders
train_loader = DataLoader(TextClassificationDataset(train), batch_size=32, shuffle=True)
test_loader = DataLoader(TextClassificationDataset(test), batch_size=32, shuffle=True)

# Defining the number of epochs
epochs = 100

# Setting the model to train mode
model.train()

# Saving of the loss values
losses = []

# Iterating through epochs
for epoch in range(epochs):
    # Initiating the total loss 
    total_loss = 0

    for batch_idx, (inputs, labels) in enumerate(train_loader):
        # Zero the gradients
        optimizer.zero_grad()  # Zero the gradients
        outputs = model(inputs)  # Forward pass

        loss = criterion(outputs, labels)  # Compute the loss
        loss.backward()  # Backpropagation
        optimizer.step()  # Update the model's parameters

        # Adding the loss to the total loss
        total_loss += loss.item()

    # Calculating the average loss
    avg_loss = total_loss / len(train_loader)

    # Appending the loss to the list containing the losses
    losses.append(avg_loss)

    # Printing the loss every n epochs
    if epoch % 20 == 0:
        print(f'Epoch: {epoch}, Loss: {avg_loss}')

Epoch: 0, Loss: 0.6951859079329055
Epoch: 20, Loss: 0.6478807757224292
Epoch: 40, Loss: 0.6398377026877882
Epoch: 60, Loss: 0.6353290403144067
Epoch: 80, Loss: 0.6312290856884758

# Setting the model to eval model
model.eval()

# List to track the test acc 
total_correct = 0
total_obs = 0

# Iterating over the test set
for batch_idx, (inputs, labels) in enumerate(test_loader):
    outputs = model(inputs)  # Forward pass

    # Getting the number of correct predictions 
    correct = ((outputs > 0.5).float() == labels).float().sum()

    # Getting the total number of predictions
    total = labels.size(0)

    # Updating the total correct and total observations
    total_correct += correct
    total_obs += total

print(f'The test accuracy is: {total_correct / total_obs}')
The test accuracy is: 0.6447750926017761

这篇文章深入探讨了 LSTM 单元的内部工作细节。某些 LSTM 层的实现可能与这里展示的有所不同，但长期和短期记忆的整体部分在绝大多数实现中都存在。

我希望读者现在对 LSTM 层有了更好的理解，并希望他们能够立即将其实施到他们的工作流程中！

特别感谢 StatQuest 提供的精彩讲解视频⁵。

[1]

名称： 推特情感分析

网址： www.kaggle.com/datasets/jp797498e/twitter-entity-sentiment-analysis

数据集许可证： creativecommons.org/publicdomain/zero/1.0/

[2]

名称： 梯度消失问题

网址： k21academy.com/datascience-blog/machine-learning/recurrent-neural-networks/#:~:text=Two%20Issues%20of%20Standard%20RNNs&text=RNNs%20suffer%20from%20the%20matter,of%20long%20data%20sequences%20difficult

[3]

名称： 长短期记忆

网址： www.bioinf.jku.at/publications/older/2604.pdf

年份： 1997

[4]

名称： 理解 LSTM 网络

网址： colah.github.io/posts/2015-08-Understanding-LSTMs/

年份： 2015

[5]

名称： 长短期记忆（LSTM），清晰解释

网址： www.youtube.com/watch?v=YCzL96nL7j0

年份： 2022

掌握模型可解释性：对部分依赖图的全面分析

原文：towardsdatascience.com/mastering-model-interpretability-a-comprehensive-look-at-partial-dependence-plots-d5e78e86051a?source=collection_archive---------14-----------------------#2023-07-07

开始你在可解释 AI 世界中的旅程。

Tiago Toledo Jr.

·

关注 发表在 Towards Data Science ·7 分钟阅读·2023 年 7 月 7 日

–

图片由 David Pupăză 提供，来源于 Unsplash

了解如何解释你的模型对于确定其是否表现异常至关重要。你对模型了解得越多，当它投入生产时，你就越不容易被其行为所惊讶。

此外，你对模型的领域知识越多，你就越能更好地向业务部门推销它。最糟糕的情况是他们意识到你实际上并不确定你在向他们销售什么。

我从未开发过一个不需要解释如何根据输入变量进行预测的模型。至少，向业务说明哪些特征有正面或负面影响是必不可少的。

你可以使用的一个工具来理解你的模型如何工作的就是部分依赖图（PDP），我们将在这篇文章中探讨它。

PDP 是什么

PDP 是一种全局可解释性方法，专注于向你展示模型的特征值如何与模型的输出相关。

这不是理解数据的方法，它只为你的模型生成见解，因此无法从中推断出目标与特征之间的因果关系。然而，它可以让你对模型进行因果推断。

这是因为该方法探测你的模型，因此你可以准确看到当特征变量改变时模型的表现。

它是如何工作的

首先，PDP 允许我们一次只研究一个或两个特征。在这篇文章中，我们将专注于单个特征分析的情况。

在模型训练后，我们生成一个探测数据集。该数据集按照以下算法创建：

• 我们选择我们感兴趣的特征的每个唯一值

• 对于每个唯一值，我们制作一份整个数据集的副本，将特征值设置为该唯一值

• 然后，我们使用模型对这个新数据集进行预测

• 最后，我们对每个唯一值的模型预测进行平均

让我们举个例子。假设我们有以下数据集：

现在，如果我们想将 PDP 应用于特征 0，我们将为特征的每个唯一值重复数据集，例如：

然后，在应用我们的模型后，我们将得到类似这样的结果：

然后，我们计算每个值的平均输出，最终得到以下数据集：

然后只是将这些数据用折线图绘制出来。

对于回归问题，计算每个特征值的平均输出是简单的。对于分类方法，我们可以使用每个类别的预测概率，然后对这些值进行平均。在这种情况下，我们将为数据集中的每个特征和类别对生成一个 PDP。

数学解释

PDP 的解释是我们在边际化一个或两个特征，以评估它们对模型预测输出的边际效应。这由以下公式给出：

其中 $f$ 是机器学习模型，$x_S$ 是我们感兴趣分析的特征集合，$x_C$ 是我们将要平均的其他特征集合。上述函数可以使用以下近似计算：

PDP 的问题

PDP 有一些我们必须注意的限制。首先，由于我们在每个特征值上平均输出，我们最终会得到一个覆盖数据集中每个值的图，即使该值仅出现一次。

因此，你可能会看到数据集中一些非常少见区域的行为，这可能不代表如果这些值更频繁时会发生的情况。因此，在查看 PDP 时，始终查看特征的分布是有帮助的，以了解哪些值更可能发生。

另一个问题发生在特征值相互抵消的情况下。例如，如果你的特征具有以下分布：

计算该特征的 PDP 时，我们将得到类似于以下的结果：

注意特征的影响绝非零，但其平均值为零。这可能会误导你认为特征是无用的，而实际上并非如此。

这种方法的另一个问题是，当我们分析的特征与我们正在平均的特征相关时。如果特征之间存在相关性，如果我们强制数据集中的每个值都具有感兴趣特征的每个值，我们将创建不现实的点。

想象一个包含降雨量和天空中云量的数据集。当我们对降雨量进行平均时，我们会得到一些说天空中有降雨却没有云的点，这是一个不可行的点。

解读 PDP

让我们看看如何分析部分依赖图。请查看下面的图像：

在 x 轴上，我们有特征 0 的值，在 y 轴上，我们有模型对每个特征值的平均输出。注意，对于小于-0.10 的值，模型输出的目标预测非常低，此后预测值上升，然后在特征值超过 0.09 之前，预测值围绕 150 波动，在特征值超过 0.09 时，预测值开始急剧上升。

因此，我们可以说特征与目标预测之间存在正相关关系，但这种相关性不是线性的。

ICE 图

ICE 图尝试解决特征值相互抵消的问题。基本上，在 ICE 图中，我们绘制模型对每个值的每个单独预测，而不仅仅是其平均值。

在 Python 中实现 PDP

让我们在 Python 中实现 PDP。为此，我们首先要导入所需的库：

import numpy as np
import matplotlib.pyplot as plt
from tqdm import tqdm
from sklearn.datasets import load_diabetes
from sklearn.ensemble import RandomForestRegressor

我们将使用来自 sklearn 的糖尿病数据集。tqdm 库将用于为我们的循环生成进度条。

现在，我们将加载数据集并拟合一个随机森林回归器：

X, y = load_diabetes(return_X_y=True)
rf = RandomForestRegressor().fit(X, y)

现在，对于我们数据集中的每个特征，我们将计算在该特征固定值的情况下模型对数据集的平均预测：

features = [0, 1, 2, 3, 4, 5, 6, 7, 8, 9]

features_averages = {}
for feature in tqdm(features):
    features_averages[feature] = ([], [])

  # For each unique value in the feature
    for feature_val in np.unique(X[:, feature]):
        features_averages[feature][0].append(feature_val)

    # We remove the feature from the dataset
        aux_X = np.delete(X, feature, axis=1)
    # We add the feature value for every row of the dataset
        aux_X = np.hstack((aux_X, np.array([feature_val for i in range(aux_X.shape[0])])[:, None]))

    # We calculate the average prediction
        features_averages[feature][1].append(np.mean(rf.predict(aux_X)))

现在，我们为每个特征绘制 PDP 图：

for feature in features_averages:
    plt.figure(figsize=(5,5))
    values = features_averages[feature][0]
    predictions = features_averages[feature][1]

    plt.plot(values, predictions)
    plt.xlabel(f'Feature: {feature}')
    plt.ylabel('Target')

例如，特征 3 的图示为：

结论

现在你有了另一个工具，可以用来提升你的工作效果，并帮助业务部门理解你展示的那个黑箱模型发生了什么。

但不要让理论消失。拿一个你正在开发的模型，将 PDP 可视化应用于它。了解模型在做什么，并在假设上更为精准。

另外，这不是唯一的解释性方法。实际上，我们还有其他方法在处理相关特征时效果更佳。请关注我的下一篇文章，其中会涵盖这些方法。

参考文献

ethen8181.github.io/machine-learning/model_selection/partial_dependence/partial_dependence.html

scikit-learn.org/stable/modules/partial_dependence.html

christophm.github.io/interpretable-ml-book/pdp.html

精通模块化编程：如何提升你的 Python 技能

原文：towardsdatascience.com/mastering-modular-programming-how-to-take-your-python-skills-to-the-next-level-ba14339e8429

编写模块化 Python 代码的最佳实践

Federico Trotta

·发表于Towards Data Science ·阅读时间 9 分钟·2023 年 4 月 10 日

–

图片由全尚烨提供，来源于Pixabay

我最近写了一篇关于Python 类的文章，其中提到类在 Python 中经常作为模块使用。

在这篇文章中，我将描述什么是模块，如何使用模块，以及为什么你应该在 Python 中使用模块。

在这篇文章的最后，如果你从未使用过模块，你将改变看法。你会改变你的编程方式，并在每次可能的情况下使用模块。我对此深信不疑。

我能理解你的感受：当我发现了模块及其强大功能时，我开始在每次可能的情况下使用它们。

我和你之间唯一的区别是我没有找到任何完整而有价值的资源来讲解你需要了解的所有关于 Python 模块的知识：这就是我创建这篇文章的原因。

我发现了一些零碎的资源。感谢这些资源和高级开发者的建议，经过几天的学习，我对 Python 模块有了全面的了解。

在阅读完这篇文章后，我相信你也会清楚地理解这个主题，并不需要其他资料去了解它。因此，如果你是 Python 新手，或者你是一个“摸索”中的 Python 爱好者但从未使用过模块，那么这篇文章绝对适合你。

为什么你应该将 Python 代码模块化

首先，让我们从一些你应该在 Python 中使用模块的理由开始讨论。我找到了 6 个理由：

1. 更好的代码组织。模块允许开发者将代码组织成可重用的单元；这些单元可以被导入并在其他程序中使用。这意味着你不再需要把你的 Python 应用程序看作是一个写满所有代码的整体程序。你会有一个主文件，在其中加载和使用模块。我们将在接下来的段落中以实际的方式看到这个理念。

2. 代码重用性。这是模块的一个优点：你创建一个独立的程序，可以在其他程序中导入和重用。可能需要稍微修改一下；但通常，你可以直接使用。这帮助我们节省了大量时间和编写代码的量：编写类似的代码发生的频率比你想象的要高。

3. 命名空间管理。如果你曾经用任何语言编码，你肯定遇到过需要命名变量的情况。随着代码在应用程序中增长，你可能需要写出具有相似名称的变量。模块化帮助我们修复名称，以便我们不必每次都重新发明轮子，因为我们将使用在模块中创建的名称，即使在主文件中也是如此。

4. 协作。当在公司工作或与朋友合作进行分布式项目（特别是使用 Git）时，模块化允许我们独立工作在模块上。这样，我们可以避免在同一文件中出现重叠问题。

5. 更容易调试。如果你调试过一个 200 行的应用程序，你就知道这有多痛苦。找到一个有几十行代码的错误可能需要几天时间。模块化避免了所有这些痛苦，因为模块彼此独立。它们代表了小程序，因此我们可以比调试一个整体应用程序更容易地调试它们。

6. 性能优化。模块可以优化你的机器性能，因为我们可以只导入我们在特定应用程序中需要的代码。

Python 中的模块是什么？

正如[1]所说：“在 Python 中，模块指的是包和库，但也包括任何可以与其他代码分离且可以独立工作的代码片段”。

基本上，每个 Python 程序都是一个模块，因为它可以独立工作。但这并不意味着我们可以导入所有模块。或者，更好地说：这并不意味着将一个模块导入到另一个文件中总是有意义的。

例如，假设你想创建一个 Python 文件，检查是否存在三个文件夹。如果不存在，它将创建这些文件夹。它可能类似于这样：

import os

# Define folders name
folders = ["fold_1", "fold_2", "fold_3"]

# Check if folders exist. If not, create folders
for folder in folders:
  if not os.path.isdir(folder):
    os.makedirs(folder)

这个文件可以保存为folders.py并运行。它可以独立工作，因此是一个模块。但让我问你一个问题：将它导入到另一个 Python 文件中是否有意义？

停下来思考一下，并回答这个问题。

好吧，答案是否定的，原因很简单：它太具体了。我的意思是：如果你需要将这三个文件夹重新创建到另一个文件夹中，你可以将文件移动到新文件夹中并运行它。

因此，当我们导入模块时，我们想要的是通用性。这就是为什么我们创建类和函数作为模块：因为它们是通用的。如果你想在其他文件中重用上述模块，你必须使它通用。例如，我们可以创建一个如下的函数：

import os

def create_folders(fold_1, fold_2, fold_3):
  # Define folders name
  folders = [fold_1, fold_2, fold_3]

  # Check if folders exist. If not, create folders
  for folder in folders:
    if not os.path.isdir(folder):
      os.makedirs(folder)

因此，这个函数可以用来检查在一个特定的文件夹中是否存在三个通用的文件夹。但这次，名称可以是任何的，因为它们作为函数的参数传递。如果它们不存在，则会被创建。所以，例如，我们可以这样调用它：

# Invoke the function
create_folders("audio", "image", "documents")

然后我们的函数会检查目录中是否存在“audio”、“image”和“documents”作为文件夹。如果不存在，它会创建这些文件夹。

这段代码适合被导入到另一个 Python 文件中，因为它是通用的：无论何时我们需要检查三个文件夹，我们都可以使用它。我们只需声明这三个文件夹的名称即可。

如何在 Python 中使用模块

从示意图来看，模块是这样工作的：

模块工作示意图。作者图片。

我们有一个通常可以称为 main.py 的“主” Python 文件。在该 Python 文件中，我们导入了两个模块。

也许你没有注意到，但在上述代码中，我们导入了 os 模块，它是一个帮助创建、管理和删除目录的 Python 模块。这意味着在管理目录时，我们不需要创建自定义的函数（或类）：我们可以使用 os 模块并将其应用于我们的具体情况。例如，我们使用了 os.makedirs 来创建文件夹。

但除了已知的模块和包，我们如何在日常编程活动中使用 Python 的模块化呢？好吧，如果你从未使用过模块化的力量，你需要“转换”你的思维；但相信我，这值得一试！

如果你不使用模块，你的 Python 应用程序看起来会像这样：

一个单体 Python 应用程序。作者图片。

所以，在一个文件中你：

• 创建函数和类。

• 调用使用变量或其他你可能需要的内容创建的函数和类。

相反，如果你使用模块，你的 Python 应用程序会变成这样：

一个模块化的 Python 应用程序。作者图片。

因此，这样，你在单独的 Python 文件中创建两个函数和两个类——这些都是模块。然后，你将它们导入到主文件中。

此时，你的文件夹组织应该是这样的：

project_name

├── main.py
│   └── packages
│       └── __init__.py
│       └── module_1.py
│       └── module_2.py
│       └── module_3.py
└──     └── module_4.py

当然，你可以创建你所需的任何数量的 packages 子文件夹。

重要的是，每个子文件夹中必须有一个 __init__.py 文件。

让我们看看它是什么，以及如何使用它。

Python 的一个实际示例

那么，让我们做一个实际的 Python 示例。我们创建一个名为 my_project 的项目，其中包含：

• 一个 main.py 文件是我们的主要文件。

• 一个名为 operations 的子文件夹，其中我们将有三个 Python 文件：__init__.py、addition.py、subtraction.py。

所以，这是我们的结构：

my_project

├── main.py
│   └── operations
│       └── __init__.py
│       └── addition.py
└──     └── subtraction.py

我们想创建两个简单的模块，用于计算两个整数的和与差。如果输入不是整数，程序将返回错误。

所以，这两个模块可以是这样的：

def sum(x: int, y: int) -> int:
        """ this function sums two integers.

        Args:
            param 1: integer                         
            param 2: integer

        Returns:
            an integer that is the sum of the two arguments
        """
        if type(x) and type(y) == int:
            print(f"the sum of these two numbers is: {x+y}")
        else:
            print("arguments are not integers!")

if __name__ == "__main__":
     sum()

和这些：

def diff(x: int, y: int) -> int:
        """ this function subtracts two integers.

        Args:
            param 1: integer                         
            param 2: integer

        Returns:
            an integer that is the difference of the two arguments
        """
        if type(x) and type(y) == int:
            print(f"the difference of these two numbers is: {x-y}")
        else:
            print("arguments are not integers!")

if __name__ == "__main__":
     diff()

**NOTE:**

If you are not familiar with this way of writing Python with type hints,
and if you don't know how ' if __name__ == "__main__"' work, you
have to read my article [here](https://medium.com/towards-data-science/python-classes-made-easy-the-definitive-guide-to-object-oriented-programming-881ed609fb6) and eveything will be clarified.

现在，我们希望将这两个模块导入到主文件中。为此，我们必须在 __init__.py 中写下以下内容：

__all__ = ["sum", "diff"]

所以，__init__.py 是一个在处理模块时必需的 Python 文件，因为它以某种方式告诉主文件在 operations 文件夹中有一些包需要导入。因此，在 __init__.py 中，我们必须声明所有包中的函数，如上所示。

现在，main.py 可以像这样编程：

from operations.addition import sum
from operations.subtraction import diff

# Sum two numbers
a = sum(1, 2)
# Subtract two numbers
b = diff(3,2)

# Print results
print(a)
print(b)

关于主文件，有几个考虑因素：

1. 由于这些模块每个只包含一个函数，我们可以写 from operations.addition import *。不过，从模块中仅导入我们使用的函数和类是一种良好的实践。

2. 如我们所见，主文件非常简洁。它仅包含导入、使用导入方法声明的变量和打印结果。

最后的提示

我从前辈那里得到的建议是为每个模块创建一个类，因为这样更容易维护。然而，这并不是必须的：这取决于你编写的代码。

在上面的示例中，我创建了两个包含一个函数的单独模块。一种替代方案是创建一个名为 Operations 的类，并将这两个函数作为其方法来编写。

这将有助于我们在一个文件中创建“类似函数”，并在主文件中进行一次导入。

所以，不要总是把建议当作必须的：接受它们，但要根据你需要创建的东西来思考，同时考虑你需要的优化。

结论

现在，告诉我：你对这种方法有多兴奋？好吧，当我理清了对模块的理解后，我开始将它们应用到我的日常工作中，相信我，情况只会变得更好。

希望我已经为你澄清了这个话题；如果没有，请：在评论中告诉我。我可以创作其他文章作为补充。

• 订阅 我的 Substack 新闻通讯 以获取更多 Python 相关内容。

• 通过 我的推荐链接加入 Medium：以 5$/月（无额外费用）解锁 Medium 上的所有内容。

• 在这里找到/联系我： https://bio.link/federicotrotta

• 觉得有用吗？ 请给我买杯 Ko-fi。

更多来自我的：

[## Python 类简易教程：面向对象编程的终极指南

通过这本全面的类参考提升你的 Python 技能

Python 循环：如何在 Python 中迭代的完整指南

利用 Python 中循环的强大功能

Python 循环：如何在 Python 中迭代的完整指南

视频制作：

[1] 如果 name == “main” 对于 Python 开发者 (视频)

精通蒙特卡洛：如何通过模拟提升机器学习模型

原文：towardsdatascience.com/mastering-monte-carlo-how-to-simulate-your-way-to-better-machine-learning-models-6b57ec4e5514

通过模拟技术提升预测算法的概率方法应用

Sydney Nye

·发表于 Towards Data Science ·26 分钟阅读·2023 年 8 月 2 日

–

“在蒙特卡洛的轮盘赌桌上”由爱德华·蒙克（1892 年）

一位科学家玩扑克牌如何永远改变了统计学游戏

在 1945 年的动荡岁月中，当世界正经历第二次世界大战的最后阶段时，一场扑克牌游戏悄然引发了计算领域的一次突破。这不仅仅是普通的游戏，而是导致蒙特卡洛方法诞生的关键。玩家？正是科学家斯坦尼斯瓦夫·乌拉姆，他当时还深陷于曼哈顿计划中。乌拉姆在从疾病中恢复的过程中，沉迷于扑克牌游戏。游戏的复杂概率引起了他的兴趣，他意识到反复模拟游戏可以提供这些概率的良好近似。这是一个灵光一现的时刻，类似于牛顿的苹果，但换成了扑克牌。随后，乌拉姆与他的同事约翰·冯·诺依曼讨论了这些想法，两人共同形式化了蒙特卡洛方法，命名来源于著名的摩纳哥蒙特卡洛赌场（在上面的爱德华·蒙克的著名画作中描绘），在那里赌注高昂，运气主宰——这正如方法本身。

快进到今天，蒙特卡罗方法已成为机器学习领域中的一张王牌，包括在强化学习、贝叶斯滤波和复杂模型优化中的应用(4)。其稳健性和多功能性确保了其持续的相关性，距其诞生已超过七十年。从乌拉姆的纸牌游戏到今天复杂的 AI 应用，蒙特卡罗方法依然是应对复杂系统中模拟和逼近力量的见证。

用概率模拟来玩转你的牌

在数据科学和机器学习的复杂世界中，蒙特卡罗模拟类似于精心计算的赌注。这种统计技术使我们能够在不确定性面前进行战略性投注，为复杂的确定性问题提供概率性的解决方案。在这篇文章中，我们将揭示蒙特卡罗模拟的神秘面纱，并探讨其在统计学和机器学习中的强大应用。

我们将深入探讨蒙特卡罗模拟背后的理论，阐明使这一技术成为强大问题解决工具的原理。我们将通过一些 Python 中的实际应用，演示如何在实践中实现蒙特卡罗模拟。

接下来，我们将探讨如何利用蒙特卡罗模拟来优化机器学习模型。我们将专注于通常具有挑战性的超参数调优任务，提供一个实用的工具包，以帮助应对这一复杂的领域。

所以，下注吧，开始吧！

理解蒙特卡罗模拟

蒙特卡罗模拟对于数学家和数据科学家来说是一项宝贵的技术。这些模拟提供了一种在广泛而复杂的可能性中导航的方法，制定有根据的假设，并逐步完善选择，直到出现最合适的解决方案。

其工作原理如下：我们生成大量随机场景，遵循某个预定义的过程，然后审视这些场景，以估计各种结果的概率。这里有一个类比来帮助理解：将每个场景视为流行的哈斯布罗棋盘游戏“Clue”中的一次回合。对于不熟悉的人来说，“Clue”是一个侦探风格的游戏，玩家在一座豪宅中移动，收集证据以推断犯罪的细节——谁、什么和哪里。每回合或提问都会排除潜在的答案，使玩家更接近揭示真实的犯罪场景。同样，蒙特卡罗研究中的每次模拟提供的见解使我们更接近解决复杂问题的方案。

在机器学习的领域，这些“场景”可以代表不同的模型配置、各种超参数集、将数据集分割为训练集和测试集的不同方式以及其他许多应用。通过评估这些场景的结果，我们可以获得对机器学习算法行为的宝贵洞察，从而做出有关其优化的明智决策。

飞镖游戏

要理解蒙特卡洛模拟，想象你在玩飞镖游戏。但不是瞄准特定目标，而是蒙上眼睛随机地将飞镖扔向一个大的正方形飞镖盘。这个正方形里有一个圆形目标。你的目标是估算圆周率π，即圆的周长与直径的比值。

听起来不可能，对吧？但这里有个窍门：圆的面积与正方形的面积之比是π/4。因此，如果你投掷大量的飞镖，落在圆内的飞镖数与总飞镖数的比率应大致为π/4。将这个比率乘以 4，你就能估算出π！

随机猜测 vs. 蒙特卡洛

为了展示蒙特卡洛模拟的威力，让我们将它与一种更简单的方法进行比较，也许是所有方法中最简单的：随机猜测。

当你运行下面的代码进行随机和蒙特卡洛两种情况时，每次都会得到一组不同的预测。这是意料之中的，因为飞镖是随机投掷的。这是蒙特卡洛模拟的一个关键特征：它们本质上是随机的。然而，尽管有这种随机性，它们在正确使用时可以提供非常准确的估算。因此，虽然你的结果不会完全像我的，但它们会讲述相同的故事。

在第一组可视化图（图 1a 到 图 1f）中，我们对π的值进行了一系列随机猜测，每次生成一个基于猜测值的圆。让我们给这个“随机性”一个正确的推动，假设虽然我们不能记住π的确切值，但我们知道它高于 2 且低于 4。从结果图中可以看到，圆的大小根据猜测值变化很大，这显示了这种方法的不准确性（这也不应让人感到惊讶）。每个图中的绿色圆圈代表单位圆，即我们试图估算的“真实”圆。蓝色圆圈则基于我们的随机猜测。

#Random Guessing of Pi

# Before running this code, make sure you have the necessary packages installed.
# You can install them using "pip install" on your terminal or "conda install" if using a conda environment

# Import necessary libraries
import random
import plotly.graph_objects as go
import numpy as np

# Number of guesses to make. Adjust this for more guesses and subsequent plots
num_guesses = 6

# Generate the coordinates for the unit circle
# We use np.linspace to generate evenly spaced numbers over the range from 0 to 2*pi.
# These represent the angles in the unit circle.
theta = np.linspace(0, 2*np.pi, 100)

# The x and y coordinates of the unit circle are then the cosine and sine of these angles, respectively.
unit_circle_x = np.cos(theta)
unit_circle_y = np.sin(theta)

# We'll make a number of guesses for the value of pi
for i in range(num_guesses):
    # Make a random guess for pi between 2 and 4
    pi_guess = random.uniform(2, 4)

    # Generate the coordinates for the circle based on the guess
    # The radius of the circle is the guessed value of pi divided by 4.
    radius = pi_guess / 4

    # The x and y coordinates of the circle are then the radius times the cosine and sine of the angles, respectively.
    circle_x = radius * np.cos(theta)
    circle_y = radius * np.sin(theta)

    # Create a scatter plot of the circle
    fig = go.Figure()

    # Add the circle to the plot
    # We use a Scatter trace with mode 'lines' to draw the circle.
    fig.add_trace(go.Scatter(
        x = circle_x,
        y = circle_y,
        mode='lines',
        line=dict(
            color='blue',
            width=3
        ),
        name='Estimated Circle'
    ))

    # Add the unit circle to the plot
    fig.add_trace(go.Scatter(
        x = unit_circle_x,
        y = unit_circle_y,
        mode='lines',
        line=dict(
            color='green',
            width=3
        ),
        name='Unit Circle'
    ))

    # Update the layout of the plot
    # We set the title to include the guessed value of pi, and adjust the size and axis ranges to properly display the circles.
    fig.update_layout(
        title=f"Fig1{chr(97 + i)}: Randomly Guessing Pi: {pi_guess}",
        width=600,
        height=600,
        xaxis=dict(
            constrain="domain",
            range=[-1, 1]
        ),
        yaxis=dict(
            scaleanchor="x",
            scaleratio=1,
            range=[-1, 1]
        )
    )

    # Display the plots
    fig.show()

图 1a-1f：π的随机估算

你可能会注意到一个奇怪的现象：在随机猜测方法中，有时接近真实π值的猜测会导致一个更远离单位圆的圆。这种明显的矛盾产生是因为我们关注的是圆的周长，而不是半径或面积。两个圆之间的视觉差距表示的是基于猜测的圆周长估算误差，而不是整个圆的误差。

在第二组可视化（图 2a至图 2f）中，我们使用蒙特卡罗方法来估算π值。我们不是进行随机猜测，而是向一个正方形上扔大量的飞镖，并计算落在正方形内切圆内的数量。由此得到的π值估算更为准确，如图中所示，圆的大小更接近实际单位圆。绿色点表示落在单位圆内的飞镖，红色点表示落在单位圆外的飞镖。

#Monte Carlo Estimation of Pi

# Import our required libraries
import random
import math
import plotly.graph_objects as go
import plotly.io as pio
import numpy as np

# We'll simulate throwing darts at a dartboard to estimate pi. Let's throw 10,000 darts.
num_darts = 10000

# To keep track of darts that land in the circle.
darts_in_circle = 0

# We'll store the coordinates of darts that fall inside and outside the circle.
x_coords_in, y_coords_in, x_coords_out, y_coords_out = [], [], [], []

# Let's generate 6 figures throughout the simulation. Therefore, we will create a new figure every 1,666 darts (10,000 divided by 6).
num_figures = 6
darts_per_figure = num_darts // num_figures

# Create a unit circle to compare our estimates against. Here, we use polar coordinates and convert to Cartesian.
theta = np.linspace(0, 2*np.pi, 100)
unit_circle_x = np.cos(theta)
unit_circle_y = np.sin(theta)

# We start throwing the darts (simulating random points inside a 1x1 square and checking if they fall inside a quarter circle).
for i in range(num_darts):
    # Generate random x, y coordinates between -1 and 1.
    x, y = random.uniform(-1, 1), random.uniform(-1, 1)

    # If a dart (point) is closer to the origin (0,0) than the distance of 1, it's inside the circle.
    if math.sqrt(x**2 + y**2) <= 1:
        darts_in_circle += 1
        x_coords_in.append(x)
        y_coords_in.append(y)
    else:
        x_coords_out.append(x)
        y_coords_out.append(y)

    # After every 1,666 darts, let's see how our estimate looks compared to the real unit circle.
    if (i + 1) % darts_per_figure == 0:
        # We estimate pi by seeing the proportion of darts that landed in the circle (out of the total number of darts).
        pi_estimate = 4 * darts_in_circle / (i + 1)

        # Now we create a circle from our estimate to compare visually with the unit circle.
        estimated_circle_radius = pi_estimate / 4
        estimated_circle_x = estimated_circle_radius * np.cos(theta)
        estimated_circle_y = estimated_circle_radius * np.sin(theta)

        # Plot the results using Plotly.
        fig = go.Figure()

        # Add the darts that landed inside and outside the circle to the plot.
        fig.add_trace(go.Scattergl(x=x_coords_in, y=y_coords_in, mode='markers', name='Darts Inside Circle', marker=dict(color='green', size=4, opacity=0.8)))
        fig.add_trace(go.Scattergl(x=x_coords_out, y=y_coords_out, mode='markers', name='Darts Outside Circle', marker=dict(color='red', size=4, opacity=0.8)))

        # Add the real unit circle and our estimated circle to the plot.
        fig.add_trace(go.Scatter(x=unit_circle_x, y=unit_circle_y, mode='lines', name='Unit Circle', line=dict(color='green', width=3)))
        fig.add_trace(go.Scatter(x=estimated_circle_x, y=estimated_circle_y, mode='lines', name='Estimated Circle', line=dict(color='blue', width=3)))

        # Customize the plot layout.
        fig.update_layout(title=f"Figure {chr(97 + (i + 1) // darts_per_figure - 1)}: Thrown Darts: {(i + 1)}, Estimated Pi: {pi_estimate}", width=600, height=600, xaxis=dict(constrain="domain", range=[-1, 1]), yaxis=dict(scaleanchor="x", scaleratio=1, range=[-1, 1]), legend=dict(yanchor="top", y=0.99, xanchor="left", x=0.01))

        # Display the plot.
        fig.show()

        # Save the plot as a PNG image file.
        pio.write_image(fig, f"fig2{chr(97 + (i + 1) // darts_per_figure - 1)}.png")

图 2a-2f：蒙特卡罗估算π值

在蒙特卡罗方法中，π的估算是基于落在圆内的“飞镖”与总投掷次数的比例。得到的π值估算用于生成一个圆。如果蒙特卡罗估算不准确，生成的圆将再次出现错误。这个估算圆与单位圆之间的差距宽度，提供了蒙特卡罗估算准确性的指示。

然而，由于蒙特卡罗方法在“飞镖”数量增加时会产生更准确的估算，因此随着更多“飞镖”的投掷，估算的圆应当会趋近于单位圆。因此，虽然两种方法在估算不准确时都显示出差距，但随着“飞镖”数量的增加，这种差距在蒙特卡罗方法中应当会更加一致地减少。

预测π值：概率的力量

蒙特卡罗模拟之所以如此强大，是因为它们能够利用随机性来解决确定性问题。通过生成大量随机场景并分析结果，我们可以估算不同结果的概率，即使对于那些难以通过解析方法解决的复杂问题也是如此。

在估计π的情况下，蒙特卡罗方法允许我们即使随机投掷飞镖，也能做出非常准确的估计。如前所述，投掷的飞镖越多，我们的估计就越准确。这展示了大数法则，这是概率论中的一个基本概念，指出从大量试验中获得的结果的平均值应该接近期望值，并且随着试验次数的增加，会趋近于期望值。让我们通过绘制投掷的飞镖数量与蒙特卡罗估计的π与真实π之间的差异的关系图来检验我们的图 2a-2f中的六个例子是否趋于真实。一般来说，我们的图表（图 2g）应该呈现负趋势。这里是实现这一点的代码：

# Calculate the differences between the real pi and the estimated pi
diff_pi = [abs(estimate - math.pi) for estimate in pi_estimates]

# Create the figure for the number of darts vs difference in pi plot (Figure 2g)
fig2g = go.Figure(data=go.Scatter(x=num_darts_thrown, y=diff_pi, mode='lines'))

# Add title and labels to the plot
fig2g.update_layout(
    title="Fig2g: Darts Thrown vs Difference in Estimated Pi",
    xaxis_title="Number of Darts Thrown",
    yaxis_title="Difference in Pi",
)

# Display the plot
fig2g.show()

# Save the plot as a png
pio.write_image(fig2g, "fig2g.png")

注意，即使只有 6 个例子，整体模式也如预期：投掷的飞镖更多（更多的场景），估计值与真实值之间的差异更小，从而预测更准确。

假设我们投掷 1,000,000 支飞镖，并进行 500 次预测。换句话说，我们将记录在 1,000,000 支飞镖模拟过程中，π的估计值与实际值之间的差异，在 500 个均匀间隔的位置上。与其生成 500 个额外的图形，不如直接跳到我们要确认的内容：是否真的如飞镖数量增加那样，π的预测值与真实π之间的差异变小。我们将使用散点图（图 2h）：

#500 Monte Carlo Scenarios; 1,000,000 darts thrown
import random
import math
import plotly.graph_objects as go
import numpy as np

# Total number of darts to throw (1M)
num_darts = 1000000
darts_in_circle = 0

# Number of scenarios to record (500)
num_scenarios = 500
darts_per_scenario = num_darts // num_scenarios

# Lists to store the data for each scenario
darts_thrown_list = []
pi_diff_list = []

# We'll throw a number of darts
for i in range(num_darts):
    # Generate random x, y coordinates between -1 and 1
    x, y = random.uniform(-1, 1), random.uniform(-1, 1)

    # Check if the dart is inside the circle
    # A dart is inside the circle if the distance from the origin (0,0) is less than or equal to 1
    if math.sqrt(x**2 + y**2) <= 1:
        darts_in_circle += 1

    # If it's time to record a scenario
    if (i + 1) % darts_per_scenario == 0:
        # Estimate pi with Monte Carlo method
        # The estimate is 4 times the number of darts in the circle divided by the total number of darts
        pi_estimate = 4 * darts_in_circle / (i + 1)

        # Record the number of darts thrown and the difference between the estimated and actual values of pi
        darts_thrown_list.append((i + 1) / 1000)  # Dividing by 1000 to display in thousands
        pi_diff_list.append(abs(pi_estimate - math.pi))

# Create a scatter plot of the data
fig = go.Figure(data=go.Scattergl(x=darts_thrown_list, y=pi_diff_list, mode='markers'))

# Update the layout of the plot
fig.update_layout(
    title="Fig2h: Difference between Estimated and Actual Pi vs. Number of Darts Thrown (in thousands)",
    xaxis_title="Number of Darts Thrown (in thousands)",
    yaxis_title="Difference between Estimated and Actual Pi",
)

# Display the plot
fig.show()
# Save the plot as a png
pio.write_image(fig2h, "fig2h.png")

蒙特卡罗模拟与超参数调优：完美组合

你可能会想，“蒙特卡罗是一个有趣的统计工具，但它如何应用于机器学习呢？”简短的回答是：有很多种方式。蒙特卡罗模拟在机器学习中的一个应用就是超参数调优。

超参数是我们（人类）在设置机器学习算法时调整的旋钮和拨盘。它们控制算法行为的方面，这些方面关键的是从数据中无法学习到的。例如，在决策树中，树的最大深度是一个超参数。在神经网络中，学习率和隐藏层的数量是超参数。

选择合适的超参数可以决定一个模型表现差异与表现优异之间的差别。但是我们如何知道选择哪些超参数呢？这就是蒙特卡罗模拟的作用所在。

传统上，机器学习从业者使用网格搜索或随机搜索等方法来调整超参数。这些方法涉及为每个超参数指定一组可能的值，然后为每一个可能的超参数组合训练和评估模型。这可能会非常耗费计算资源和时间，尤其是当有许多超参数需要调整或每个超参数可以取的值范围很大时。

蒙特卡洛模拟提供了一个更高效的替代方案。我们可以根据某种概率分布从超参数空间中随机采样，而不是穷尽所有可能的超参数组合。这使我们能够更有效地探索超参数空间，并更快地找到好的超参数组合。

在下一部分，我们将使用一个真实数据集来演示如何在实践中使用蒙特卡洛模拟进行超参数调整。让我们开始吧！

蒙特卡洛模拟用于超参数调整

实验的核心：心脏病数据集

在机器学习的世界里，数据是驱动我们模型的命脉。为了探索蒙特卡洛模拟在超参数调整中的应用，让我们来看一个与心脏息息相关的数据集——字面上的意思。心脏病数据集（CC BY 4.0）来自 UCI 机器学习库，是一些患者的医学记录，其中有些人患有心脏病。

数据集包含 14 个属性，包括年龄、性别、胸痛类型、静息血压、胆固醇水平、空腹血糖等。目标变量是心脏病的存在，使这是一个二分类任务。由于包含了分类和数值特征，这是一个有趣的数据集，用于演示超参数调整。

首先，让我们查看一下我们的数据集，以了解我们将要处理的内容——总是一个好的开始。

 #Load and view first few rows of dataset

# Import necessary libraries
import pandas as pd
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler, OneHotEncoder
from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import GridSearchCV
from sklearn.metrics import roc_auc_score
import numpy as np
import plotly.graph_objects as go

# Load the dataset
# The dataset is available at the UCI Machine Learning Repository
# It's a dataset about heart disease and includes various patient measurements
url = "https://archive.ics.uci.edu/ml/machine-learning-databases/heart-disease/processed.cleveland.data"

# Define the column names for the dataframe
column_names = ["age", "sex", "cp", "trestbps", "chol", "fbs", "restecg", "thalach", "exang", "oldpeak", "slope", "ca", "thal", "target"]

# Load the dataset into a pandas dataframe
# We specify the column names and also tell pandas to treat '?' as NaN
df = pd.read_csv(url, names=column_names, na_values="?")

# Print the first few rows of the dataframe
# This gives us a quick overview of the data
print(df.head())

这展示了数据集中所有列的前四个值。如果你加载了正确的 csv 文件并且像我一样命名了你的列，你的输出将会类似于图 3。

图 3：数据集中前 4 行的数据显示

设定脉搏：数据预处理

在我们可以使用心脏病数据集进行超参数调整之前，我们需要对数据进行预处理。这涉及几个步骤：

1. 处理缺失值：数据集中有些记录缺失了值。我们需要决定如何处理这些缺失值，可以通过删除记录、填补缺失值或其他方法。

2. 编码分类变量：许多机器学习算法要求输入数据为数值型。我们需要将分类变量转换为数值格式。

3. 规范化数值特征：机器学习算法在数值特征处于类似尺度时通常表现更好。我们将应用规范化来调整这些特征的尺度。

首先处理缺失值。在我们的心脏病数据集中，‘ca’和‘thal’列中有一些缺失值。我们将用各自列的中位数填补这些缺失值。这是一种处理缺失数据的常见策略，因为它不会严重影响数据的分布。

接下来，我们将对分类变量进行编码。在我们的数据集中，‘cp’，‘restecg’，‘slope’，‘ca’和‘thal’列是分类的。我们将使用标签编码将这些分类变量转换为数值型的。标签编码将列中的每个唯一类别分配给不同的整数。

最后，我们将规范化数值特征。规范化调整数值特征的尺度，使它们都落在类似的范围内。这可以帮助提高许多机器学习算法的性能。我们将使用标准化来进行规范化，将数据转化为均值为 0，标准差为 1。

这是执行所有这些预处理步骤的 Python 代码：

# Preprocess

# Import necessary libraries
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import LabelEncoder

# Identify missing values in the dataset
# This will print the number of missing values in each column
print(df.isnull().sum())

# Fill missing values with the median of the column
# The SimpleImputer class from sklearn provides basic strategies for imputing missing values
# We're using the 'median' strategy, which replaces missing values with the median of each column
imputer = SimpleImputer(strategy='median')

# Apply the imputer to the dataframe
# The result is a new dataframe where missing values have been filled in
df_filled = pd.DataFrame(imputer.fit_transform(df), columns=df.columns)

# Print the first few rows of the filled dataframe
# This gives us a quick check to make sure the imputation worked correctly
print(df_filled.head())

# Identify categorical variables in the dataset
# These are variables that contain non-numerical data
categorical_vars = df_filled.select_dtypes(include='object').columns

# Encode categorical variables
# The LabelEncoder class from sklearn converts each unique string into a unique integer
encoder = LabelEncoder()
for var in categorical_vars:
    df_filled[var] = encoder.fit_transform(df_filled[var])

# Normalize numerical features
# The StandardScaler class from sklearn standardizes features by removing the mean and scaling to unit variance
scaler = StandardScaler()

# Apply the scaler to the dataframe
# The result is a new dataframe where numerical features have been normalized
df_normalized = pd.DataFrame(scaler.fit_transform(df_filled), columns=df_filled.columns)

# Print the first few rows of the normalized dataframe
# This gives us a quick check to make sure the normalization worked correctly
print(df_normalized.head())

第一个打印语句显示了原始数据集中每列的缺失值数量。在我们的情况下，‘ca’和‘thal’列中有一些缺失值。

第二个打印语句显示了填补缺失值后的数据集前几行。如前所述，我们使用了每列的中位数来填补缺失值。

第三个打印语句显示了对分类变量进行编码后数据集的前几行。经过这一步骤，我们的数据集中的所有变量都是数值型的。

最后的打印语句显示了在规范化数值特征后数据集的前几行，其中数据的均值为 0，标准差为 1。经过这一步骤，我们数据集中的所有数值特征都在类似的尺度上。检查您的输出是否类似于图 4：

图 4：预处理打印语句输出

运行这段代码后，我们得到了一个已经过预处理的数据集，准备好进行建模。

实现一个基本的机器学习模型

现在我们已经预处理了数据，准备实施一个基本的机器学习模型。这将作为我们的基准模型，之后我们将尝试通过超参数调优来改进它。

我们将使用一个简单的逻辑回归模型来完成这个任务。请注意，尽管它被称为“回归”，但它实际上是处理二分类问题（如我们在心脏病数据集中遇到的）的最受欢迎的算法之一。这是一个线性模型，用于预测正类的概率。

在训练我们的模型后，我们将使用两个常见指标来评估其性能：准确度和 ROC-AUC。准确度是所有预测中正确预测的比例，而 ROC-AUC（接收者操作特征曲线—曲线下面积）衡量真实正率和假正率之间的权衡。

那么这与蒙特卡洛模拟有什么关系呢？好吧，像逻辑回归这样的机器学习模型有几个可以调整的超参数，以提高性能。然而，找到最佳的超参数组合就像在大海捞针。这就是蒙特卡洛模拟的用武之地。通过随机采样不同的超参数组合并评估其性能，我们可以估算出良好超参数的概率分布，并对使用哪些最佳参数做出有根据的猜测，类似于我们在掷飞镖练习中挑选更好π值的方式。

这是实现并评估我们新预处理数据的基本逻辑回归模型的 Python 代码：

# Logistic Regression Model - Baseline

# Import necessary libraries
from sklearn.model_selection import train_test_split
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import accuracy_score, roc_auc_score

# Replace the 'target' column in the normalized DataFrame with the original 'target' column
# This is done because the 'target' column was also normalized, which is not what we want
df_normalized['target'] = df['target']

# Binarize the 'target' column
# This is done because the original 'target' column contains values from 0 to 4
# We want to simplify the problem to a binary classification problem: heart disease or no heart disease
df_normalized['target'] = df_normalized['target'].apply(lambda x: 1 if x > 0 else 0)

# Split the data into training and test sets
# The 'target' column is our label, so we drop it from our features (X)
# We use a test size of 20%, meaning 80% of the data will be used for training and 20% for testing
X = df_normalized.drop('target', axis=1)
y = df_normalized['target']
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# Implement a basic logistic regression model
# Logistic Regression is a simple yet powerful linear model for binary classification problems
model = LogisticRegression()
model.fit(X_train, y_train)

# Make predictions on the test set
# The model has been trained, so we can now use it to make predictions on unseen data
y_pred = model.predict(X_test)

# Evaluate the model
# We use accuracy (the proportion of correct predictions) and ROC-AUC (a measure of how well the model distinguishes between classes) as our metrics
accuracy = accuracy_score(y_test, y_pred)
roc_auc = roc_auc_score(y_test, y_pred)

# Print the performance metrics
# These give us an indication of how well our model is performing
print("Baseline Model " + f'Accuracy: {accuracy}')
print("Baseline Model " + f'ROC-AUC: {roc_auc}')

我们的基本逻辑回归模型具有 0.885 的准确度和 0.884 的 ROC-AUC 分数，为我们提供了一个坚实的基准，以便进一步改进。这些指标表明我们的模型在区分有心脏病和无心脏病的患者方面表现相当出色。让我们看看是否能进一步提高。

使用网格搜索进行超参数调整

在机器学习中，模型的性能通常可以通过调整其超参数来提高。超参数是从数据中未学习到的参数，而是在学习过程开始之前设置的。例如，在逻辑回归中，正则化强度‘C’和惩罚类型‘l1’或‘l2’都是超参数。

让我们使用网格搜索对逻辑回归模型进行超参数调整。我们将调整‘C’和‘penalty’超参数，并使用 ROC-AUC 作为我们的评分指标。让我们看看能否超越基准模型的表现。

现在，让我们开始这一部分的 Python 代码。

# Grid Search

# Import necessary libraries
from sklearn.model_selection import GridSearchCV

# Define the hyperparameters and their values
# 'C' is the inverse of regularization strength (smaller values specify stronger regularization)
# 'penalty' specifies the norm used in the penalization (l1 or l2)
hyperparameters = {'C': [0.001, 0.01, 0.1, 1, 10, 100, 1000], 
                   'penalty': ['l1', 'l2']}

# Implement grid search
# GridSearchCV is a method used to tune our model's hyperparameters
# We pass our model, the hyperparameters to tune, and the number of folds for cross-validation
# We're using ROC-AUC as our scoring metric
grid_search = GridSearchCV(LogisticRegression(), hyperparameters, cv=5, scoring='roc_auc')
grid_search.fit(X_train, y_train)

# Get the best hyperparameters
# GridSearchCV has found the best hyperparameters for our model, so we print them out
best_params = grid_search.best_params_
print(f'Best hyperparameters: {best_params}')

# Evaluate the best model
# GridSearchCV also gives us the best model, so we can use it to make predictions and evaluate its performance
best_model = grid_search.best_estimator_
y_pred_best = best_model.predict(X_test)
accuracy_best = accuracy_score(y_test, y_pred_best)
roc_auc_best = roc_auc_score(y_test, y_pred_best)

# Print the performance metrics of the best model
# These give us an indication of how well our model is performing after hyperparameter tuning
print("Grid Search Method " + f'Accuracy of the best model: {accuracy_best}')
print("Grid Search Method " + f'ROC-AUC of the best model: {roc_auc_best}')

经过网格搜索发现最佳超参数为{‘C’: 0.1, ‘penalty’: ‘l2’}，我们最佳模型的准确度为 0.852，ROC-AUC 分数为 0.853。有趣的是，这一表现稍低于我们的基准模型。这可能是因为基准模型的超参数已经非常适合这个数据集，或者可能是由于训练-测试划分中的随机性所致。不管怎样，这都提醒我们更复杂的模型和技术并不总是更好。

不过，你可能也注意到，我们的网格搜索只探索了相对较少的超参数组合。在实际操作中，超参数及其潜在值的数量可能会大得多，从而使网格搜索计算上昂贵甚至不可行。

这就是蒙特卡罗方法的作用。让我们看看这种更有指导性的方法是否能改善原始基线或基于网格搜索的模型的表现：

#Monte Carlo

# Import necessary libraries
from sklearn.metrics import accuracy_score, roc_auc_score
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import train_test_split
import numpy as np

# Split the data into training and test sets
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# Define the range of hyperparameters
# 'C' is the inverse of regularization strength (smaller values specify stronger regularization)
# 'penalty' specifies the norm used in the penalization (l1 or l2)
C_range = np.logspace(-3, 3, 7)
penalty_options = ['l1', 'l2']

# Initialize variables to store the best score and hyperparameters
best_score = 0
best_hyperparams = None

# Perform the Monte Carlo simulation
# We're going to perform 1000 iterations. You can play with this number to see how the performance changes.
# Remember the Law of Large Numbers!
for _ in range(1000):  

    # Randomly select hyperparameters from the defined range
    C = np.random.choice(C_range)
    penalty = np.random.choice(penalty_options)

    # Create and evaluate the model with these hyperparameters
    # We're using 'liblinear' solver as it supports both L1 and L2 regularization
    model = LogisticRegression(C=C, penalty=penalty, solver='liblinear')
    model.fit(X_train, y_train)
    y_pred = model.predict(X_test)

    # Calculate the accuracy and ROC-AUC
    accuracy = accuracy_score(y_test, y_pred)
    roc_auc = roc_auc_score(y_test, y_pred)

    # If this model's ROC-AUC is the best so far, store its score and hyperparameters
    if roc_auc > best_score:
        best_score = roc_auc
        best_hyperparams = {'C': C, 'penalty': penalty}

# Print the best score and hyperparameters
print("Monte Carlo Method " + f'Best ROC-AUC: {best_score}')
print("Monte Carlo Method " + f'Best hyperparameters: {best_hyperparams}')

# Train the model with the best hyperparameters
best_model = LogisticRegression(**best_hyperparams, solver='liblinear')
best_model.fit(X_train, y_train)

# Make predictions on the test set
y_pred = best_model.predict(X_test)

# Calculate and print the accuracy of the best model
accuracy = accuracy_score(y_test, y_pred)
print("Monte Carlo Method " + f'Accuracy of the best model: {accuracy}')

在蒙特卡罗方法中，我们发现最佳 ROC-AUC 分数为 0.9014，最佳超参数为{‘C’: 0.1, ‘penalty’: ‘l1’}。最佳模型的准确率为 0.9016。

两种技术的故事：网格搜索与蒙特卡罗

网格搜索通过评估超参数空间中所有可能组合在训练集上的性能，并根据预定义的指标（例如，准确率、ROC-AUC 等）选择产生最佳平均性能的组合，从而选择“最佳”超参数。这涉及将训练数据分成几个子集（或“折叠”，我们在代码中设置为 5），在一些折叠上训练模型，在其余折叠上评估模型，然后在所有折叠中平均性能。这被称为交叉验证，有助于通过确保模型在不同训练数据子集上的平均表现良好来减少过拟合。

相比之下，蒙特卡罗方法不会对不同训练数据子集进行任何平均。它随机选择超参数，在整个训练集上评估模型，然后选择在测试集上产生最佳性能的超参数。这种方法不使用任何交叉验证，因此不会对不同训练数据子集的性能进行平均。

在上述实验中，尽管网格搜索方法评估了蒙特卡罗方法选择的组合，但由于该组合在交叉验证过程中未能在不同训练数据子集上产生最佳平均性能，因此未被选为“最佳”超参数。然而，蒙特卡罗方法选择的组合恰好在本实验中使用的特定测试集上表现更好。这表明所选超参数在特定测试集的特征上表现良好，即使它们在不同训练数据子集上的平均表现不是最好。

这突显了网格搜索方法中通过对不同训练数据子集进行平均引入的偏差，与蒙特卡罗方法中通过在整个训练集上评估模型并根据单个测试集选择超参数引入的方差之间的权衡。

我鼓励你动手尝试 Python 代码，观察不同因素如何影响性能。你还可以比较两种方法在不同超参数空间下的计算时间，以了解它们的效率。这个练习旨在展示这些方法的动态，它们各有优缺点，并突出“最佳”方法可能取决于你的具体场景、计算资源和模型要求。因此，随意实验，祝学习愉快！

结论

起源于纸牌游戏的蒙特卡罗方法，无疑重新塑造了计算数学和数据科学的格局。它的力量在于其简洁性和多功能性，使我们能够相对轻松地处理复杂的高维问题。从通过掷飞镖估计π的值到在机器学习模型中调整超参数，蒙特卡罗模拟已经证明了其在数据科学武器库中的宝贵价值。

在这篇文章中，我们从蒙特卡罗方法的起源出发，了解了它的理论基础，并探讨了其在机器学习中的实际应用。我们看到了它如何用来优化机器学习模型，并通过使用真实数据集进行超参数调整的实践探索。我们还与其他方法进行了比较，展示了它的效率和有效性。

但蒙特卡罗的故事远未结束。随着我们不断推动机器学习和数据科学的边界，蒙特卡罗方法无疑将继续发挥重要作用。无论我们是在开发复杂的 AI 应用、理解复杂数据，还是只是玩一局纸牌游戏，蒙特卡罗方法都证明了模拟和近似在解决复杂问题中的强大力量。

当我们向前迈进时，让我们花一点时间欣赏这个方法的美妙——一个起源于简单纸牌游戏的方法，却能够驱动世界上最先进的计算。蒙特卡罗方法确实是一场高风险的机会与复杂性的游戏，到目前为止，它似乎总能获胜。因此，继续洗牌，继续打牌，并记住——在数据科学的游戏中，蒙特卡罗可能就是你手中的王牌。

离别感想

恭喜你成功到达终点！我们已经探索了概率的世界，挣扎于复杂的模型，并对蒙特卡罗模拟的强大能力有了新的认识。我们见证了它们如何将复杂的问题简化为可管理的组件，甚至优化了机器学习任务的超参数。

如果你像我一样喜欢深入探讨机器学习问题的复杂性，请关注我在 Medium 和 LinkedIn 上。让我们一起在人工智能的迷宫中，一次解决一个巧妙的问题。

直到我们的下一个统计冒险，继续探索，继续学习，继续模拟！在你数据科学和机器学习的旅程中，愿好运常伴你左右。

注意：除非另有说明，所有图片均由作者提供。