1.1. 线性模型#

以下是一组用于回归的方法，其中目标值预期是特征的线性组合。用数学符号表示，如果 \(\hat{y}\) 是预测值。

\[\hat{y}(w, x) = w_0 + w_1 x_1 + ... + w_p x_p\]

在本模块中，我们将向量 \(w = (w_1, ..., w_p)\) 称为 coef_，将 \(w_0\) 称为 intercept_。

要使用广义线性模型执行分类，请参阅逻辑回归。

1.1.1. 普通最小二乘法#

LinearRegression 通过系数 \(w = (w_1, ..., w_p)\) 拟合线性模型，以最小化数据集中观测目标与线性近似预测目标之间的残差平方和。数学上，它解决的是以下形式的问题：

\[\min_{w} || X w - y||_2^2\]
../_images/sphx_glr_plot_ols_ridge_001.png

LinearRegression 的 fit 方法接受参数 X、y、sample_weight，并将线性模型的系数 \(w\) 存储在其 coef_ 和 intercept_ 属性中。

>>> from sklearn import linear_model
>>> reg = linear_model.LinearRegression()
>>> reg.fit([[0, 0], [1, 1], [2, 2]], [0, 1, 2])
LinearRegression()
>>> reg.coef_
array([0.5, 0.5])
>>> reg.intercept_
0.0

普通最小二乘法的系数估计依赖于特征的独立性。当特征相关且设计矩阵 \(X\) 的某些列近似线性相关时，设计矩阵将变得接近奇异，结果是最小二乘估计对观测目标中的随机误差高度敏感，产生较大的方差。这种多重共线性情况可能出现，例如，在没有实验设计的情况下收集数据时。

示例

1.1.1.1. 非负最小二乘法#

可以限制所有系数为非负数，这在系数代表某些物理或自然非负量（例如，频率计数或商品价格）时可能很有用。LinearRegression 接受一个布尔型参数 positive：当设置为 True 时，将应用非负最小二乘法。

示例

1.1.1.2. 普通最小二乘法的复杂度#

最小二乘解是使用 \(X\) 的奇异值分解计算的。如果 \(X\) 是一个形状为 (n_samples, n_features) 的矩阵，则该方法的成本为 \(O(n_{\text{samples}} n_{\text{features}}^2)\)，假设 \(n_{\text{samples}} \geq n_{\text{features}}\)。

1.1.2. 岭回归与分类#

1.1.2.1. 回归#

Ridge 回归通过对系数大小施加惩罚来解决普通最小二乘法的一些问题。岭系数最小化一个惩罚残差平方和：

\[\min_{w} || X w - y||_2^2 + \alpha ||w||_2^2\]

复杂性参数 \(\alpha \geq 0\) 控制收缩量：\(\alpha\) 的值越大，收缩量越大，因此系数对共线性更具鲁棒性。

../_images/sphx_glr_plot_ridge_path_001.png

与其他线性模型一样，Ridge 将在其 fit 方法中接受数组 X、y，并将线性模型的系数 \(w\) 存储在其 coef_ 成员中。

>>> from sklearn import linear_model
>>> reg = linear_model.Ridge(alpha=.5)
>>> reg.fit([[0, 0], [0, 0], [1, 1]], [0, .1, 1])
Ridge(alpha=0.5)
>>> reg.coef_
array([0.34545455, 0.34545455])
>>> reg.intercept_
np.float64(0.13636)

请注意，类 Ridge 允许用户通过设置 solver="auto" 来自动选择求解器。当指定此选项时，Ridge 将在 "lbfgs"、"cholesky" 和 "sparse_cg" 求解器之间进行选择。Ridge 将从上到下检查下表中所示的条件。如果条件为真，则选择相应的求解器。

求解器

条件

‘lbfgs’

指定了 positive=True 选项。

‘cholesky’

输入数组 X 不是稀疏的。

‘sparse_cg’

以上条件均未满足。

示例

1.1.2.2. 分类#

Ridge 回归器有一个分类器变体：RidgeClassifier。此分类器首先将二元目标转换为 {-1, 1}，然后将问题视为回归任务，优化与上述相同的目标。预测类别对应于回归器预测的符号。对于多类别分类，问题被视为多输出回归，预测类别对应于具有最高值的输出。

使用（惩罚的）最小二乘损失来拟合分类模型而不是更传统的逻辑损失或铰链损失可能看起来有些疑问。然而，在实践中，所有这些模型在准确性或精确度/召回率方面都可以产生相似的交叉验证分数，而 RidgeClassifier 使用的惩罚最小二乘损失允许选择非常不同的数值求解器，具有不同的计算性能特征。

当类别数量很多时，RidgeClassifier 可以比（例如）LogisticRegression 快得多，因为它只需计算一次投影矩阵 \((X^T X)^{-1} X^T\)。

此分类器有时被称为具有线性核的最小二乘支持向量机。

示例

1.1.2.3. 岭回归复杂度#

此方法的复杂度与普通最小二乘法相同。

1.1.2.4. 设置正则化参数：留一法交叉验证#

RidgeCV 和 RidgeClassifierCV 实现了岭回归/分类，并内置了 alpha 参数的交叉验证。它们的工作方式与 GridSearchCV 相同，只是默认使用高效的留一法交叉验证。当使用默认的交叉验证时，由于计算留一法误差所使用的公式，alpha 不能为 0。详情请参阅 [RL2007]。

使用示例

>>> import numpy as np
>>> from sklearn import linear_model
>>> reg = linear_model.RidgeCV(alphas=np.logspace(-6, 6, 13))
>>> reg.fit([[0, 0], [0, 0], [1, 1]], [0, .1, 1])
RidgeCV(alphas=array([1.e-06, 1.e-05, 1.e-04, 1.e-03, 1.e-02, 1.e-01, 1.e+00, 1.e+01,
      1.e+02, 1.e+03, 1.e+04, 1.e+05, 1.e+06]))
>>> reg.alpha_
np.float64(0.01)

指定 cv 属性的值将触发与 GridSearchCV 的交叉验证使用，例如 cv=10 用于 10 折交叉验证，而不是留一法交叉验证。

参考文献#
[RL2007]

“正则化最小二乘法笔记”，Rifkin & Lippert (技术报告, 课程幻灯片)。

1.1.3. Lasso#

Lasso 是一种估计稀疏系数的线性模型。在某些情况下它很有用，因为它倾向于选择具有较少非零系数的解，从而有效地减少给定解所依赖的特征数量。因此，Lasso 及其变体是压缩感知领域的基础。在某些条件下，它可以恢复非零系数的精确集合（参见压缩感知：使用 L1 先验（Lasso）进行层析重建）。

数学上，它由一个线性模型和附加的正则化项组成。要最小化的目标函数是：

\[\min_{w} { \frac{1}{2n_{\text{samples}}} ||X w - y||_2 ^ 2 + \alpha ||w||_1}\]

