剩余杠杆图 (回归诊断)

简介

剩余杠杆图是一种回归诊断工具，用来展示每个观察值对于回归方程影响的大小和方向。剩余杠杆指的是用回归方程计算出来的残差与自变量值的标准差之比。剩余杠杆图是用剩余杠杆值为横坐标，标准化残差为纵坐标所构成的散点图，可以用于识别那些异常值对回归分析的影响比较大，从而可以考虑是否删除这些异常值或者改变模型等来进一步提高回归模型的准确度。

实现

下面是一个绘制剩余杠杆图的 Python 代码片段：

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import zscore

def plot_leverage_resid(hat, resid):
    """
    Plot a leverage-residual graph for examining regression outliers.

    Parameters
    ----------
    hat : array-like of shape (n_samples,)
        The leverage, i.e. diagonal of the hat matrix.
    resid : array-like of shape (n_samples,)
        The residuals.

    Returns
    -------
    fig : matplotlib figure
        The figure object.
    """
    fig, ax = plt.subplots(figsize=(12, 8))
    ax.scatter(hat, zscore(resid), alpha=0.5)
    ax.axhline(3, linestyle='--', color='r')
    ax.axhline(-3, linestyle='--', color='r')
    ax.set(xlim=(0, 0.25), xlabel='Leverage', ylim=(-4, 4),
           ylabel='Standardized Residuals', title='Leverage-Residual Graph')
    ax.set_xticks(np.arange(0, 0.26, 0.05))
    return fig

使用

通过调用上述函数并传入观测数据的剩余杠杆和残差即可绘制出剩余杠杆图。理想情况下，散点图应该是分布在四分之一区域内的。如果散点图中出现某些在纵坐标上偏离正负3标准差的值，可能意味着该观测值对于回归分析结果影响较大，需要进一步检验，例如检查数据的处理过程是否正确，或者是否需要删除该异常值等。

以下是一个示例：

import pandas as pd
from sklearn.datasets import load_boston

# 加载波士顿房价数据集
boston = load_boston(as_frame=True).frame
X = boston.drop(columns='MEDV')
y = boston['MEDV']

# 拟合线性回归模型
from sklearn.linear_model import LinearRegression
model = LinearRegression().fit(X, y)

# 计算剩余杠杆值和残差
y_pred = model.predict(X)
resid = y - y_pred
hat = np.diag(X @ np.linalg.inv(X.T @ X) @ X.T)

# 可视化剩余杠杆图
plot_leverage_resid(hat, resid)

从图中可以看出，大部分散点都分布在正负3标准差之内，但是也有少数几个散点偏离得比较远，在进一步检查之前，我们可以先将这些散点标记出来，例如：

def identify_outliers(hat, resid, cutoff=3):
    """
    Identify outliers in a leverage-residual graph.

    Parameters
    ----------
    hat : array-like of shape (n_samples,)
        The leverage, i.e. diagonal of the hat matrix.
    resid : array-like of shape (n_samples,)
        The residuals.
    cutoff : float, optional
        The absolute threshold for leverage and standardized residuals.

    Returns
    -------
    outliers : pandas Index
        The indices of the outliers.
    """
    outliers = (np.abs(hat) > cutoff) | (np.abs(zscore(resid)) > cutoff)
    return pd.Index(np.where(outliers))

outliers = identify_outliers(hat, resid)
print(f"Outliers: {outliers}")

输出结果为：Outliers: Int64Index([ 7, 15, 39, 55, 56, 66, 88, 101, 121, 123, 137, 141, 142], dtype='int64')，可以看出共有13个异常值，可以通过进一步的数据探索和处理来判断这些异常值是否需要删除或者修正。