稳健协方差估计与马氏距离的相关性

简介

在本实验中，我们将探讨在高斯分布数据上使用具有马氏距离的稳健协方差估计。马氏距离是衡量一个点与一个分布之间距离的指标。它被定义为一个点与分布均值之间的距离，并通过分布协方差矩阵的逆进行缩放。对于高斯分布的数据，马氏距离可用于计算一个观测值到分布众数的距离。我们将在计算一个受污染数据集的马氏距离时，比较协方差的稳健估计器——最小协方差行列式（MCD）估计器与标准协方差最大似然估计器（MLE）的性能。

虚拟机使用提示

虚拟机启动完成后，点击左上角切换到“笔记本”标签页，以访问 Jupyter Notebook 进行练习。

有时，你可能需要等待几秒钟让 Jupyter Notebook 完成加载。由于 Jupyter Notebook 的限制，操作验证无法自动化。

如果你在学习过程中遇到问题，随时向 Labby 提问。课程结束后提供反馈，我们会及时为你解决问题。

Skills Graph

生成数据

首先，我们生成一个包含 125 个样本和 2 个特征的数据集。两个特征均呈高斯分布，均值为 0。然而，特征 1 的标准差等于 2，特征 2 的标准差等于 1。接下来，我们用高斯异常值样本替换 25 个样本，其中特征 1 的标准差等于 1，特征 2 的标准差等于 7。

import numpy as np

## for consistent results
np.random.seed(7)

n_samples = 125
n_outliers = 25
n_features = 2

## generate Gaussian data of shape (125, 2)
gen_cov = np.eye(n_features)
gen_cov[0, 0] = 2.0
X = np.dot(np.random.randn(n_samples, n_features), gen_cov)
## add some outliers
outliers_cov = np.eye(n_features)
outliers_cov[np.arange(1, n_features), np.arange(1, n_features)] = 7.0
X[-n_outliers:] = np.dot(np.random.randn(n_outliers, n_features), outliers_cov)

将 MCD 和 MLE 协方差估计器应用于数据

我们将把基于 MCD 和 MLE 的协方差估计器应用于我们的数据，并打印估计出的协方差矩阵。请注意，基于 MLE 的估计器得出的特征 2 的估计方差（7.5）比 MCD 稳健估计器（1.2）高得多。这表明基于 MCD 的稳健估计器对异常值样本的抵抗力要强得多，这些异常值样本被设计为在特征 2 上具有大得多的方差。

from sklearn.covariance import EmpiricalCovariance, MinCovDet

## fit a MCD robust estimator to data
robust_cov = MinCovDet().fit(X)
## fit a MLE estimator to data
emp_cov = EmpiricalCovariance().fit(X)
print(
    "Estimated covariance matrix:\nMCD (Robust):\n{}\nMLE:\n{}".format(
        robust_cov.covariance_, emp_cov.covariance_
    )
)

绘制马氏距离的等高线

我们将绘制两种方法计算出的马氏距离的等高线。请注意，基于稳健 MCD 的马氏距离能更好地拟合内点（黑色点），而基于 MLE 的距离受异常值（红色点）的影响更大。

import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(10, 5))
## Plot data set
inlier_plot = ax.scatter(X[:, 0], X[:, 1], color="black", label="inliers")
outlier_plot = ax.scatter(
    X[:, 0][-n_outliers:], X[:, 1][-n_outliers:], color="red", label="outliers"
)
ax.set_xlim(ax.get_xlim()[0], 10.0)
ax.set_title("Mahalanobis distances of a contaminated data set")

## Create meshgrid of feature 1 and feature 2 values
xx, yy = np.meshgrid(
    np.linspace(plt.xlim()[0], plt.xlim()[1], 100),
    np.linspace(plt.ylim()[0], plt.ylim()[1], 100),
)
zz = np.c_[xx.ravel(), yy.ravel()]
## Calculate the MLE based Mahalanobis distances of the meshgrid
mahal_emp_cov = emp_cov.mahalanobis(zz)
mahal_emp_cov = mahal_emp_cov.reshape(xx.shape)
emp_cov_contour = plt.contour(
    xx, yy, np.sqrt(mahal_emp_cov), cmap=plt.cm.PuBu_r, linestyles="dashed"
)
## Calculate the MCD based Mahalanobis distances
mahal_robust_cov = robust_cov.mahalanobis(zz)
mahal_robust_cov = mahal_robust_cov.reshape(xx.shape)
robust_contour = ax.contour(
    xx, yy, np.sqrt(mahal_robust_cov), cmap=plt.cm.YlOrBr_r, linestyles="dotted"
)

## Add legend
ax.legend(
    [
        emp_cov_contour.collections[1],
        robust_contour.collections[1],
        inlier_plot,
        outlier_plot,
    ],
    ["MLE dist", "MCD dist", "inliers", "outliers"],
    loc="upper right",
    borderaxespad=0,
)

plt.show()

比较 MLE 和 MCD 马氏距离

我们将突出基于 MCD 的马氏距离区分异常值的能力。我们对马氏距离取立方根，得到近似正态分布。然后，我们用箱线图绘制内点和异常值样本的值。对于基于稳健 MCD 的马氏距离，异常值样本的分布与内点样本的分布分得更开。

fig, (ax1, ax2) = plt.subplots(1, 2)
plt.subplots_adjust(wspace=0.6)

## Calculate cubic root of MLE Mahalanobis distances for samples
emp_mahal = emp_cov.mahalanobis(X - np.mean(X, 0)) ** (0.33)
## Plot boxplots
ax1.boxplot([emp_mahal[:-n_outliers], emp_mahal[-n_outliers:]], widths=0.25)
## Plot individual samples
ax1.plot(
    np.full(n_samples - n_outliers, 1.26),
    emp_mahal[:-n_outliers],
    "+k",
    markeredgewidth=1,
)
ax1.plot(np.full(n_outliers, 2.26), emp_mahal[-n_outliers:], "+k", markeredgewidth=1)
ax1.axes.set_xticklabels(("inliers", "outliers"), size=15)
ax1.set_ylabel(r"$\sqrt[3]{\rm{(Mahal. dist.)}}$", size=16)
ax1.set_title("Using non-robust estimates\n(Maximum Likelihood)")

## Calculate cubic root of MCD Mahalanobis distances for samples
robust_mahal = robust_cov.mahalanobis(X - robust_cov.location_) ** (0.33)
## Plot boxplots
ax2.boxplot([robust_mahal[:-n_outliers], robust_mahal[-n_outliers:]], widths=0.25)
## Plot individual samples
ax2.plot(
    np.full(n_samples - n_outliers, 1.26),
    robust_mahal[:-n_outliers],
    "+k",
    markeredgewidth=1,
)
ax2.plot(np.full(n_outliers, 2.26), robust_mahal[-n_outliers:], "+k", markeredgewidth=1)
ax2.axes.set_xticklabels(("inliers", "outliers"), size=15)
ax2.set_ylabel(r"$\sqrt[3]{\rm{(Mahal. dist.)}}$", size=16)
ax2.set_title("Using robust estimates\n(Minimum Covariance Determinant)")

plt.show()

总结

在本实验中，我们学习了如何在高斯分布数据上使用马氏距离进行稳健协方差估计。我们比较了基于 MCD 和 MLE 的协方差估计器在计算受污染数据集的马氏距离时的性能。我们观察到，基于 MCD 的稳健估计器对异常值样本的抵抗力要强得多，并且基于 MCD 的马氏距离在区分异常值方面表现更好。