因此，Lasso 估计通过添加 \(\alpha ||w||_1\) 来解决最小二乘惩罚的最小化问题，其中 \(\alpha\) 是一个常数，\(||w||_1\) 是系数向量的 \(\ell_1\) 范数。

类 Lasso 中的实现使用坐标下降作为拟合系数的算法。有关另一种实现，请参阅最小角回归。

>>> from sklearn import linear_model
>>> reg = linear_model.Lasso(alpha=0.1)
>>> reg.fit([[0, 0], [1, 1]], [0, 1])
Lasso(alpha=0.1)
>>> reg.predict([[1, 1]])
array([0.8])

函数 lasso_path 对于较低级别的任务很有用，因为它沿所有可能值的完整路径计算系数。

示例

注意

使用 Lasso 进行特征选择

由于 Lasso 回归产生稀疏模型，因此可用于执行特征选择，详情请参阅基于 L1 的特征选择。

参考文献#

以下两篇参考文献解释了 scikit-learn 坐标下降求解器中使用的迭代，以及用于收敛控制的对偶间隙计算。

  • “基于坐标下降的广义线性模型正则化路径”，Friedman, Hastie & Tibshirani, J Stat Softw, 2010 (论文)。

  • “用于大规模 L1 正则化最小二乘的内点法”，S. J. Kim, K. Koh, M. Lustig, S. Boyd and D. Gorinevsky，发表于 IEEE Journal of Selected Topics in Signal Processing, 2007 (论文)

1.1.3.1. 设置正则化参数#

alpha 参数控制估计系数的稀疏程度。

1.1.3.1.1. 使用交叉验证#

scikit-learn 提供了通过交叉验证设置 Lasso alpha 参数的对象：LassoCV 和 LassoLarsCV。LassoLarsCV 基于下文解释的最小角回归算法。

对于具有许多共线性特征的高维数据集，通常首选 LassoCV。然而，LassoLarsCV 的优点是能够探索更多相关的 alpha 参数值，并且如果样本数量与特征数量相比非常小，它通常比 LassoCV 更快。

lasso_cv_1 lasso_cv_2

1.1.3.1.2. 基于信息准则的模型选择#

或者，估计器 LassoLarsIC 建议使用赤池信息准则（AIC）和贝叶斯信息准则（BIC）。这是一种计算成本更低的替代方法，用于找到最佳 alpha 值，因为正则化路径只计算一次，而不是在使用 k 折交叉验证时计算 k+1 次。

实际上，这些准则是在样本内训练集上计算的。简而言之，它们通过 Lasso 模型的灵活性来惩罚其过于乐观的分数（参见下文“数学细节”部分）。

然而，这些准则需要对解的自由度进行适当的估计，并且是针对大样本（渐近结果）推导出来的，并假设正在研究的候选模型是正确的模型。当问题条件不佳时（例如，特征多于样本），它们也容易失效。

../_images/sphx_glr_plot_lasso_lars_ic_001.png

示例

1.1.3.1.3. AIC 和 BIC 准则#

AIC（以及 BIC）的定义在文献中可能有所不同。在本节中，我们将提供更多关于 scikit-learn 中计算准则的信息。

数学细节#

AIC 准则定义为：

\[AIC = -2 \log(\hat{L}) + 2 d\]

其中 \(\hat{L}\) 是模型的最大似然，\(d\) 是参数数量（在上一节中也称为自由度）。

BIC 的定义将常数 \(2\) 替换为 \(\log(N)\)：

\[BIC = -2 \log(\hat{L}) + \log(N) d\]

其中 \(N\) 是样本数量。

对于线性高斯模型，最大对数似然定义为：

\[\log(\hat{L}) = - \frac{n}{2} \log(2 \pi) - \frac{n}{2} \ln(\sigma^2) - \frac{\sum_{i=1}^{n} (y_i - \hat{y}_i)^2}{2\sigma^2}\]

其中 \(\sigma^2\) 是噪声方差的估计，\(y_i\) 和 \(\hat{y}_i\) 分别是真实目标和预测目标，\(n\) 是样本数量。

将最大对数似然代入 AIC 公式得到：

\[AIC = n \log(2 \pi \sigma^2) + \frac{\sum_{i=1}^{n} (y_i - \hat{y}_i)^2}{\sigma^2} + 2 d\]

上述表达式的第一项有时会被舍弃，因为当提供 \(\sigma^2\) 时，它是一个常数。此外，有时会指出 AIC 等同于 \(C_p\) 统计量 [12]。然而，严格来说，它仅在某些常数和乘法因子上等同。

最后，我们上面提到 \(\sigma^2\) 是噪声方差的估计。在 LassoLarsIC 中，当未提供参数 noise_variance（默认）时，噪声方差通过无偏估计器 [13] 进行估计，定义为：

\[\sigma^2 = \frac{\sum_{i=1}^{n} (y_i - \hat{y}_i)^2}{n - p}\]

其中 \(p\) 是特征数量，\(\hat{y}_i\) 是使用普通最小二乘回归预测的目标。请注意，此公式仅在 n_samples > n_features 时有效。

参考文献

1.1.3.1.4. 与 SVM 正则化参数的比较#

alpha 和 SVM 的正则化参数 C 之间的等价关系由 alpha = 1 / C 或 alpha = 1 / (n_samples * C) 给出，具体取决于估计器和模型优化的精确目标函数。

1.1.4. 多任务 Lasso#

MultiTaskLasso 是一种线性模型，用于联合估计多个回归问题（y 是一个形状为 (n_samples, n_tasks) 的二维数组）的稀疏系数。约束是所有回归问题（也称为任务）选择的特征相同。

下图比较了通过简单 Lasso 或 MultiTaskLasso 获得的系数矩阵 W 中非零项的位置。Lasso 估计产生分散的非零项，而 MultiTaskLasso 的非零项是完整的列。

multi_task_lasso_1 multi_task_lasso_2

拟合时间序列模型，要求任何活动特征始终保持活动。

示例

数学细节#

数学上，它由一个线性模型组成，该模型通过混合 \(\ell_1\) \(\ell_2\) 范数进行正则化训练。要最小化的目标函数是：

\[\min_{W} { \frac{1}{2n_{\text{samples}}} ||X W - Y||_{\text{Fro}} ^ 2 + \alpha ||W||_{21}}\]

其中 \(\text{Fro}\) 表示 Frobenius 范数：

\[||A||_{\text{Fro}} = \sqrt{\sum_{ij} a_{ij}^2}\]

且 \(\ell_1\) \(\ell_2\) 表示为：

\[||A||_{2 1} = \sum_i \sqrt{\sum_j a_{ij}^2}.\]

类 MultiTaskLasso 中的实现使用坐标下降作为拟合系数的算法。

1.1.5. Elastic-Net#

ElasticNet 是一种线性回归模型，它通过 \(\ell_1\) 和 \(\ell_2\) 范数对系数进行正则化训练。这种组合允许学习稀疏模型（其中很少有非零权重，如 Lasso），同时仍然保持 Ridge 的正则化特性。我们使用 l1_ratio 参数控制 \(\ell_1\) 和 \(\ell_2\) 的凸组合。

