[人工智能]用于比较分类算法的统计测试

用于比较分类算法的统计测试

回顾开创性的论文和实施，以发现适合您的数据的最佳选择

对大多数数据科学家来说，比较预测方法以确定哪种方法应该用于手头的任务是一项日常活动。通常，我们会有一个分类模型库，并使用交叉验证来验证它们，以确定哪一个是最好的。

然而，另一个目标不是比较分类器，而是比较学习算法本身。其想法是：给定这个任务（数据），哪种学习算法（KNN、SVM、随机森林等）将在大小为D的数据集上生成更准确的分类器？

正如我们将看到的，这里介绍的每种方法都有一些优点和缺点。然而，使用两个比例测试的第一直觉可能会导致一些非常糟糕的结果。

为了更多地了解我们如何比较这些算法，并提高我们的统计学知识，今天我将解释和实现用于比较监督分类学习算法的近似统计测试[1]中的方法，这是一篇关于该领域的开创性论文。

在接下来的部分中，我将描述每个测试，讨论其优缺点，实现它们，然后将结果与可用的实现进行比较。

这篇文章的笔记本可以在 Kaggle 和 Github 上找到。

初始代码设置

对于本文中的代码，我们将使用两种分类算法：KNN和随机森林，在 sklearn 软件包上免费提供的UCI机器学习库中的葡萄酒数据集[2]上预测葡萄酒质量。为此，我们将导入一些必需的库，并实例化算法：

# Importing the required libs
import numpy as np
import pandas as pd

from tqdm import tqdm
from scipy.stats import norm, chi2
from scipy.stats import t as t_dist
from sklearn.datasets import load_wine
from sklearn.metrics import accuracy_score
from sklearn.ensemble import RandomForestClassifier
from sklearn.neighbors import KNeighborsClassifier
from sklearn.model_selection import train_test_split, KFold

# Libs implementations
from mlxtend.evaluate import mcnemar
from mlxtend.evaluate import mcnemar_table
from mlxtend.evaluate import paired_ttest_5x2cv
from mlxtend.evaluate import proportion_difference
from mlxtend.evaluate import paired_ttest_kfold_cv
from mlxtend.evaluate import paired_ttest_resampled# Getting the wine data from sklearn
X, y = load_wine(return_X_y = True)# Instantiating the classification algorithms
rf = RandomForestClassifier(random_state=42)
knn = KNeighborsClassifier(n_neighbors=1)# For holdout cases
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.20, random_state=42)

双比例试验(Two Proportions Test)

比较两个比例是一个非常古老的问题，统计学有一个经典的假设检验来解决这个问题：给定两个总体的两个比例，无效假设是比例之间的差异的平均值等于零。

我们可以用以下统计数据来计算：双比例检验统计量(Two proportions test statistic)
$$z=\frac{p_A?p_B}{\sqrt{2p(1?p)/n}}$$
看起来很简单，对吧？我们只需要得到我们算法的命中率（准确度）并进行比较。然而，这个测试有一个重要的假设：样本的独立性。

大家可以很快猜到，这里的样本不是独立的，因为两种算法的测试集和训练集都是相同的。所以这个假设是错误的。

这种方法还有两个问题：

1. 它不考虑测试集的方差。如果我们改变它，我们可能会有非常不同的结果
2. 它并不代表整个数据集，而是代表了一组较小的数据集，这些数据集被选择用于训练

要使用此测试，可以使用以下代码：

# First we fit the classification algorithms
rf.fit(X_train, y_train)
knn.fit(X_train, y_train)# Generate the predictions
rf_y = rf.predict(X_test)
knn_y = knn.predict(X_test)# Calculate the accuracy
acc1 = accuracy_score(y_test, rf_y)
acc2 = accuracy_score(y_test, knn_y)# Run the test
print("Proportions Z-Test")
z, p = proportion_difference(acc1, acc2, n_1=len(y_test))
print(f"z statistic: {z}, p-value: {p}\n")

