逻辑回归LR与广义线性模型GLM

逻辑回归（Logistic Regression, LR）和广义线性模型（Generalized Linear Model, GLM）是统计建模和机器学习中的重要工具。它们常用于分类问题和回归分析。以下是关于逻辑回归和广义线性模型开发实践的一些专业知识：

1. 广义线性模型（GLM）概述

广义线性模型是一个扩展的线性回归模型，它允许响应变量（因变量）具有不同的分布族，而不仅仅是正态分布。GLM包括三个主要组成部分：

逻辑回归是一种用于解决二分类问题的对数线性模型。它的核心思想是将逻辑函数（sigmoid函数）应用于线性回归的预测值，从而将原始的实数输出转换为概率形式。这样，我们就可以通过设定一个阈值来将概率转换为最终的分类结果。

在实践中，逻辑回归具有简单、高效和易于解释等优点。然而，它也有一些局限性，例如对异常值的敏感性以及在多分类问题中的表现不佳。

• 线性预测器：表示为 η=Xβ\eta = X\betaη=Xβ，其中 XXX 是输入特征矩阵，β\betaβ 是回归系数。
• 链接函数：将线性预测器 η\etaη 映射到响应变量的均值。例如，逻辑回归使用逻辑函数作为链接函数。
• 分布族：响应变量的概率分布，如二项分布（逻辑回归）、泊松分布（泊松回归）、正态分布（线性回归）等。

2. 逻辑回归（LR）

逻辑回归是广义线性模型的一个特例，主要用于二分类问题。它的核心思想是通过逻辑函数将线性预测器的输出映射到一个介于0和1之间的概率值。

逻辑回归模型的主要组成部分：

• 线性预测器：η=Xβ\eta = X\betaη=Xβ
• 逻辑链接函数：p=11+e−(Xβ)p = \frac{1}{1 + e^{-(X\beta)}}p=1+e−(Xβ)1​，其中 ppp 是事件发生的概率。
• 损失函数：逻辑回归使用的是交叉熵损失函数，也称为对数损失（log-loss）。

逻辑回归模型的开发实践：

1. 数据预处理：

   • 特征标准化：标准化特征有助于加快收敛速度，特别是在梯度下降算法中。
   • 处理不平衡数据：对于类别不平衡的数据，可以使用过采样、欠采样或生成对抗样本等方法进行处理。

2. 模型拟合：

   • 选择优化算法：逻辑回归可以使用梯度下降法、拟牛顿法（如BFGS）或其他优化算法进行拟合。
   • 正则化：使用L1（Lasso）或L2（Ridge）正则化来防止过拟合。

3. 模型评估：

   • 评估指标：常用指标包括准确率、精确率、召回率、F1分数、ROC曲线和AUC值。
   • 交叉验证：使用k折交叉验证评估模型的稳定性和泛化能力。

4. 超参数调优：

   • 选择正则化参数：通过交叉验证选择最佳的正则化参数。

除了逻辑回归，GLM还包括其他回归模型，如：

• 泊松回归：用于计数数据，响应变量服从泊松分布。适用于事件发生的次数建模。

  • 链接函数：对数链接函数，即 log⁡(λ)=Xβ\log(\lambda) = X\betalog(λ)=Xβ，其中 λ\lambdaλ 是事件的发生率。

• 负二项回归：用于过度离散的计数数据，响应变量服从负二项分布。

• Gamma回归：用于正值连续响应变量，适合处理响应变量为非负连续数据。

1. 选择合适的分布和链接函数：

   • 根据响应变量的类型选择合适的分布族和链接函数。

2. 数据预处理：

   • 处理缺失值：根据情况选择插补或删除缺失值。
   • 特征选择：选择对响应变量有显著影响的特征。

3. 模型拟合与评估：

   • 使用最大似然估计（MLE）进行参数估计。
   • 对于泊松回归和Gamma回归，可以使用伪R²、AIC（赤池信息准则）和BIC（贝叶斯信息准则）等指标评估模型的拟合优度。

