python - 为什么我会收到此 DCPError？

Question

我正在尝试使用 CVXPY 优化二元投资组合向量，使其大于基准。

import cvxpy as cp
import numpy as np

# Generate a random non-trivial quadratic program.

n = 10 # number of options

np.random.seed(1)
mu = np.random.randn(n) # expected means
var_covar = np.random.randn(n,n) # variance-covariance matrix
var_covar = var_covar.T.dot(var_covar) # cont'd
bench_cov = np.random.randn(n) # n-length vector of cov(benchmark, returns)

lamd = 0.01 # risk tolerance

# Define and solve the CVXPY problem.

x = cp.Variable(n, boolean=True)

prob = cp.Problem(cp.Maximize(mu.T@x + lamd * (cp.quad_form(x, var_covar) - (2 * bench_cov.T@x))), [cp.sum(x) == 4])

prob.solve()

我使用 CVXPY 版本 1.1.0a0(直接从 github 下载)遇到此错误:

DCPError: Problem does not follow DCP rules. Specifically:

The objective is not DCP, even though each sub-expression is.

You are trying to maximize a function that is convex.

据我所知，最大化凸函数是非常困难的，但我从一篇论文中得到了这个方程。我想我一定做错了什么，因为我是二次规划和 CVXPY 的新手。

谢谢!

Answer 1

您的模型的问题是 max x'Qx 是非凸的。由于我们有二进制变量x，我们可以使用一个技巧。

定义

y(i,j) = x(i)*x(j)

作为额外的二进制变量。然后我们就可以这样写

sum((i,j), x(i)*Q(i,j)*x(j)) 

作为

sum((i,j), y(i,j)*Q(i,j)) 

二进制乘法y(i,j) = x(i)*x(j)可以线性化为:

 y(i,j) <= x(i)
 y(i,j) <= x(j)
 y(i,j) >= x(i)+x(j)-1

通过这种重新表述，我们有了一个完全线性的模型。它是一个 MIP，因为我们有二进制变量。

我们可以在 CVXPY 中执行此操作:

import numpy as np
import cvxpy as cp

# Generate a random non-trivial quadratic program.

n = 10 # number of options

np.random.seed(1)
mu = np.random.randn(n) # expected means
var_covar = np.random.randn(n,n) # variance-covariance matrix
var_covar = var_covar.T.dot(var_covar) # cont'd
bench_cov = np.random.randn(n) # n-length vector of cov(benchmark, returns)

lamd = 0.01 # risk tolerance

e = np.ones((1,n))

x = cp.Variable((n,1), "x", boolean=True)
y = cp.Variable((n,n), "y", boolean=True)


prob = cp.Problem(cp.Maximize(mu.T@x + lamd * (cp.sum(cp.multiply(y,var_covar)) -2*bench_cov.T@x) ),
                  [y <= x@e, y <= (x@e).T, y >= x@e + (x@e).T - e.T@e, cp.sum(x)==4 ])

prob.solve(solver=cp.ECOS_BB)
print("status",prob.status)
print("obj",prob.value)
print("x",x.value)

这给出了结果:

status optimal
obj 4.765120794509871
x [[1.00000000e+00]
 [3.52931931e-10]
 [3.80644178e-10]
 [2.53300872e-10]
 [9.99999999e-01]
 [1.79871537e-10]
 [1.00000000e+00]
 [3.46298454e-10]
 [9.99999999e-01]
 [1.00172269e-09]]

注释:

• 我们鼓励您使用比 ECOS_BB 更好的 MIP 求解器。对于这个模型，它给出了正确的结果，但它有点像一个玩具求解器，并且已知会在更困难的数据集上给出问题。
• 我不了解该模型的经济学原理。我们在这里将风险最大化。根据此模型的结果做出投资决策可能并不谨慎。
• 请注意，一些高端求解器(例如 Cplex 和 Gurobi)会自动执行此重新公式化。但是 CVXPY 不允许您将非凸模型传递给求解器。