在这里，我们只是在保持测试集上拟合算法，并对结果精度进行测试。

重抽样配对t检验(Resampled Paired t-test)

为了解释测试集的方差，可以使用重采样配对t检验。在这个测试中，我们将设置一些试验（例如30次），并使用一个保持测试集来测量每个试验中每个算法的准确性。

然后，如果我们假设 $p_i=p{A}_i?p{B}_i$，对于每个试验 $i$ 是正态分布的，我们可以应用配对学生的 $t$ 检验：配对t检验统计量(Paired t-test statistic)
$$t=\frac{\hat{p}\sqrt{n}}{\sqrt{\frac{{\sum }_{i=1}^n{\left(p^i?\hat{p}\right)}^2}{n?1}}}$$
因为在每次试验中，我们都会改变我们的测试集，它的方差会被考虑进去，从而改善了之前测试中的一个问题。然而，我们仍然面临一些问题：

1. $p_i$ 的正态分布不成立，因为比例不是独立的，因为它们是在同一测试集上计算的
2. 在每个试验中，训练集和测试集之间都有重叠，因此 $p_i$ 不是独立的
3. 它要求我们的算法被拟合多次，如果拟合时间太长，这可能是禁止的

对于该实现，我们将定义并创建一个函数，该函数将接收 $p_{i}s$ 作为参数：

def paired_t_test(p):
    p_hat = np.mean(p)
    n = len(p)
    den = np.sqrt(sum([(diff - p_hat)**2 for diff in p]) / (n - 1))
    t = (p_hat * (n**(1/2))) / den
    
    p_value = t_dist.sf(t, n-1)*2
    
    return t, p_value

在这个函数上，我们只是根据下面的等式，从测试中创建 t统计量(t statistic )。然后，我们使用scipy中的学生t分布找出测试的p值(p-value)，然后返回统计数据和p值(p-value)。

运行重采样t检验的代码如下：

n_tests = 30

p_ = []
rng = np.random.RandomState(42)
for i in range(n_tests):
    randint = rng.randint(low=0, high=32767)
    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.20, random_state=randint)
    rf.fit(X_train, y_train)
    knn.fit(X_train, y_train)

    acc1 = accuracy_score(y_test, rf.predict(X_test))
    acc2 = accuracy_score(y_test, knn.predict(X_test))
    p_.append(acc1 - acc2)
    
print("Paired t-test Resampled")
t, p = paired_t_test(p_)
print(f"t statistic: {t}, p-value: {p}\n")

在这里，我们将迭代次数定义为30次，对于每一次迭代，我们分割数据，拟合分类器，然后计算精度之间的差异。我们将这个值保存在一个列表中，然后用这个列表调用上面定义的函数。

Mlxtend库已经实现了此测试，因此还可以使用：

print("Paired t-test Resampled")
t, p = paired_ttest_resampled(estimator1=rf, estimator2=knn, X=X,y=y, random_seed=42, num_rounds=30, test_size=0.2)
print(f"t statistic: {t}, p-value: {p}\n")

请注意，我们必须通过函数的整个训练集，因为它将在自身内部创建拆分。如果给定相同的种子，可以验证结果是否相同。

交叉验证配对t检验(Cross-Validated Paired t-test)

此方法的结构与上述方法相同。但是，我们不会在每次试验中使用保持测试集，而是使用K倍进行交叉验证。

这将消除测试集重叠的问题，因为现在每个样本都将根据不同的数据进行测试。

然而，我们仍然存在训练数据重叠的问题。在10倍交叉验证中，每轮培训都与其他人分享80%的培训数据。

此外，我们仍然存在时间问题，因为我们必须多次调整分类器。然而，由于通常进行10倍或5倍交叉验证，这比我们之前做的30个测试小得多，因此通常比重取样配对t检验小。