当存在多个相互关联的特征时，Elastic-Net 非常有用。Lasso 可能会随机选择其中一个，而 Elastic-Net 则可能选择两者。

Lasso 和 Ridge 之间权衡的一个实际优势是，它允许 Elastic-Net 继承 Ridge 在旋转下的一些稳定性。

在这种情况下，要最小化的目标函数是：

\[\min_{w} { \frac{1}{2n_{\text{samples}}} ||X w - y||_2 ^ 2 + \alpha \rho ||w||_1 + \frac{\alpha(1-\rho)}{2} ||w||_2 ^ 2}\]
../_images/sphx_glr_plot_lasso_lasso_lars_elasticnet_path_002.png

类 ElasticNetCV 可用于通过交叉验证设置参数 alpha (\(\alpha\)) 和 l1_ratio (\(\rho\))。

示例

参考文献#

以下两篇参考文献解释了 scikit-learn 坐标下降求解器中使用的迭代，以及用于收敛控制的对偶间隙计算。

  • “基于坐标下降的广义线性模型正则化路径”，Friedman, Hastie & Tibshirani, J Stat Softw, 2010 (论文)。

  • “用于大规模 L1 正则化最小二乘的内点法”，S. J. Kim, K. Koh, M. Lustig, S. Boyd and D. Gorinevsky，发表于 IEEE Journal of Selected Topics in Signal Processing, 2007 (论文)

1.1.6. 多任务 Elastic-Net#

MultiTaskElasticNet 是一种 Elastic-Net 模型，用于联合估计多个回归问题（Y 是一个形状为 (n_samples, n_tasks) 的二维数组）的稀疏系数。约束是所有回归问题（也称为任务）选择的特征相同。

数学上，它由一个线性模型组成，该模型通过混合 \(\ell_1\) \(\ell_2\) 范数和 \(\ell_2\) 范数进行正则化训练。要最小化的目标函数是：

\[\min_{W} { \frac{1}{2n_{\text{samples}}} ||X W - Y||_{\text{Fro}}^2 + \alpha \rho ||W||_{2 1} + \frac{\alpha(1-\rho)}{2} ||W||_{\text{Fro}}^2}\]

类 MultiTaskElasticNet 中的实现使用坐标下降作为拟合系数的算法。

类 MultiTaskElasticNetCV 可用于通过交叉验证设置参数 alpha (\(\alpha\)) 和 l1_ratio (\(\rho\))。

1.1.7. 最小角回归#

最小角回归（LARS）是一种用于高维数据的回归算法，由 Bradley Efron、Trevor Hastie、Iain Johnstone 和 Robert Tibshirani 开发。LARS 类似于向前逐步回归。在每一步中，它都会找到与目标最相关的特征。当有多个特征具有相同的相关性时，它不会沿着相同的特征继续，而是沿着特征之间等角的方向前进。

LARS 的优点包括：

  • 在特征数量显著大于样本数量的情况下，它在数值上是高效的。

  • 它的计算速度与向前选择一样快，并且具有与普通最小二乘法相同的复杂度。

  • 它产生完整的 piecewise 线性解路径，这在交叉验证或类似的模型调优尝试中很有用。

  • 如果两个特征与目标的关联度几乎相同，那么它们的系数应该以大致相同的速率增加。因此，该算法的行为符合直觉预期，并且也更稳定。

  • 它可以很容易地修改以产生其他估计器（如 Lasso）的解。

LARS 方法的缺点包括：

  • 由于 LARS 基于对残差的迭代重新拟合，它似乎对噪声的影响特别敏感。Weisberg 在 Efron 等人 (2004) 的《统计年鉴》文章的讨论部分详细讨论了这个问题。

LARS 模型可以通过估计器 Lars，或其低级实现 lars_path 或 lars_path_gram 来使用。

1.1.8. LARS Lasso#

LassoLars 是一个使用 LARS 算法实现的 lasso 模型，与基于坐标下降的实现不同，它产生了精确解，该解是其系数范数的 piecewise 线性函数。

../_images/sphx_glr_plot_lasso_lasso_lars_elasticnet_path_001.png 
>>> from sklearn import linear_model
>>> reg = linear_model.LassoLars(alpha=.1)
>>> reg.fit([[0, 0], [1, 1]], [0, 1])
LassoLars(alpha=0.1)
>>> reg.coef_
array([0.6, 0.        ])

示例

LARS 算法几乎免费提供了沿正则化参数的完整系数路径，因此常见的操作是使用函数 lars_path 或 lars_path_gram 来检索该路径。

数学公式#

该算法类似于向前逐步回归，但不是在每一步都包含特征，而是沿着与每个特征的残差相关性等角的方向增加估计系数。

LARS 解不是给出向量结果，而是由一条曲线组成，该曲线表示参数向量 \(\ell_1\) 范数每个值的解。完整的系数路径存储在形状为 '(n_features, max_features + 1)' 的数组 coef_path_ 中。第一列始终为零。

参考文献

1.1.9. 正交匹配追踪 (OMP)#

OrthogonalMatchingPursuit 和 orthogonal_mp 实现了 OMP 算法，用于在对非零系数数量（即 \(\ell_0\) 伪范数）施加约束的情况下近似线性模型的拟合。

作为一种像最小角回归一样的向前特征选择方法，正交匹配追踪可以近似具有固定数量非零元素的最佳解向量：

\[\underset{w}{\operatorname{arg\,min\,}} ||y - Xw||_2^2 \text{ subject to } ||w||_0 \leq n_{\text{nonzero_coefs}}\]

或者，正交匹配追踪可以针对特定的误差，而不是特定的非零系数数量。这可以表示为：

\[\underset{w}{\operatorname{arg\,min\,}} ||w||_0 \text{ subject to } ||y-Xw||_2^2 \leq \text{tol}\]

OMP 基于一种贪婪算法，该算法在每一步都包含与当前残差高度相关的原子。它类似于更简单的匹配追踪（MP）方法，但更优越之处在于，在每次迭代中，使用对先前选择的字典元素的正交投影重新计算残差。

示例

参考文献#

1.1.10. 贝叶斯回归#

贝叶斯回归技术可用于在估计过程中包含正则化参数：正则化参数不是严格设置，而是根据手头数据进行调整。

这可以通过引入模型超参数的无信息先验来完成。岭回归与分类中使用的 \(\ell_{2}\) 正则化等同于在系数 \(w\) 上带有精度 \(\lambda^{-1}\) 的高斯先验下找到最大后验估计。与其手动设置 lambda，不如将其视为一个随机变量，从数据中进行估计。

为了获得一个完全概率模型，输出 \(y\) 被假定为围绕 \(X w\) 呈高斯分布：

\[p(y|X,w,\alpha) = \mathcal{N}(y|X w,\alpha^{-1})\]

其中 \(\alpha\) 再次被视为要从数据中估计的随机变量。

贝叶斯回归的优点是：

  • 它能适应手头的数据。

  • 它可用于在估计过程中包含正则化参数。