4. 诊断与调整：

   • 诊断残差：检查残差的分布，识别模型假设是否合理。
   • 模型调整：根据诊断结果进行模型调整，如变换响应变量或选择不同的链接函数。

• 多重共线性：在特征间存在高度相关性时，可能导致回归系数不稳定。可以通过主成分分析（PCA）或其他降维技术来缓解。
• 异常值：对模型有较大影响的异常值需要特别处理，如使用鲁棒回归方法。
• 数据不平衡：逻辑回归和其他GLM在处理数据不平衡问题时可能表现不佳，考虑使用加权损失函数或其他处理技术。

适用问题类型：逻辑回归适用于二分类问题；而GLM可以应用于多分类、有序分类和回归问题。
模型假设：逻辑回归假设因变量服从二项分布；而GLM的假设更为灵活，可以根据具体问题选择不同的分布形式。

IRLS 优势
收敛速度：IRLS通常在少数几次迭代后就能快速收敛，特别适用于样本量较大时。它在许多实际问题中表现出良好的收敛性。

简单易行：IRLS的实现相对简单，基于最小二乘法的迭代加权更新，易于理解和实现。

鲁棒性：IRLS在处理某些类型的响应变量（如二项分布、泊松分布等）时表现良好，特别适用于这些分布的特定性质。

广泛适用性：IRLS不仅适用于广义线性回归模型，也可以用于其他一些需要迭代加权的统计模型。

隐语的GLM是基于秘密共享mpc协议实现的，计算的算子基本都可以化简为加法和乘法，因此列举了这两类算子的秘密贡献计算逻辑，加法比较容易理解，乘法需要借助beaver三元组辅助实现。

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
 
from functools import partial
from sklearn.datasets import fetch_openml
from sklearn.metrics import (
    mean_absolute_error,
    mean_squared_error,
    mean_tweedie_deviance,
)
 
 
def load_mtpl2(n_samples=None):
    """获取法国机动车第三方责任保险索赔数据集。
    参数
    ----------
    n_samples: int, 默认值=None
      要选择的样本数量（为了更快的运行时间）。完整数据集有 678013 个样本。
    """
    # 从 https://www.openml.org/d/41214 获取 freMTPL2freq 数据集
    df_freq = fetch_openml(data_id=41214, as_frame=True).data
    df_freq["IDpol"] = df_freq["IDpol"].astype(int)
    df_freq.set_index("IDpol", inplace=True)
 
    # 从 https://www.openml.org/d/41215 获取 freMTPL2sev 数据集
    df_sev = fetch_openml(data_id=41215, as_frame=True).data
 
    # 按照相同的 ID 累计索赔金额
    df_sev = df_sev.groupby("IDpol").sum()
 
    df = df_freq.join(df_sev, how="left")
    df["ClaimAmount"].fillna(0, inplace=True)
 
    # 去掉字符串字段的引号
    for column_name in df.columns[df.dtypes.values == object]:
        df[column_name] = df[column_name].str.strip("'")
    return df.iloc[:n_samples]
 
 
