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\]

LinearRegression 将在其fit方法中接收数组X，y，并将线性模型的系数\(w\)存储在其coef_成员中。

>>> 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])

普通最小二乘法的系数估计依赖于特征的独立性。当特征相关且设计矩阵\(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\)的值越大，收缩量越大，因此系数对共线性更稳健。

与其他线性模型一样，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_
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\)。

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

示例

• 绘制 Ridge 系数作为正则化的函数

• 使用稀疏特征对文本文档进行分类

• 线性模型系数解释中的常见陷阱

1.1.2.3. Ridge 复杂度

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

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_
0.01

指定cv属性的值将触发使用带有GridSearchCV的交叉验证，例如cv=10表示10倍交叉验证，而不是留一交叉验证。

参考文献

[RL2007]

“关于正则化最小二乘法的笔记”，Rifkin & Lippert（技术报告，课程幻灯片）。

1.1.3. Lasso

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 对于较低级别的任务很有用，因为它计算了所有可能值的完整路径上的系数。

示例

• 用于稀疏信号的基于 L1 的模型

• 压缩感知：使用 L1 先验的断层扫描重建 (Lasso)

• 线性模型系数解释中的常见陷阱

注意

使用 Lasso 进行特征选择

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

参考文献

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

• “通过坐标下降进行广义线性模型的正则化路径”，Friedman、Hastie 和 Tibshirani，J Stat Softw，2010 (论文）。

• “一种用于大规模 L1 正则化最小二乘法的内点法”，S. J. Kim、K. Koh、M. Lustig、S. Boyd 和 D. Gorinevsky，IEEE 信号处理精选主题期刊，2007 年 (论文)

1.1.3.1. 设置正则化参数

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

1.1.3.1.1. 使用交叉验证

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

对于具有许多共线性特征的高维数据集，LassoCV 通常更可取。但是，LassoLarsCV 的优势在于它可以探索更多相关的 alpha 参数值，并且如果样本数量与特征数量相比非常小，它通常比 LassoCV 更快。

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

或者，估计器 LassoLarsIC 建议使用 Akaike 信息准则 (AIC) 和 Bayes 信息准则 (BIC)。它是一种计算成本更低的替代方法，用于找到 alpha 的最佳值，因为正则化路径仅计算一次，而不是使用 k 折交叉验证时计算 k+1 次。

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

然而，此类标准需要对解的自由度进行适当的估计，是针对大样本（渐近结果）导出的，并假设正确的模型是正在研究的候选模型。当问题条件很差时（例如，特征数量多于样本数量），它们也往往会失效。

示例

• Lasso 模型选择：AIC-BIC/交叉验证

• 通过信息准则进行 Lasso 模型选择

1.1.3.1.3. AIC 和 BIC 准则

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

数学细节

AIC 准则定义为

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

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

BIC 的定义是用 \(\log(N)\) 替换常数 \(2\)。

\[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 时有效。

参考文献

[12]

Zou, Hui, Trevor Hastie, and Robert Tibshirani. “On the degrees of freedom of the lasso.” The Annals of Statistics 35.5 (2007): 2173-2192.

[13]

Cherkassky, Vladimir, and Yunqian Ma. “Comparison of model selection for regression.” Neural computation 15.7 (2003): 1691-1714.

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 的非零值是完整的列。

拟合时间序列模型，强制任何活动特征在所有时间点都处于活动状态。

示例

• 使用多任务 Lasso 进行联合特征选择

数学细节

在数学上，它由使用混合 \(\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}\]

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

示例

• 用于稀疏信号的基于 L1 的模型

• Lasso、Lasso-LARS 和 Elastic Net 轨迹

• 使用预计算的 Gram 矩阵和加权样本拟合 Elastic Net

参考文献

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

• “通过坐标下降进行广义线性模型的正则化路径”，Friedman、Hastie 和 Tibshirani，J Stat Softw，2010 (论文）。

• “一种用于大规模 L1 正则化最小二乘法的内点法”，S. J. Kim、K. Koh、M. Lustig、S. Boyd 和 D. Gorinevsky，IEEE 信号处理精选主题期刊，2007 年 (论文)

1.1.6. 多任务 Elastic-Net

MultiTaskElasticNet 是一个弹性网络模型，它可以联合估计多个回归问题的稀疏系数：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 的优点在于：

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

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

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

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

• 它很容易修改以产生其他估计量的解，例如 Lasso。

LARS 方法的缺点包括：

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

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

1.1.8. LARS Lasso

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

>>> 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.        ])