贝叶斯回归的缺点包括：

  • 模型的推断可能非常耗时。

参考文献#

  • 关于贝叶斯方法的一个很好的介绍可以在 C. Bishop 的著作《模式识别与机器学习》中找到。

  • 原始算法详见 Radford M. Neal 的著作 神经网络的贝叶斯学习。

1.1.10.1. 贝叶斯岭回归#

BayesianRidge 估计如上所述的回归问题的概率模型。系数 \(w\) 的先验由球面高斯分布给出：

\[p(w|\lambda) = \mathcal{N}(w|0,\lambda^{-1}\mathbf{I}_{p})\]

\(\alpha\) 和 \(\lambda\) 的先验被选择为伽马分布，这是高斯分布精度的共轭先验。由此产生的模型称为贝叶斯岭回归，它与经典 Ridge 相似。

参数 \(w\)、\(\alpha\) 和 \(\lambda\) 在模型拟合过程中联合估计，正则化参数 \(\alpha\) 和 \(\lambda\) 通过最大化对数边缘似然来估计。scikit-learn 的实现基于 (Tipping, 2001) 附录 A 中描述的算法，其中参数 \(\alpha\) 和 \(\lambda\) 的更新按照 (MacKay, 1992) 的建议进行。最大化过程的初始值可以通过超参数 alpha_init 和 lambda_init 设置。

还有四个超参数，伽马先验分布的 \(\alpha_1\)、\(\alpha_2\)、\(\lambda_1\) 和 \(\lambda_2\)，它们分别作用于 \(\alpha\) 和 \(\lambda\)。这些通常选择为无信息的。默认情况下 \(\alpha_1 = \alpha_2 = \lambda_1 = \lambda_2 = 10^{-6}\)。

贝叶斯岭回归用于回归

>>> from sklearn import linear_model
>>> X = [[0., 0.], [1., 1.], [2., 2.], [3., 3.]]
>>> Y = [0., 1., 2., 3.]
>>> reg = linear_model.BayesianRidge()
>>> reg.fit(X, Y)
BayesianRidge()

拟合后，模型可用于预测新值

>>> reg.predict([[1, 0.]])
array([0.50000013])

可以访问模型的系数 \(w\)

>>> reg.coef_
array([0.49999993, 0.49999993])

由于贝叶斯框架，找到的权重与普通最小二乘法找到的权重略有不同。然而，贝叶斯岭回归对病态问题更具鲁棒性。

示例

参考文献#

1.1.10.2. 自动相关性确定 - ARD#

自动相关性确定（在 ARDRegression 中实现）是一种线性模型，与贝叶斯岭回归非常相似，但会导致更稀疏的系数 \(w\) [1] [2]。

ARDRegression 对 \(w\) 提出了一种不同的先验：它放弃了球面高斯分布，转而使用中心椭圆高斯分布。这意味着每个系数 \(w_{i}\) 本身可以从一个以零为中心、精度为 \(\lambda_{i}\) 的高斯分布中抽取：

\[p(w|\lambda) = \mathcal{N}(w|0,A^{-1})\]

其中 \(A\) 是一个正定对角矩阵，且 \(\text{diag}(A) = \lambda = \{\lambda_{1},...,\lambda_{p}\}\)。

与贝叶斯岭回归不同，\(w_{i}\) 的每个坐标都有自己的标准差 \(\frac{1}{\lambda_i}\)。所有 \(\lambda_i\) 的先验被选择为由超参数 \(\lambda_1\) 和 \(\lambda_2\) 给出的相同伽马分布。

ARD 在文献中也被称为稀疏贝叶斯学习和相关向量机 [3] [4]。有关 ARD 和贝叶斯岭回归之间的详细比较，请参阅下面的示例。

示例

参考文献

1.1.11. 逻辑回归#

逻辑回归在 LogisticRegression 中实现。尽管其名称如此，但在 scikit-learn/ML 术语中，它被实现为用于分类而不是回归的线性模型。逻辑回归在文献中也称为 logit 回归、最大熵分类 (MaxEnt) 或对数线性分类器。在该模型中，描述单个试验可能结果的概率使用逻辑函数进行建模。

此实现可以拟合二元、一对多或多项逻辑回归，并可选地使用 \(\ell_1\)、\(\ell_2\) 或 Elastic-Net 正则化。

注意

正则化

默认应用正则化，这在机器学习中很常见，但在统计学中不常见。正则化的另一个优点是它提高了数值稳定性。不进行正则化相当于将 C 设置为一个非常高的值。

注意

作为广义线性模型 (GLM) 特殊情况的逻辑回归

逻辑回归是广义线性模型的一种特殊情况，它具有二项式/伯努利条件分布和 Logit 链接。逻辑回归的数值输出（即预测概率）可以通过对其应用阈值（默认为 0.5）来用作分类器。scikit-learn 中就是这样实现的，因此它期望一个分类目标，从而使逻辑回归成为一个分类器。

示例

1.1.11.1. 二元情况#

为了方便表示，我们假设数据点 \(i\) 的目标 \(y_i\) 取值在集合 \(\{0, 1\}\) 中。一旦拟合完成，LogisticRegression 的 predict_proba 方法将正类的概率 \(P(y_i=1|X_i)\) 预测为：

\[\hat{p}(X_i) = \operatorname{expit}(X_i w + w_0) = \frac{1}{1 + \exp(-X_i w - w_0)}.\]

作为一个优化问题，带正则化项 \(r(w)\) 的二元逻辑回归最小化以下成本函数：

(1)#\[\min_{w} \frac{1}{S}\sum_{i=1}^n s_i \left(-y_i \log(\hat{p}(X_i)) - (1 - y_i) \log(1 - \hat{p}(X_i))\right) + \frac{r(w)}{S C}\,,\]

其中 \({s_i}\) 对应于用户分配给特定训练样本的权重（向量 \(s\) 由类别权重和样本权重的逐元素乘积形成），且总和 \(S = \sum_{i=1}^n s_i\)。

我们目前通过 penalty 参数为正则化项 \(r(w)\) 提供了四种选择：

惩罚项

\(r(w)\)

无

\(0\)

\(\ell_1\)

\(\|w\|_1\)

\(\ell_2\)

\(\frac{1}{2}\|w\|_2^2 = \frac{1}{2}w^T w\)

ElasticNet

\(\frac{1 - \rho}{2}w^T w + \rho \|w\|_1\)

对于 ElasticNet，\(\rho\)（对应于 l1_ratio 参数）控制 \(\ell_1\) 正则化与 \(\ell_2\) 正则化的强度。当 \(\rho = 1\) 时，Elastic-Net 等价于 \(\ell_1\)；当 \(\rho=0\) 时，等价于 \(\ell_2\)。

请注意，类别权重和样本权重的尺度将影响优化问题。例如，将样本权重乘以常数 \(b>0\) 等同于将（逆）正则化强度 C 乘以 \(b\)。

1.1.11.2. 多项式情况#

