初学者快速示例#

在这里,我们提供了一个适合初学者的示例,介绍了skscope在线性回归中用于特征/变量选择的基本用法。

介绍#

让我们考虑一个数据集,其中包含\(n=100\)个独立观察值的有趣响应变量和\(p=10\)个预测变量:

from sklearn.datasets import make_regression

## generate data
n, p, s = 100, 10, 3
X, y, true_coefs = make_regression(n_samples=n, n_features=p, n_informative=s, coef=True, random_state=0)

print("X shape:", X.shape)
>>> X shape: (100, 10)

print("Y shape:", y.shape)
>>> Y shape: (100,)

在这个场景中,我们的目标是使用可用的预测变量来解释响应变量的变化。然而,我们不知道由true_coefs表示的真实参数。我们只知道在\(p=10\)个预测变量中,只有\(s=3\)个变量实际上影响响应变量。此外,我们假设预测变量和响应变量之间存在线性关系:

\[y=X \theta^{*} +\epsilon.\]

其中

  • \(y \in R^n\)\(X \in R^{n\times p}\) 分别表示响应变量和预测变量的观测值;

  • \(\epsilon = (\epsilon_1, \ldots, \epsilon_n)\) 是一个由 \(n\) 个独立同分布的零均值随机噪声组成的向量;

  • \(\theta^{*}\) 是一个需要估计的未知系数向量。

使用skscope,用户可以有效地估计\(\theta^{*}\)并通过变量选择识别有影响力的预测变量。下面展示了实现这一目标的逐步过程。

逐步过程#

稀疏约束优化视角#

让我们首先将变量选择问题表述为一个优化问题:

\[\arg\min_{\theta \in R^p} ||y-X \theta||^{2} s.t. ||\theta||_0 \leq s,\]

其中:

  • \(\theta\) 是一个需要优化的系数向量。

  • \(||y-X \theta||^{2}\) 是量化系数向量 \(\theta\) 质量的目标函数。

  • \(||\theta||_0\) 计算 \(\theta\) 中非零元素的数量。稀疏性约束 \(||\theta||_0 \leq s\) 反映了我们的先验知识,即有影响力的预测变量的数量小于或等于 \(s\)

这个优化问题的直观解释是找到一个小的预测变量子集,从而产生具有最理想预测准确性的线性模型。通过这样做,我们有效地利用了观测数据中的信息和我们的先验知识。

通过skscope#解决SCO

上述优化问题是一个稀疏约束优化(SCO)问题。由于skscope是为通用SCO设计的,它可以轻松解决这个问题。使用skscope的美妙之处在于其简单性——只要目标函数\(||y-X \theta||^{2}\)被正确实现,它就能解决SCO问题。对于当前的例子,我们将目标函数编程为一个名为objective_function(coefs)的Python函数:

import jax.numpy as jnp

def objective_function(coefs):
    return jnp.linalg.norm(y - X @ coefs)

请注意,我们在第一行导入了所需的jax库[*]。该库提供了一个与numpy兼容的多维数组API,并支持自动微分。由于jax.numpy的使用与numpy非常相似,在这里,我们可以将jax.numpy理解为等同于numpy[]。通过上面的代码片段,objective_function的实现与数学函数\(||y-X \theta||^{2}\)完全一致。

接下来,我们将objective_function输入到skscope中的一个求解器中,以获得SCO问题的实际解决方案。下面,我们导入ScopeSolver []并正确配置它:

from skscope import ScopeSolver

scope_solver = ScopeSolver(
    dimensionality=p,  ## there are p parameters
    sparsity=s,        ## we want to select s variables
)

在上述配置中,我们设置了:

  • dimensionality: 参数的数量,必须指定;

  • sparsity: 期望的稀疏度水平。

在上面的例子中,\(\theta\) 是参数,所以 dimensionality\(p\) 并且 sparsity\(s\)

定义了scope_solverobjective_function后,我们可以使用一个命令来解决SCO问题:

scope_solver.solve(objective_function)

solve 方法是 skscope 中求解器的主要方法。它将目标函数作为优化目标,并指示算法执行优化过程。

检索解决方案#

由于skscope中的求解器带有提取结果的必要功能,我们只需一行代码即可获得所需的结果。例如,如果我们对检索有效变量感兴趣,我们可以使用get_support方法:

import numpy as np

estimated_support_set = scope_solver.get_support()
print("Estimated effective predictors:", estimated_support_set)
>>> Estimated effective predictors: [3 5 6]
print("True effective predictors:", np.nonzero(true_coefs)[0])
>>> True effective predictors: [3 5 6]

我们可以观察到,估计的有效预测变量与真实变量相匹配,这表明了skscope中求解器的准确性。

此外,我们可能对回归系数感兴趣:

  • get_estimated_params 方法返回优化后的系数。

est_coefs = scope_solver.get_estimated_params()
print("Estimated coefficient:", np.around(est_coefs, 2))
>>> Estimated coefficient: [ 0.    0.    0.   82.19  0.   88.31 70.05  0.    0.    0.  ]
print("True coefficient:", np.around(true_coefs, 2))
>>> True coefficient: [ 0.    0.    0.   82.19  0.   88.31 70.05  0.    0.    0.  ]

在输出中,我们可以观察到估计的系数与真实系数非常接近。

进一步阅读#

脚注#