对于这个测试，我们将使用我们在上一个测试中定义的相同函数，因为它将是一个t测试，我们只需要改变我们在循环的每个迭代上划分数据的方式（当然，循环上的迭代次数也是如此）。使用sklearn库中的KFold，我们有：

p_ = []

kf = KFold(n_splits=10, shuffle=True, random_state=42)
for train_index, test_index in kf.split(X):
    X_train, X_test, y_train, y_test = X[train_index], X[test_index], y[train_index], y[test_index]
    rf.fit(X_train, y_train)
    knn.fit(X_train, y_train)

    acc1 = accuracy_score(y_test, rf.predict(X_test))
    acc2 = accuracy_score(y_test, knn.predict(X_test))
    p_.append(acc1 - acc2)

print("Cross Validated Paired t-test")
t, p = paired_t_test(p_)
print(f"t statistic: {t}, p-value: {p}\n")

Mlxtend还为该测试提供了一个实现：

t, p = paired_ttest_kfold_cv(estimator1=rf, estimator2=knn, X=X, y=y, random_seed=42, shuffle=True, cv=10)
print(f"t statistic: {t}, p-value: {p}\n"

麦克内马尔试验

这种测试的优点是，我们的每个算法只需要一次拟合。我们使用一个保持集来适应其中的每一个，然后我们创建一个列联表：
McNemar’s Contingency Table(麦克内马尔列联表)

然后，我们陈述两种算法应具有相同错误率的零假设，并使用以下统计数据进行卡方检验：麦克内马尔检验统计量(McNemar’s test statistic)
$${\chi }^2=\frac{(|b?c|?1{)}^2}{b+c}$$
在纸质基准测试中，与其他测试相比，该测试在错误率方面排名第二，仅输给我们将在下面看到的5x2交叉验证测试。因此，作者说，如果你负担不起交叉验证，就应该使用这种方法。

然而，这方面仍然存在问题：

1. 该测试没有考虑列车组的变化，因为我们只适合一次
2. 因为我们只拟合一次，所以我们不考虑内部算法的变化
3. 我们使用比原版更小的一套。然而，请注意，这里的每个测试都会受到影响

为了实施测试，我们将创建一个专用功能：

def mcnemar_test(y_true, y_1, y_2):
    b = sum(np.logical_and((knn_y != y_test),(rf_y == y_test)))
    c = sum(np.logical_and((knn_y == y_test),(rf_y != y_test)))
    
    c_ = (np.abs(b - c) - 1)**2 / (b + c)
    
    p_value = chi2.sf(c_, 1)
    return c_, p_value

在这里，我们只是从列联表中计算数值，看看模型在哪里有正确或错误的答案。然后，我们查看卡方分布，找到我们计算的统计数据的p值。

由于这使用了一个保持集，以下步骤很简单：

print("McNemar's test")
chi2_, p = mcnemar_test(y_test, rf_y, knn_y)
print(f"chi2 statistic: {chi2_}, p-value: {p}\n")

此外，还可以使用Mlxtend库中的实现：

print("McNemar's test")
table = mcnemar_table(y_target=y_test, y_model1=rf_y, y_model2=knn_y)
chi2_, p = mcnemar(ary=table, corrected=True)
print(f"chi2 statistic: {chi2_}, p-value: {p}\n")

5x2交叉验证试验(5x2 Cross-Validation Test)

根据作者的基准测试，这项测试被认为是这5项测试中最好的一项。

其想法是进行5次2倍交叉验证，生成10种不同的估计。然后我们定义以下比例：5x2交叉验证测试参数(5x2 Cross-Validation test statistic)
$$?\begin{array}{c}{p}^1=p_A^1?p_B^1\\ p^2=p_A^2?p_B^2\\ \hat{p}=\left(p^1+p^2\right)/2\\ s^2={\left(p^1?\hat{p}\right)}^2+\left(p^2?\hat{p}\right)\end{array}$$
最后是统计数据：5x2交叉验证检验统计(5x2 Cross-Validation test statistic)
$$t=\frac{p_1^{(1)}}{\sqrt{\frac{1}{5}{\sum }_{i=1}^5{s}_i^2}}$$

这里最大的缺点是，我们必须多次适应算法。

这篇文章对这个方法有更全面的描述和演绎，所以我建议阅读它以获得全面的理解。

最后，为了实现它，我们将创建一个函数：

def five_two_statistic(p1, p2):
    p1 = np.array(p1)
    p2 = np.array(p2)
    p_hat = (p1 + p2) / 2
    s = (p1 - p_hat)**2 + (p2 - p_hat)**2
    t = p1[0] / np.sqrt(1/5. * sum(s))
    
    # p_value = t_dist.sf(t, 5)*2
    p_value = t_dist.sf(np.abs(t), 5) * 2.  # 这里要取绝对值哦
    
    return t, p_value

请注意，我们只是创建计算统计数据所需的值，然后一如既往地查看分布以找到 $p$ 值。

然后，我们进行五次双重交叉验证：

p_1 = []
p_2 = []

rng = np.random.RandomState(42)
for i in range(5):
    randint = rng.randint(low=0, high=32767)
    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.50, random_state=randint)

    rf.fit(X_train, y_train)
    knn.fit(X_train, y_train)
    acc1 = accuracy_score(y_test, rf.predict(X_test))
    acc2 = accuracy_score(y_test, knn.predict(X_test))
    p_1.append(acc1 - acc2)

    rf.fit(X_test, y_test)
    knn.fit(X_test, y_test)
    acc1 = accuracy_score(y_train, rf.predict(X_train))
    acc2 = accuracy_score(y_train, knn.predict(X_train))
    p_2.append(acc1 - acc2)# Running the test