二元情况可以扩展到 \(K\) 个类别，从而得到多项逻辑回归，另请参阅对数线性模型。

注意

可以使用仅 \(K-1\) 个权重向量来参数化 \(K\) 类分类模型，利用所有类别概率之和必须为一的事实，使一个类别概率完全由其他类别概率决定。我们特意选择使用 \(K\) 个权重向量来过度参数化模型，以方便实现并保持关于类别顺序的对称归纳偏差，参见 [16]。当使用正则化时，这种影响变得尤为重要。过度参数化的选择对于非惩罚模型可能是有害的，因为那样解可能不是唯一的，如 [16] 所示。

数学细节#

设 \(y_i \in \{1, \ldots, K\}\) 是观测 \(i\) 的标签（序数）编码目标变量。与单个系数向量不同，我们现在有一个系数矩阵 \(W\)，其中每个行向量 \(W_k\) 对应于类别 \(k\)。我们旨在通过 predict_proba 预测类别概率 \(P(y_i=k|X_i)\)，如下：

\[\hat{p}_k(X_i) = \frac{\exp(X_i W_k + W_{0, k})}{\sum_{l=0}^{K-1} \exp(X_i W_l + W_{0, l})}.\]

优化目标变为

\[\min_W -\frac{1}{S}\sum_{i=1}^n \sum_{k=0}^{K-1} s_{ik} [y_i = k] \log(\hat{p}_k(X_i)) + \frac{r(W)}{S C}\,,\]

其中 \([P]\) 表示艾弗森括号（Iverson bracket），如果 \(P\) 为假，则其值为 \(0\)，否则为 \(1\)。

同样，\(s_{ik}\) 是用户分配的权重（样本权重与类别权重的乘积），其总和为 \(S = \sum_{i=1}^n \sum_{k=0}^{K-1} s_{ik}\)。

我们目前通过 penalty 参数提供四种正则化项 \(r(W)\) 的选择，其中 \(m\) 是特征的数量

惩罚项

\(r(W)\)

无

\(0\)

\(\ell_1\)

\(\|W\|_{1,1} = \sum_{i=1}^m\sum_{j=1}^{K}|W_{i,j}|\)

\(\ell_2\)

\(\frac{1}{2}\|W\|_F^2 = \frac{1}{2}\sum_{i=1}^m\sum_{j=1}^{K} W_{i,j}^2\)

ElasticNet

\(\frac{1 - \rho}{2}\|W\|_F^2 + \rho \|W\|_{1,1}\)

1.1.11.3. 求解器#

类 LogisticRegression 中实现的求解器有 “lbfgs”、“liblinear”、“newton-cg”、“newton-cholesky”、“sag” 和 “saga”

下表总结了每个求解器支持的惩罚项和多项式多类别

求解器

惩罚项

‘lbfgs’

‘liblinear’

‘newton-cg’

‘newton-cholesky’

‘sag’

‘saga’

L2 惩罚

是

是

是

是

是

是

L1 惩罚

否

是

否

否

否

是

Elastic-Net (L1 + L2)

否

否

否

否

否

是

无惩罚 (‘none’)

是

否

是

是

是

是

多类别支持

多项式多类别

是

否

是

是

是

是

行为

惩罚截距 (不好)

否

是

否

否

否

否

对于大型数据集更快

否

否

否

否

是

是

对未缩放的数据集具有鲁棒性

是

是

是

是

否

否

“lbfgs” 求解器因其鲁棒性而被默认使用。对于 n_samples >> n_features 的情况，“newton-cholesky” 是一个不错的选择，并且可以达到高精度（极小的 tol 值）。对于大型数据集，“saga” 求解器通常更快（比“lbfgs”快），特别是对于低精度（高 tol）。对于大型数据集，您也可以考虑使用 SGDClassifier 并设置 loss="log_loss"，它可能更快，但需要更多的调优。

1.1.11.3.1. 求解器之间的差异#

当 fit_intercept=False 且拟合的 coef_（或）要预测的数据为零时，使用 solver=liblinear 的 LogisticRegression 或 LinearSVC 与直接使用外部 liblinear 库之间获得的分数可能存在差异。这是因为对于 decision_function 为零的样本，LogisticRegression 和 LinearSVC 预测负类别，而 liblinear 预测正类别。请注意，一个 fit_intercept=False 且有许多 decision_function 为零的样本的模型，很可能是欠拟合的坏模型，建议您设置 fit_intercept=True 并增加 intercept_scaling。

求解器详情#

  • “liblinear” 求解器使用坐标下降 (CD) 算法，并依赖于随 scikit-learn 提供的优秀 C++ LIBLINEAR 库。然而，liblinear 中实现的 CD 算法无法学习真正的多项式（多类别）模型；相反，优化问题被分解为“一对多”的方式，因此为所有类别训练了单独的二元分类器。这在内部发生，因此使用此求解器的 LogisticRegression 实例表现为多类别分类器。对于 \(\ell_1\) 正则化，sklearn.svm.l1_min_c 允许计算 C 的下限，以便获得非“空”（所有特征权重为零）模型。

  • “lbfgs”、“newton-cg”和“sag”求解器仅支持 \(\ell_2\) 正则化或无正则化，并且对于某些高维数据被发现收敛更快。将这些求解器的 multi_class 设置为“multinomial”会学习一个真正的多项式逻辑回归模型[5]，这意味着其概率估计应该比默认的“一对多”设置更好地校准。

  • “sag” 求解器使用随机平均梯度下降 [6]。当样本数量和特征数量都很大时，它比其他求解器对于大型数据集更快。

  • “saga” 求解器 [7] 是“sag”的一个变体，它也支持非平滑的 penalty="l1"。因此，它是稀疏多项式逻辑回归的首选求解器。它也是唯一支持 penalty="elasticnet" 的求解器。

  • “lbfgs” 是一种优化算法，它近似于 Broyden–Fletcher–Goldfarb–Shanno 算法 [8]，属于拟牛顿方法。因此，它可以处理各种不同的训练数据，是默认的求解器。然而，其性能在缩放不良的数据集和具有稀有类别的独热编码分类特征的数据集上会受到影响。

  • “newton-cholesky” 求解器是一个精确的牛顿求解器，它计算 Hessian 矩阵并求解所得的线性系统。对于 n_samples >> n_features 的情况，它是一个非常好的选择，并且可以达到高精度（极小的 tol 值），但也有一些缺点：仅支持 \(\ell_2\) 正则化。此外，由于 Hessian 矩阵是显式计算的，内存使用量与 n_features 以及 n_classes 呈二次方关系。

有关其中一些求解器的比较，请参阅 [9]。

参考文献

注意

使用稀疏逻辑回归进行特征选择

带有 \(\ell_1\) 惩罚的逻辑回归会产生稀疏模型，因此可以用于执行特征选择，如 基于 L1 的特征选择 中所述。

注意

P 值估计

在无惩罚回归的情况下，可以获得系数的 P 值和置信区间。statsmodels 包 原生支持此功能。在 sklearn 中，也可以使用自举法代替。