def plot_obs_pred(
    df,
    feature,
    weight,
    observed,
    predicted,
    y_label=None,
    title=None,
    ax=None,
    fill_legend=False,
):
    """绘制观察值和预测值 - 按特征级别聚合。
    参数
    ----------
    df : DataFrame
        输入数据
    feature: str
        要绘制的特征列名
    weight : str
        df 中权重或曝光量的列名
    observed : str
        df 中观察目标的列名
    predicted : DataFrame
        具有与 df 相同索引的预测目标数据框
    fill_legend : bool, 默认值=False
        是否显示填充区域的图例
    """
    # 按特征级别聚合观察变量和预测变量
    df_ = df.loc[:, [feature, weight]].copy()
    df_["observed"] = df[observed] * df[weight]
    df_["predicted"] = predicted * df[weight]
    df_ = (
        df_.groupby([feature])[[weight, "observed", "predicted"]]
        .sum()
        .assign(observed=lambda x: x["observed"] / x[weight])
        .assign(predicted=lambda x: x["predicted"] / x[weight])
    )
 
    ax = df_.loc[:, ["observed", "predicted"]].plot(style=".", ax=ax)
    y_max = df_.loc[:, ["observed", "predicted"]].values.max() * 0.8
    p2 = ax.fill_between(
        df_.index,
        0,
        y_max * df_[weight] / df_[weight].values.max(),
        color="g",
        alpha=0.1,
    )
    if fill_legend:
        ax.legend([p2], ["{} distribution".format(feature)])
    ax.set(
        ylabel=y_label if y_label is not None else None,
        title=title if title is not None else "Train: Observed vs Predicted",
    )
 
 
def score_estimator(
    estimator,
    X_train,
    X_test,
    df_train,
    df_test,
    target,
    weights,
    tweedie_powers=None,
):
    """使用不同的指标评估估计器在训练集和测试集上的表现"""
 
    metrics = [
        ("D² explained", None),  # 如果存在，使用默认评分器
        ("mean abs. error", mean_absolute_error),
        ("mean squared error", mean_squared_error),
    ]
    if tweedie_powers:
        metrics += [
            (
                "mean Tweedie dev p={:.4f}".format(power),
                partial(mean_tweedie_deviance, power=power),
            )
            for power in tweedie_powers
        ]
 
    res = []
    for subset_label, X, df in [
        ("train", X_train, df_train),
        ("test", X_test, df_test),
    ]:
        y, _weights = df[target], df[weights]
        for score_label, metric in metrics:
            if isinstance(estimator, tuple) and len(estimator) == 2:
                # 对由频率模型和严重性模型乘积组成的模型进行评分
                est_freq, est_sev = estimator
                y_pred = est_freq.predict(X) * est_sev.predict(X)
            else:
                y_pred = estimator.predict(X)
 
            if metric is None:
                if not hasattr(estimator, "score"):
                    continue
                score = estimator.score(X, y, sample_weight=_weights)
            else:
                score = metric(y, y_pred, sample_weight=_weights)
 
            res.append({"subset": subset_label, "metric": score_label, "score": score})
 
    res = (
        pd.DataFrame(res)
        .set_index(["metric", "subset"])
        .score.unstack(-1)
        .round(4)
        .loc[:, ["train", "test"]]
    )
    return res

除了IRLS，拟牛顿法（Quasi-Newton Methods）也是经常被使用的二阶优化器。L-BFGS（Limited-memory BFGS）：适用于大规模问题的BFGS变种，利用有限内存近似海森矩阵

计算效率：对于大规模数据集，逻辑回归通常更快，因为它的计算过程相对简单；而GLM可能需要更多的计算资源和时间来拟合模型。
解释性：逻辑回归的输出可以直接解释为概率，因此更容易解释；而GLM的输出解释可能需要更多的专业知识或统计技巧。
综上所述，逻辑回归和广义线性模型都是强大的对数线性模型工具。在选择使用哪一个时，我们需要考虑问题的类型、数据的特性以及计算资源的限制等因素。在实践中，我们可以根据具体的需求和限制来选择最合适的模型。例如，对于二分类问题，如果数据的分布符合二项分布的假设并且计算资源有限，逻辑回归可能是一个更好的选择。而对于多分类或有序分类问题，或者当因变量的分布不符合任何预先设定的形式时，广义线性模型可能更加合适。

此外，对于需要同时处理分类和回归问题的情况，我们可以考虑使用一种称为“广义线性混合模型”的方法。这种方法结合了GLM的灵活性和混合模型的强大功能，允许我们在同一框架内处理各种类型的数据和问题。这将是我们在未来深入研究的一个重要方向。