示例

• Lasso、Lasso-LARS 和 Elastic Net 轨迹

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

数学公式

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

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

参考文献

• 原始算法在 Hastie 等人的论文 最小角回归 中有详细介绍。

1.1.9. 正交匹配追踪 (OMP)

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

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

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

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

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

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

示例

• 正交匹配追踪

参考文献

• https://www.cs.technion.ac.il/~ronrubin/Publications/KSVD-OMP-v2.pdf

• 使用时频字典的匹配追踪，S. G. Mallat，Z. Zhang，

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 的《神经网络的贝叶斯学习》Bayesian learning for neural networks一书中有详细描述。

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])

由于贝叶斯框架，找到的权重与普通最小二乘法找到的权重略有不同。但是，贝叶斯岭回归对不适定问题更稳健。

示例

• 使用贝叶斯岭回归进行曲线拟合

参考文献

• Christopher M. Bishop：《模式识别与机器学习》，2006 年，第 3.3 节。

• David J. C. MacKay，贝叶斯插值，1992 年。

• Michael E. Tipping，稀疏贝叶斯学习和相关向量机，2001 年。

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]

Christopher M. Bishop: 模式识别与机器学习，第 7.2.1 章

[2]

David Wipf 和 Srikantan Nagarajan: 自动相关性确定的一种新视角

[3]

Michael E. Tipping: 稀疏贝叶斯学习和相关向量机

[4]

Tristan Fletcher: 相关向量机详解

1.1.11. 逻辑回归

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

此实现可以拟合二元、一对多或多项逻辑回归，并可以选择\(\ell_1\)、\(\ell_2\)或弹性网络正则化。

注意

正则化

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

注意

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

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

示例

• L1 惩罚和逻辑回归中的稀疏性

• L1-逻辑回归的正则化路径

• 多项式和一对多逻辑回归的决策边界

• 在 20newgroups 上进行多类别稀疏逻辑回归

• 使用多项式逻辑回归 + L1 进行 MNIST 分类

• 绘制分类概率

1.1.11.1. 二元情况

为方便起见，我们假设目标\(y_i\)对于数据点\(i\)取值于集合\(\{0, 1\}\)。拟合后，predict_proba方法LogisticRegression预测正类概率\(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)\)的四种选择

penalty	\(r(w)\)
None	\(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\)时，Elastic-Net 等效于\(\ell_2\)。

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

1.1.11.2. 多项式情况

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

注意

可以使用仅\(K-1\)个权重向量来参数化\(K\)类分类模型，通过利用所有类概率必须加起来等于1的事实，留下一个类概率完全由其他类概率确定。我们故意选择使用\(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]\)表示艾弗森括号，如果\(P\)为假，则其值为\(0\)，否则其值为\(1\)。

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

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