LogisticRegressionCV 实现了带有内置交叉验证支持的逻辑回归，以根据 scoring 属性找到最优的 C 和 l1_ratio 参数。“newton-cg”、“sag”、“saga”和“lbfgs”求解器被发现对于高维密集数据更快，这是由于热启动（请参阅 术语表）。

1.1.12. 广义线性模型#

广义线性模型 (GLM) 以两种方式扩展了线性模型[10]。首先，预测值 \(\hat{y}\) 通过逆链接函数 \(h\) 与输入变量 \(X\) 的线性组合关联，如下所示：

\[\hat{y}(w, X) = h(Xw).\]

其次，平方损失函数被指数族分布（或更准确地说，可再生的指数离散模型 (EDM)[11]）的单位偏差 \(d\) 所取代。

最小化问题变为

\[\min_{w} \frac{1}{2 n_{\text{samples}}} \sum_i d(y_i, \hat{y}_i) + \frac{\alpha}{2} ||w||_2^2,\]

其中 \(\alpha\) 是 L2 正则化惩罚。当提供样本权重时，平均值变为加权平均值。

下表列出了一些特定的 EDM 及其单位偏差

分布

目标域

单位偏差 \(d(y, \hat{y})\)

正态分布

\(y \in (-\infty, \infty)\)

\((y-\hat{y})^2\)

伯努利分布

\(y \in \{0, 1\}\)

\(2({y}\log\frac{y}{\hat{y}}+({1}-{y})\log\frac{{1}-{y}}{{1}-\hat{y}})\)

分类分布

\(y \in \{0, 1, ..., k\}\)

\(2\sum_{i \in \{0, 1, ..., k\}} I(y = i) y_\text{i}\log\frac{I(y = i)}{\hat{I(y = i)}}\)

泊松分布

\(y \in [0, \infty)\)

\(2(y\log\frac{y}{\hat{y}}-y+\hat{y})\)

伽马分布

\(y \in (0, \infty)\)

\(2(\log\frac{\hat{y}}{y}+\frac{y}{\hat{y}}-1)\)

逆高斯分布

\(y \in (0, \infty)\)

\(\frac{(y-\hat{y})^2}{y\hat{y}^2}\)

这些分布的概率密度函数 (PDF) 如下图所示，

../_images/poisson_gamma_tweedie_distributions.png

服从泊松分布、Tweedie (power=1.5) 和伽马分布的随机变量 Y 的 PDF，具有不同的均值 (\(\mu\))。注意泊松分布和 Tweedie (power=1.5) 分布在 \(Y=0\) 处存在点质量，而伽马分布没有，因为它具有严格的正目标域。#

伯努利分布是一种离散概率分布，模拟伯努利试验——一个只有两个互斥结果的事件。分类分布是伯努利分布的推广，用于分类随机变量。伯努利分布中的随机变量有两种可能的结果，而分类随机变量可以取 K 个可能类别中的一个，每个类别的概率单独指定。

分布的选择取决于手头的问题

  • 如果目标值 \(y\) 是计数（非负整数值）或相对频率（非负），您可以使用带有对数链接的泊松分布。

  • 如果目标值为正且偏斜，您可以尝试使用带有对数链接的伽马分布。

  • 如果目标值的尾部比伽马分布更重，您可以尝试逆高斯分布（甚至 Tweedie 族中更高方差的幂）。

  • 如果目标值 \(y\) 是概率，您可以使用伯努利分布。带有 logit 链接的伯努利分布可用于二元分类。带有 softmax 链接的分类分布可用于多类别分类。

用例示例#

  • 农业/天气建模：每年降雨事件次数（泊松），每次事件降雨量（伽马），每年总降雨量（Tweedie / 复合泊松伽马）。

  • 风险建模/保险保单定价：每年每保单持有人索赔事件次数（泊松），每次事件成本（伽马），每年每保单持有人总成本（Tweedie / 复合泊松伽马）。

  • 信用违约：贷款无法偿还的概率（伯努利）。

  • 欺诈检测：现金转账等金融交易是欺诈交易的概率（伯努利）。

  • 预测性维护：每年生产中断事件次数（泊松），中断持续时间（伽马），每年总中断时间（Tweedie / 复合泊松伽马）。

  • 药物测试：在系列试验中治愈患者的概率或患者出现副作用的概率（伯努利）。

  • 新闻分类：将新闻文章分类为商业新闻、政治和娱乐新闻三类（分类）。

参考文献

1.1.12.1. 用法#

TweedieRegressor 实现了 Tweedie 分布的广义线性模型，可以使用适当的 power 参数来建模上述任何分布。特别地，

  • power = 0：正态分布。在这种情况下，诸如 Ridge、ElasticNet 等特定估计器通常更合适。

  • power = 1：泊松分布。PoissonRegressor 为方便而暴露。然而，它严格等同于 TweedieRegressor(power=1, link='log')。

  • power = 2：伽马分布。GammaRegressor 为方便而暴露。然而，它严格等同于 TweedieRegressor(power=2, link='log')。

  • power = 3：逆高斯分布。

链接函数由 link 参数决定。

使用示例

>>> from sklearn.linear_model import TweedieRegressor
>>> reg = TweedieRegressor(power=1, alpha=0.5, link='log')
>>> reg.fit([[0, 0], [0, 1], [2, 2]], [0, 1, 2])
TweedieRegressor(alpha=0.5, link='log', power=1)
>>> reg.coef_
array([0.2463, 0.4337])
>>> reg.intercept_
np.float64(-0.7638)

示例

实际考虑#

特征矩阵 X 在拟合之前应进行标准化。这确保了惩罚对特征一视同仁。

由于线性预测器 \(Xw\) 可以为负值，而泊松、伽马和逆高斯分布不支持负值，因此需要应用一个保证非负性的逆链接函数。例如，当 link='log' 时，逆链接函数变为 \(h(Xw)=\exp(Xw)\)。

如果您想建模相对频率，即每暴露（时间、体积等）的计数，您可以通过使用泊松分布并将 \(y=\frac{\mathrm{counts}}{\mathrm{exposure}}\) 作为目标值以及 \(\mathrm{exposure}\) 作为样本权重来实现。具体示例请参阅 保险索赔的 Tweedie 回归。

在对 TweedieRegressor 的 power 参数进行交叉验证时，建议指定一个显式的 scoring 函数，因为默认评分器 TweedieRegressor.score 本身是 power 的一个函数。

1.1.13. 随机梯度下降 - SGD#

随机梯度下降是一种简单但非常有效的线性模型拟合方法。当样本数量（和特征数量）非常大时，它尤其有用。partial_fit 方法允许在线/核外学习。

类 SGDClassifier 和 SGDRegressor 提供了使用不同（凸）损失函数和不同惩罚来拟合分类和回归线性模型的功能。例如，当 loss="log" 时，SGDClassifier 拟合逻辑回归模型，而当 loss="hinge" 时，它拟合线性支持向量机 (SVM)。

