solve.QP 的数值问题，R 中的二次规划函数

Question

假设我有以下二次规划问题，其中的数字非常大：

> require(quadprog)
> Q <- diag(c(34890781627, 34890781627, 34890781627, 34890781627, 34890781627))
> c <- c(133013236723, 29459621018, 31362634710, 24032348447, 23332381578)
> ( A <- t(kronecker(diag(1,5),as.matrix(c(1,-1)))) )
     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
[1,]    1   -1    0    0    0    0    0    0    0     0
[2,]    0    0    1   -1    0    0    0    0    0     0
[3,]    0    0    0    0    1   -1    0    0    0     0
[4,]    0    0    0    0    0    0    1   -1    0     0
[5,]    0    0    0    0    0    0    0    0    1    -1
> ( e <- rep(c(0,1),5) )
 [1] 0 1 0 1 0 1 0 1 0 1

约束是Ax >= -e：我将所有变量都约束在[0,1]内。我猜是因为这里涉及的数字太大了，所以会造成一些数值问题。如果我运行

> solve.QP(Q,c,A,bvec=-e)$solution

我会得到Error in solve.QP(Q, c, A, bvec = -e) : constraints are inconsistent, no solution!

当我尝试将所有内容除以一个大数以缩小这些疯狂的数字时，该函数能够生成一些输出，但不正确：

> 1E4*solve.QP(1E-8*Q,1E-4*c,A,bvec=-1E-4*e)$solution
[1] 1 1 1 1 1

我认为正确答案是

> w <- (1/diag(Q))*c
> w[w>1] <- 1
> w[w<0] <- 0
> w
[1] 1.0000000 0.8443382 0.8988803 0.6887879 0.6687263

我真的不需要在这里 QP 因为 Q 矩阵是对角线的，但我只是想用这个作为说明并获得一些关于处理与 solve.QP() 相关的数值问题的建议.假设我有一个密集的 Q 矩阵，我应该怎么做才能得到正确的结果？

提前致谢！

Answer 1

抱歉各位，我刚刚意识到我在尝试缩小矩阵 Q 和向量 c 时犯了一个错误。

min x^T Q x /2 - c^T x 等同于 min a * (x^T Q x /2 - c^T x) 任何正数 a。所以我需要做的就是将 Q 和 c 乘以相同的小数，如下所示：

> solve.QP(1E-8*Q,1E-8*c,A,bvec=-e)$solution
[1] 1.0000000 0.8443382 0.8988803 0.6887879 0.6687263

这次它给出了正确的结果。之前我想改造x = ay，但是我的逻辑不对。正确的做法是：

min x^T Q x /2 - c^T x 等同于 min y^T (a^2*Q) y /2 - a*c^T y 任何正数 a。并且约束也会相应地变为Ay >= -e/a。所以你也可以做

1E-4*solve.QP(1E-8*Q,1E-4*c,A,bvec=-1E4*e)$solution
[1] 1.0000000 0.8443382 0.8988803 0.6887879 0.6687263

剩余问题

一般来说，我们什么时候应该考虑这个数值问题，最佳收缩尺寸是多少a？或者有什么方法可以更好的解决这个问题？

solve.QP 的数值问题，R 中的二次规划函数

A numerical issue with solve.QP, the quadratic programming function in R

r

quadratic-programming

剩余问题