penalty	\(r(W)\)
None	\(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惩罚	否	是		否	否	否	是
弹性网络 (L1 + L2)	否	否		否	否	否	是
无惩罚 (‘none’)	是	否		是	是	是	是
多类别支持
多项式多类别	是	否		是	否	是	是
行为
惩罚截距 (不好)	否	是		否	否	否	否
对于大型数据集更快	否	否		否	否	是	是
对未缩放的数据集鲁棒	是	是		是	是	否	否

默认情况下使用“lbfgs”求解器，因为它具有鲁棒性。对于大型数据集，“saga”求解器通常更快。对于大型数据集，您也可以考虑使用SGDClassifier，其中loss="log_loss"，这可能更快，但需要更多调整。

1.1.11.3.1. 求解器之间的差异

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

求解器详情

• 求解器“liblinear”使用坐标下降 (CD) 算法，并依赖于优秀的 C++ LIBLINEAR 库，该库随 scikit-learn 一起提供。但是，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，这是一个非常好的选择，但它有一些缺点：仅支持 \(\ell_2\) 正则化。此外，由于显式计算 Hessian 矩阵，内存使用量对 n_features 和 n_classes 都有二次依赖性。因此，多类别情况下只实现了“一对多”方案。

有关这些求解器的一些比较，请参见 [9]。

参考文献

[5]

Christopher M. Bishop：模式识别和机器学习，第 4.3.4 章

[6]

Mark Schmidt、Nicolas Le Roux 和 Francis Bach：使用随机平均梯度最小化有限和。

[7]

Aaron Defazio、Francis Bach、Simon Lacoste-Julien：SAGA：一种支持非强凸组合目标的快速增量梯度方法。

[8]

https://en.wikipedia.org/wiki/Broyden%E2%80%93Fletcher%E2%80%93Goldfarb%E2%80%93Shanno_algorithm

[9]

Thomas P. Minka “逻辑回归的数值优化器比较”

[16] (1,2)

Simon、Noah、J. Friedman 和 T. Hastie。“一种用于组惩罚多响应和多项式回归的分块下降算法。”

注意

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

具有 \(\ell_1\) 惩罚的逻辑回归会产生稀疏模型，因此可用于执行特征选择，详见 基于 L1 的特征选择。

注意

P 值估计

在无惩罚回归的情况下，可以获得系数的 p 值和置信区间。statsmodels 包 原生支持这一点。在 sklearn 中，也可以使用 bootstrapping。

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).\]

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

最小化问题变为

\[\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) 在下图中进行了说明，

服从泊松分布、Tweedie（幂=1.5）分布和伽马分布的不同均值 (\(\mu\)) 的随机变量 Y 的 PDF。观察泊松分布和 Tweedie（幂=1.5）分布在 \(Y=0\) 处的点质量，但伽马分布没有，因为伽马分布具有严格为正的目标域。

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

分布的选择取决于具体的问题。

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

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

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

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

用例示例

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

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

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

• 欺诈检测：金融交易（如现金转账）属于欺诈交易的概率（伯努利分布）。

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

• 医学药物测试：在一组试验中治愈患者的概率或患者出现副作用的概率（伯努利分布）。

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

参考文献

[10]

McCullagh, Peter; Nelder, John (1989). Generalized Linear Models, Second Edition. Boca Raton: Chapman and Hall/CRC. ISBN 0-412-31760-5.

[11]

Jørgensen, B. (1992). The theory of exponential dispersion models and analysis of deviance. Monografias de matemática, no. 51. 另见 指数弥散模型。

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_
-0.7638...

示例

• 泊松回归和非正态损失

• 保险索赔的 Tweedie 回归

实际注意事项

特征矩阵 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）。

参考文献

• “在线被动攻击算法” K. Crammer, O. Dekel, J. Keshat, S. Shalev-Shwartz, Y. Singer - JMLR 7 (2006)

1.1.16. 稳健回归：异常值和建模误差

稳健回归旨在拟合存在数据损坏的回归模型：异常值或模型误差。

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

处理受异常值损坏的数据时，需要记住以下几点

• X或y中的异常值?

  y方向的异常值	X方向的异常值

• 异常值的比例与误差幅度

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

  小的异常值	大的异常值

稳健拟合的一个重要概念是崩溃点：数据中可以偏离多少比例，拟合才开始错过内点数据。

请注意，通常情况下，高维设置（大的n_features）中的稳健拟合非常困难。这里的稳健模型可能在这种情况下不起作用。

权衡：哪个估计器？

Scikit-learn 提供了 3 种稳健回归估计器：RANSAC、Theil Sen 和HuberRegressor。