您可以参考专门的 随机梯度下降 文档部分获取更多详细信息。

1.1.14. 感知机#

Perceptron 是另一种适用于大规模学习的简单分类算法。默认情况下：

  • 它不需要学习率。

  • 它未被正则化（惩罚）。

  • 它只在出现错误时更新模型。

最后一个特性意味着感知机比使用铰链损失的 SGD 训练速度略快，并且生成的模型更稀疏。

实际上，Perceptron 是 SGDClassifier 类的一个封装，它使用感知机损失和恒定学习率。有关更多详细信息，请参阅 SGD 过程的 数学部分。

1.1.15. 被动攻击算法#

被动攻击算法是一系列用于大规模学习的算法。它们与感知机类似，不需要学习率。然而，与感知机不同的是，它们包含一个正则化参数 C。

对于分类，PassiveAggressiveClassifier 可与 loss='hinge' (PA-I) 或 loss='squared_hinge' (PA-II) 一起使用。对于回归，PassiveAggressiveRegressor 可与 loss='epsilon_insensitive' (PA-I) 或 loss='squared_epsilon_insensitive' (PA-II) 一起使用。

参考文献#

1.1.16. 鲁棒回归：异常值和建模误差#

鲁棒回归旨在在存在损坏数据（异常值或模型误差）的情况下拟合回归模型。

../_images/sphx_glr_plot_theilsen_001.png

1.1.16.1. 不同场景和有用概念#

处理被异常值污染的数据时需要注意以下几点：

  • X 或 y 中的异常值?

    y 方向的异常值

    X 方向的异常值

    y_outliers

    X_outliers

  • 异常值比例与误差幅度

    异常点的数量很重要，但它们的异常程度也同样重要。

    小异常值

    大异常值

    y_outliers

    large_y_outliers

鲁棒拟合的一个重要概念是**崩溃点**：拟合开始偏离内点数据时，异常数据的比例。

请注意，通常在高维设置（大 n_features）中进行鲁棒拟合非常困难。此处的鲁棒模型可能不适用于这些设置。

1.1.16.2. RANSAC: 随机样本一致性#

RANSAC (RANdom SAmple Consensus) 从完整数据集的内点随机子集中拟合模型。

RANSAC 是一种非确定性算法，只能以一定的概率产生合理的结果，该概率取决于迭代次数（参见 max_trials 参数）。它通常用于线性和非线性回归问题，在摄影测量计算机视觉领域尤为流行。

该算法将完整的输入样本数据分为一组内点（可能受到噪声影响）和异常值（例如由错误测量或对数据无效假设引起）。然后仅根据确定的内点估计结果模型。

../_images/sphx_glr_plot_ransac_001.png

示例

算法详情#

每次迭代执行以下步骤：

  1. 从原始数据中选择 min_samples 个随机样本，并检查数据是否有效（参见 is_data_valid）。

  2. 将模型拟合到随机子集（estimator.fit），并检查估计模型是否有效（参见 is_model_valid）。

  3. 通过计算与估计模型的残差（estimator.predict(X) - y）将所有数据分类为内点或异常值——所有绝对残差小于或等于 residual_threshold 的数据样本都被视为内点。

  4. 如果内点样本数量最大，则将拟合模型保存为最佳模型。如果当前估计模型具有相同数量的内点，则只有当它具有更好的分数时才被视为最佳模型。

这些步骤将执行最多 max_trials 次，或者直到满足其中一个特殊停止条件（参见 stop_n_inliers 和 stop_score）。最终模型是使用先前确定的最佳模型的所有内点样本（共识集）估计的。

is_data_valid 和 is_model_valid 函数允许识别和拒绝随机子样本的退化组合。如果不需要估计模型来识别退化情况，则应使用 is_data_valid，因为它在拟合模型之前被调用，从而带来更好的计算性能。

参考文献#

1.1.16.3. Theil-Sen 估计器：基于广义中位数估计器#

TheilSenRegressor 估计器使用多维中位数的一种推广。因此，它对多元异常值具有鲁棒性。但请注意，估计器的鲁棒性随问题维度快速下降。在高维情况下，它会失去其鲁棒性，并且不会比普通最小二乘法更好。

示例

理论考量#

TheilSenRegressor 在渐近效率和作为无偏估计器方面与 普通最小二乘法 (OLS) 相当。与 OLS 相比，Theil-Sen 是一种非参数方法，这意味着它不对数据的基础分布做任何假设。由于 Theil-Sen 是一种基于中位数的估计器，因此它对损坏的数据（即异常值）更具鲁棒性。在单变量设置中，对于简单的线性回归，Theil-Sen 的崩溃点约为 29.3%，这意味着它可以容忍高达 29.3% 的任意损坏数据。

../_images/sphx_glr_plot_theilsen_001.png

scikit-learn 中 TheilSenRegressor 的实现遵循对多元线性回归模型的推广 [14]，该模型使用空间中位数，它是中位数在多维中的推广 [15]。

在时间和空间复杂度方面，Theil-Sen 的缩放比例为

\[\binom{n_{\text{samples}}}{n_{\text{subsamples}}}\]

这使得它无法穷举地应用于样本和特征数量庞大的问题。因此，可以通过仅考虑所有可能组合的随机子集来选择子群体的大小，以限制时间和空间复杂度。

参考文献

另请参阅 Wikipedia 页面

1.1.16.4. Huber 回归#

HuberRegressor 与 Ridge 不同，因为它对由 epsilon 参数定义为异常值的样本应用线性损失。如果样本的绝对误差小于阈值 epsilon，则该样本被分类为内点。它与 TheilSenRegressor 和 RANSACRegressor 的区别在于它不忽略异常值的影响，而是赋予它们较小的权重。

../_images/sphx_glr_plot_huber_vs_ridge_001.png

示例

数学细节#

HuberRegressor 最小化

\[\min_{w, \sigma} {\sum_{i=1}^n\left(\sigma + H_{\epsilon}\left(\frac{X_{i}w - y_{i}}{\sigma}\right)\sigma\right) + \alpha {||w||_2}^2}\]

其中损失函数由下式给出

\[\begin{split}H_{\epsilon}(z) = \begin{cases} z^2, & \text {if } |z| < \epsilon, \\ 2\epsilon|z| - \epsilon^2, & \text{otherwise} \end{cases}\end{split}\]

建议将参数 epsilon 设置为 1.35 以达到 95% 的统计效率。

参考文献

  • Peter J. Huber, Elvezio M. Ronchetti: Robust Statistics, Concomitant scale estimates, p. 172.

HuberRegressor 与使用损失设置为 huber 的 SGDRegressor 的区别在于以下几个方面。

  • HuberRegressor 具有尺度不变性。一旦设置了 epsilon，即使以不同的值缩放 X 和 y，也能产生与以前相同的对异常值的鲁棒性。而 SGDRegressor 则在缩放 X 和 y 后需要重新设置 epsilon。

  • HuberRegressor 在样本数量较少的数据上使用效率更高，而 SGDRegressor 需要在训练数据上进行多次迭代才能达到相同的鲁棒性。

