线性规划与迭代重新加权最小二乘运行时间

Question

我试图运行最小绝对偏差回归（L1 回归）。我通过设计一个线性程序并使用 CVXOPT（Python）求解来做到这一点。特征数量（包括常量）为 13 个，样本数量约为 100K。花了好几个小时才完成。 然后我意识到我可以用

statsmodels' QuantReg

运行

q=0.5

。花了几秒钟。 我读到他们使用迭代重新加权最小二乘法（IRLS）。 为此目的，IRLS 是否比线性编程更快，或者线性优化问题的不同表述，以及使用像 Gurobi 这样更好的求解器，会大大减少线性程序的运行时间？ 我问的原因是因为我想尝试一个损失函数，该函数涉及两个具有不同权重的不同绝对函数的求和。我相信我也可以使用 IRLS 来解决这个问题，但我不确定，而且我对线性规划更满意。

我用于线性规划的公式如下：

CVXOPT 的结果代码如下：

def l1_fit(x,y):
'''

:param x: feature matrix
:param y: label vector
:return: LAD coefficients

The problem I construct is:
    min c'*(v|A)
    subject to:
        (-I|x)*(v|A)<=y
        (-I|-x)*(v|A)<=-y
'''
c = matrix(np.concatenate((np.ones(y.shape[0]),np.zeros(x.shape[1]))).astype(float).tolist(),
           (y.shape[0]+x.shape[1],1),'d')
G = scipy_sparse_to_spmatrix(np.transpose(hstack((vstack((-identity(y.shape[0]),np.transpose(x))),
                          vstack((-identity(y.shape[0]),np.transpose(-x))))).astype(float)))
h = matrix(np.concatenate((y,-y)).astype(float).tolist(),(2*y.shape[0],1),'d')
sol = lp(c, G, h)
return sol[1][y.shape[0]:]

Answer 1

如果您可以访问 Gurobi 或 CPlex，那么您应该遵循 Erwin Kalvelagen - 线性规划和 LAD 回归中的分析。

它表明线性规划问题的最佳表述如下：

但是，我不确定将模型转换为 LP 问题是否可行。 我生成了一个模型矩阵 $ oldsymbol{A} \in \mathbb{R}^{100,000 imes 15}$.

我使用

Convex.jl

（相当于 Julia 的

CVXPY

）和 ECOS Solver（应该在

CVXPY

中可用）以其直接形式解决了这个问题。
问题已在

50 [Sec]

下解决。所以快了大约 1-2 个数量级。

此外，我还实现了一个基于“迭代重新加权最小二乘法”（IRLS）的求解器。我可以在不到一秒内解决问题。

速度的原因是 IRLS 最消耗的操作（求解线性系统）是在

A' * A

上完成的，在您的情况下，它非常小。

我在 Julia 中实现了如下：

最小绝对偏差 (LAD) 线拟合/回归
• 。
通过稳健拟合来拟合数据样本
• 。
具有 P 范数正则化的最小二乘线性回归，其中 1≤P≤2
• 。

该代码是在 Julia 中使用我编写的相对高效的 IRLS 代码实现的（无分配）。

---

您可以在 Python 中使用

NumPy
、

SciPy

和

Numba

来非常接近地模仿它。所以解决你的问题几乎需要一个小时的编码时间。

代码可以在我的

StackOverflow GitHub Repository

上找到（查看 StackOverflow\Q64422417 文件夹）。