3、互信息计算：离散变量互信息计算、连续变量互信息估计（核密度估计、k近邻法）、互信息标准化

互信息（Mutual Information, MI）是衡量两个变量之间相互依赖程度的非参数度量。与皮尔逊相关系数仅捕捉线性关系不同，互信息能够检测任意类型的统计依赖关系（包括非线性、非单调关系），因此在MRM特征筛选中具有不可替代的地位。本节从离散变量互信息的精确计算出发，逐步深入到连续变量的估计方法，最后介绍标准化技术以消除量纲影响。

3.1 离散变量互信息计算

当两个变量均为离散类别型时，互信息可以直接通过联合概率分布与边缘概率分布精确计算。其核心公式为：

I(X;Y) = Σ Σ p(x,y) * log( p(x,y) / (p(x) * p(y)) )

其中求和遍历X和Y的所有取值组合。实际计算时，通常使用频率估计概率：

以下Python代码展示了离散互信息的直接实现：

import numpy as np
from collections import Counter

def discrete_mi(x, y):
    """计算两个离散变量的互信息"""
    n = len(x)
    # 联合分布
    joint = Counter(zip(x, y))
    # 边缘分布
    px = Counter(x)
    py = Counter(y)
    
    mi = 0.0
    for (xi, yj), n_ij in joint.items():
        p_xy = n_ij / n
        p_x = px[xi] / n
        p_y = py[yj] / n
        mi += p_xy * np.log(p_xy / (p_x * p_y))
    return mi

3.2 连续变量互信息估计

当变量为连续型时，无法直接计算联合概率分布，必须通过估计方法近似。常用的两类方法为核密度估计法和k近邻法。

3.2.1 核密度估计法（KDE）

核密度估计法通过平滑核函数（如高斯核）来估计连续变量的概率密度函数，进而计算互信息。其流程如下：

1. 联合密度估计：对二维样本点(x_i, y_i)使用二维核密度估计，得到\hat{p}(x,y)。
2. 边缘密度估计：分别对x和y进行一维核密度估计，得到\hat{p}(x)和\hat{p}(y)。
3. 数值积分：在定义域网格上离散化，计算I = Σ\hat{p}(x,y)·log(\hat{p}(x,y)/(\hat{p}(x)\hat{p}(y)))·Δx·Δy。

KDE法的优点是理论成熟、估计平滑；缺点是计算复杂度高（O(N²)），且在高维数据中表现不佳（维度灾难）。

3.2.2 k近邻法（kNN）

k近邻法（也称为KSG估计器，由Kraskov、Stögbauer和Grassberger提出）是目前最广泛使用的连续互信息估计方法。它不直接估计密度，而是基于样本点邻域的距离信息来估计互信息。

KSG估计器的核心公式为：

I(X;Y) = ψ(k) + ψ(N) - ⟨ψ(n_x+1) + ψ(n_y+1)⟩

其中ψ为digamma函数，N为样本总数，k为近邻数，n_x和n_y分别为在联合空间中第k近邻距离投影到X和Y边缘空间后包含的样本点数。

以下为KSG估计器的Python实现（基于Chebyshev距离）：

import numpy as np
from scipy.special import digamma
from sklearn.neighbors import NearestNeighbors

def ksg_mi(x, y, k=5):
    """KSG互信息估计器（k近邻法）"""
    n = len(x)
    data = np.column_stack([x, y])
    
    # 联合空间中的k近邻距离
    nn_joint = NearestNeighbors(metric='chebyshev', n_neighbors=k+1)
    nn_joint.fit(data)
    distances, _ = nn_joint.kneighbors(data)
    eps = distances[:, k]  # 第k近邻距离
    
    # 在X边缘空间中计数
    nn_x = NearestNeighbors(metric='chebyshev', n_neighbors=n)
    nn_x.fit(x.reshape(-1, 1))
    n_x = np.array([np.sum(nn_x.radius_neighbors([x[i]], eps[i])[0]) 
                    for i in range(n)])
    
    # 在Y边缘空间中计数
    nn_y = NearestNeighbors(metric='chebyshev', n_neighbors=n)
    nn_y.fit(y.reshape(-1, 1))
    n_y = np.array([np.sum(nn_y.radius_neighbors([y[i]], eps[i])[0]) 
                    for i in range(n)])
    
    # KSG公式
    mi = digamma(k) + digamma(n) - np.mean(digamma(n_x + 1) + digamma(n_y + 1))
    return mi

3.3 互信息标准化

原始互信息的取值范围为[0, +∞)，且受变量熵值影响——熵大的变量天然具有更高的互信息上界。这使得不同特征对之间的互信息值难以直接比较。标准化（归一化）技术将互信息映射到固定区间（通常为[0,1]），实现特征重要性排序的公平性。

常用的标准化方法包括：

以下为SU的计算实现（基于离散互信息或KSG估计）：

from scipy.stats import entropy

def symmetric_uncertainty(x, y, mi_func=ksg_mi, **kwargs):
    """计算对称不确定性（SU）"""
    mi = mi_func(x, y, **kwargs)
    # 计算边缘熵
    h_x = entropy(np.histogram(x, bins='auto')[0] / len(x))
    h_y = entropy(np.histogram(y, bins='auto')[0] / len(y))
    # 防止除零
    if h_x + h_y == 0:
        return 0.0
    return 2.0 * mi / (h_x + h_y)

3.4 方法选择与实战建议

在实际MRM特征筛选项目中，变量类型组合决定了互信息计算方法的选择：