请注意，此估计器与 R 中鲁棒回归的实现 不同，因为 R 实现基于残差大于某个阈值的程度对每个样本赋予权重，从而进行加权最小二乘实现。

1.1.17. 分位数回归#

分位数回归估计以 \(X\) 为条件的 \(y\) 的中位数或其他分位数，而普通最小二乘法 (OLS) 估计条件均值。

如果您对预测区间而不是点预测感兴趣，分位数回归可能会很有用。有时，预测区间是基于预测误差服从零均值和恒定方差的正态分布的假设来计算的。分位数回归即使对于具有非恒定（但可预测）方差或非正态分布的误差，也能提供合理的预测区间。

../_images/sphx_glr_plot_quantile_regression_002.png

基于最小化弹子球损失，条件分位数也可以通过线性模型以外的模型来估计。例如，如果将其参数 loss 设置为 "quantile" 且参数 alpha 设置为应预测的分位数，则 GradientBoostingRegressor 可以预测条件分位数。请参阅 梯度提升回归的预测区间 中的示例。

大多数分位数回归的实现都基于线性规划问题。当前的实现基于 scipy.optimize.linprog。

示例

数学细节#

作为线性模型，QuantileRegressor 为第 \(q\) 个分位数（\(q \in (0, 1)\)）给出线性预测 \(\hat{y}(w, X) = Xw\)。然后通过以下最小化问题找到权重或系数 \(w\)：

\[\min_{w} {\frac{1}{n_{\text{samples}}} \sum_i PB_q(y_i - X_i w) + \alpha ||w||_1}.\]

这包括弹子球损失（也称为线性损失），另请参阅 mean_pinball_loss，

\[\begin{split}PB_q(t) = q \max(t, 0) + (1 - q) \max(-t, 0) = \begin{cases} q t, & t > 0, \\ 0, & t = 0, \\ (q-1) t, & t < 0 \end{cases}\end{split}\]

以及由参数 alpha 控制的 L1 惩罚，类似于 Lasso。

由于弹子球损失仅在残差上呈线性，分位数回归比基于平方误差的均值估计对异常值更具鲁棒性。介于两者之间的是 HuberRegressor。

参考文献#

1.1.18. 多项式回归：用基函数扩展线性模型#

机器学习中一个常见的模式是使用在数据的非线性函数上训练的线性模型。这种方法保持了线性模型通常较快的性能，同时允许它们拟合更广泛的数据。

数学细节#

例如，通过从系数构建**多项式特征**，可以扩展简单的线性回归。在标准线性回归情况下，对于二维数据，您可能有一个模型如下所示：

\[\hat{y}(w, x) = w_0 + w_1 x_1 + w_2 x_2\]

如果我们要拟合一个抛物面而不是平面到数据，我们可以将特征组合成二阶多项式，使模型看起来像这样：

\[\hat{y}(w, x) = w_0 + w_1 x_1 + w_2 x_2 + w_3 x_1 x_2 + w_4 x_1^2 + w_5 x_2^2\]

（有时令人惊讶的）观察结果是，这**仍然是一个线性模型**：要理解这一点，想象创建一个新的特征集

\[z = [x_1, x_2, x_1 x_2, x_1^2, x_2^2]\]

通过对数据的重新标记，我们的问题可以写成

\[\hat{y}(w, z) = w_0 + w_1 z_1 + w_2 z_2 + w_3 z_3 + w_4 z_4 + w_5 z_5\]

我们看到，由此产生的**多项式回归**属于我们上面考虑的线性模型类别（即模型在 \(w\) 中是线性的），并且可以使用相同的技术解决。通过在由这些基函数构建的更高维空间中考虑线性拟合，模型具有灵活性，可以拟合更广泛的数据范围。

这是一个将此思想应用于一维数据的示例，使用不同程度的多项式特征

../_images/sphx_glr_plot_polynomial_interpolation_001.png

此图是使用 PolynomialFeatures 转换器创建的，该转换器将输入数据矩阵转换为给定度数的新数据矩阵。它可以按如下方式使用：

>>> from sklearn.preprocessing import PolynomialFeatures
>>> import numpy as np
>>> X = np.arange(6).reshape(3, 2)
>>> X
array([[0, 1],
       [2, 3],
       [4, 5]])
>>> poly = PolynomialFeatures(degree=2)
>>> poly.fit_transform(X)
array([[ 1.,  0.,  1.,  0.,  0.,  1.],
       [ 1.,  2.,  3.,  4.,  6.,  9.],
       [ 1.,  4.,  5., 16., 20., 25.]])

X 的特征已从 \([x_1, x_2]\) 转换为 \([1, x_1, x_2, x_1^2, x_1 x_2, x_2^2]\)，现在可以在任何线性模型中使用。

这种预处理可以使用 Pipeline 工具进行简化。可以创建一个表示简单多项式回归的单一对象并按如下方式使用：

>>> from sklearn.preprocessing import PolynomialFeatures
>>> from sklearn.linear_model import LinearRegression
>>> from sklearn.pipeline import Pipeline
>>> import numpy as np
>>> model = Pipeline([('poly', PolynomialFeatures(degree=3)),
...                   ('linear', LinearRegression(fit_intercept=False))])
>>> # fit to an order-3 polynomial data
>>> x = np.arange(5)
>>> y = 3 - 2 * x + x ** 2 - x ** 3
>>> model = model.fit(x[:, np.newaxis], y)
>>> model.named_steps['linear'].coef_
array([ 3., -2.,  1., -1.])

在多项式特征上训练的线性模型能够精确地恢复输入多项式系数。

在某些情况下，不需要包含任何单个特征的更高次幂，而只需包含最多将 \(d\) 个不同特征相乘的所谓*交互特征*。这些可以通过 PolynomialFeatures 并设置 interaction_only=True 获得。

例如，当处理布尔特征时，对于所有 \(n\)，\(x_i^n = x_i\)，因此无用；但 \(x_i x_j\) 表示两个布尔值的合取。通过这种方式，我们可以使用线性分类器解决 XOR 问题：

>>> from sklearn.linear_model import Perceptron
>>> from sklearn.preprocessing import PolynomialFeatures
>>> import numpy as np
>>> X = np.array([[0, 0], [0, 1], [1, 0], [1, 1]])
>>> y = X[:, 0] ^ X[:, 1]
>>> y
array([0, 1, 1, 0])
>>> X = PolynomialFeatures(interaction_only=True).fit_transform(X).astype(int)
>>> X
array([[1, 0, 0, 0],
       [1, 0, 1, 0],
       [1, 1, 0, 0],
       [1, 1, 1, 1]])
>>> clf = Perceptron(fit_intercept=False, max_iter=10, tol=None,
...                  shuffle=False).fit(X, y)

并且分类器“预测”是完美的

>>> clf.predict(X)
array([0, 1, 1, 0])
>>> clf.score(X, y)
1.0