print("5x2 CV Paired t-test")     
t, p = five_two_statistic(p_1, p_2)
print(f"t statistic: {t}, p-value: {p}\n")

我们还有Mlxtend实现：

print("5x2 CV Paired t-test")
t, p = paired_ttest_5x2cv(estimator1=rf, estimator2=knn, X=X, y=y, random_seed=42)
print(f"t statistic: {t}, p-value: {p}\n")

我们在这两种实现上得到了与预期相同的结果。

结论

重要的是要注意，在这种情况下没有银弹。这里提出的每个测试都有一些优点和缺点，它们都是近似值。但是，请注意，在任何情况下都不应使用双比例测试。

考虑到所需的时间预算，我们可以应用所有这些测试并比较它们的结果，以便更好地评估是否应该使用一类算法。

另一方面，如果感兴趣的算法可以计算其破碎系数（例如SVM、MLP或决策树），则这些测试可以与统计学习理论的结果一起使用，以确定应该使用哪种算法。但这是另一篇帖子的讨论。

我强烈建议阅读这篇论文，看看基准是如何构建的，它很容易阅读，而且信息丰富。

[1] Thomas G. Dietterich, Approximate Statistical Tests for Comparing Supervised Classification Learning Algorithms (1998), Neural Computation 1998; 10 (7): 1895–1923

[2] Lichman, M. (2013). UCI Machine Learning Repository [https://archive.ics.uci.edu/ml]. Irvine, CA: University of California, School of Information and Computer Science.

reference

@article{Toledo2022Jan,
author = {Toledo, Jr., T.},
title = {{Statistical Tests for Comparing Classification Algorithms}},
journaltitle = {Medium},
year = {2022},
month = {1},
date = {2022-01-04},
urldate = {2022-03-11},
publisher = {Towards Data Science},
language = {english},
hyphenation = {english},
url = {https://towardsdatascience.com/statistical-tests-for-comparing-classification-algorithms-ac1804e79bb7},
abstract = {{Comparing prediction methods to define which one should be used for the task at hand is a daily activity for most data scientists. Usually, one will have a pool of classification models and will…}}
}