• HuberRegressor 应该比 RANSAC 和 Theil Sen 快，除非样本数量非常大，即 n_samples >> n_features。这是因为 RANSAC 和 Theil Sen 拟合的是数据的较小子集。但是，对于默认参数，Theil Sen 和 RANSAC 可能不如 HuberRegressor 稳健。

• RANSAC 算法比 Theil-Sen 算法更快，并且随着样本数量的增加，其可扩展性也更好。

• RANSAC 算法能更好地处理 y 方向上的大型异常值（最常见的情况）。

• Theil-Sen 算法能更好地处理 X 方向上的中等大小的异常值，但在高维情况下，这种特性会消失。

如有疑问，请使用 RANSAC 算法。

1.1.16.2. RANSAC：随机采样一致性

RANSAC（RANdom SAmple Consensus，随机采样一致性）算法从完整数据集的内点随机子集中拟合模型。

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

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

示例

• 使用 RANSAC 进行稳健线性模型估计

• 稳健线性估计器拟合

算法细节

每次迭代执行以下步骤

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，因为它在拟合模型之前调用，从而提高计算性能。

参考文献

• https://en.wikipedia.org/wiki/RANSAC

• “Random Sample Consensus: A Paradigm for Model Fitting with Applications to Image Analysis and Automated Cartography” Martin A. Fischler 和 Robert C. Bolles - SRI International (1981)

• “Performance Evaluation of RANSAC Family” Sunglok Choi，Taemin Kim 和 Wonpil Yu - BMVC (2009)

1.1.16.3. Theil-Sen 估计器：基于广义中值的估计器

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

示例

• Theil-Sen 回归

• 稳健线性估计器拟合

理论考虑

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

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

就时间和空间复杂度而言，Theil-Sen 的扩展比例为：

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

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

参考文献

[14]

Dang Xin，Peng Hanxiang，Wang Xueqin 和 Zhang Heping：多元线性回归模型中的 Theil-Sen 估计量。

[15]

20. Kärkkäinen 和 S. Äyrämö：关于稳健数据挖掘中空间中位数的计算。

另请参见 维基百科页面

1.1.16.4. Huber 回归

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

示例

• HuberRegressor 与 Ridge 在具有强烈异常值的 dataset 上的比较

数学细节

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: 稳健统计，伴随尺度估计，第 172 页

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

• HuberRegressor 是尺度不变的。一旦设置了 epsilon，将 X 和 y 按不同值缩小或放大，将产生与之前相同的对异常值的稳健性。相比之下，在 SGDRegressor 中，当 X 和 y 缩放时，必须再次设置 epsilon。

• HuberRegressor 应用于样本数量较少的数据应该更高效，而 SGDRegressor 需要多次遍历训练数据才能产生相同的稳健性。

请注意，此估计器与 R 中的稳健回归实现不同 (https://stats.oarc.ucla.edu/r/dae/robust-regression/)，因为 R 实现使用加权最小二乘实现，根据残差超过某个阈值的程度为每个样本赋予权重。

1.1.17. 分位数回归

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

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

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

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

示例

• 分位数回归

数学细节

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

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

这包括了Pinball损失函数（也称为线性损失函数），参见 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。

由于Pinball损失函数仅在线性残差中是线性的，因此分位数回归比基于平方误差的均值估计对异常值更稳健。介于两者之间的是 HuberRegressor。

参考文献

• Koenker, R., & Bassett Jr, G. (1978). Regression quantiles. Econometrica: journal of the Econometric Society, 33-50.

• Portnoy, S., & Koenker, R. (1997). The Gaussian hare and the Laplacian tortoise: computability of squared-error versus absolute-error estimators. Statistical Science, 12, 279-300.

• Koenker, R. (2005). Quantile Regression. Cambridge University Press.

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\) 中是线性的），并且可以使用相同的技术求解。通过考虑使用这些基函数构建的更高维空间内的线性拟合，模型具有灵活性，可以拟合更广泛的数据范围。

以下是用不同次数的多项式特征将此思想应用于一维数据的示例

此图是使用 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