什么是CVXPY?¶
CVXPY 是一个用于凸优化问题的 Python 嵌入式建模语言。它会自动将问题转换为标准形式,调用求解器,并解包结果。
下面的代码在CVXPY中解决了一个简单的优化问题:
import cvxpy as cp
# Create two scalar optimization variables.
x = cp.Variable()
y = cp.Variable()
# Create two constraints.
constraints = [x + y == 1,
x - y >= 1]
# Form objective.
obj = cp.Minimize((x - y)**2)
# Form and solve problem.
prob = cp.Problem(obj, constraints)
prob.solve() # Returns the optimal value.
print("status:", prob.status)
print("optimal value", prob.value)
print("optimal var", x.value, y.value)
status: optimal
optimal value 0.999999999761
optimal var 1.00000000001 -1.19961841702e-11
状态,由求解方法赋予“最优”值,告诉我们问题已成功解决。最优值(这里基本上是1)是在满足约束条件的所有变量选择中目标函数的最小值。最后打印的内容给出了实现最优目标的x和y的值(基本上是1和0)。
prob.solve() 返回最优值并更新 prob.status、
prob.value,以及问题中所有变量的 value 字段。
改变问题¶
Problems 是不可变的,这意味着它们在创建后不能被更改。要更改目标或约束,请创建一个新的问题。
# Replace the objective.
prob2 = cp.Problem(cp.Maximize(x + y), prob.constraints)
print("optimal value", prob2.solve())
# Replace the constraint (x + y == 1).
constraints = [x + y <= 3] + prob2.constraints[1:]
prob3 = cp.Problem(prob2.objective, constraints)
print("optimal value", prob3.solve())
optimal value 1.0
optimal value 3.00000000006
不可行和无界问题¶
如果问题不可行或无界,状态字段将分别设置为“infeasible”或“unbounded”。问题变量的值字段不会被更新。
import cvxpy as cp
x = cp.Variable()
# An infeasible problem.
prob = cp.Problem(cp.Minimize(x), [x >= 1, x <= 0])
prob.solve()
print("status:", prob.status)
print("optimal value", prob.value)
# An unbounded problem.
prob = cp.Problem(cp.Minimize(x))
prob.solve()
print("status:", prob.status)
print("optimal value", prob.value)
status: infeasible
optimal value inf
status: unbounded
optimal value -inf
请注意,对于最小化问题,如果不可行,则最优值为inf,如果无界,则为-inf。对于最大化问题,情况则相反。
其他问题状态¶
如果由CVXPY调用的求解器解决了问题,但达到的精度低于预期,问题状态将指示达到的较低精度。指示较低精度的状态有
“optimal_inaccurate”
“无界不准确”
“不可行_不准确”
问题变量根据找到的解的类型(即最优、无界或不可行)进行常规更新。
如果求解器完全无法解决问题,CVXPY 会抛出一个 SolverError 异常。
如果发生这种情况,您应该尝试使用其他求解器。详情请参阅
Solvers 的讨论。
CVXPY 提供了以下常量作为不同状态字符串的别名:
OPTIMALINFEASIBLEUNBOUNDEDOPTIMAL_INACCURATEINFEASIBLE_INACCURATEUNBOUNDED_INACCURATEINFEASIBLE_OR_UNBOUNDED
要测试问题是否成功解决,您可以使用
prob.status == OPTIMAL
状态 INFEASIBLE_OR_UNBOUNDED 很少见。它用于当求解器能够确定问题要么是不可行的,要么是无界的,但无法确定具体是哪种情况时。您可以通过重新解决问题来确定确切的状态,其中您将目标函数设置为常数(例如,objective = cp.Minimize(0))。如果新问题以状态码 INFEASIBLE_OR_UNBOUNDED 解决,则原始问题是不可行的。如果新问题以状态 OPTIMAL 解决,则原始问题是无界的。
向量和矩阵¶
Variables 可以是标量、向量或矩阵,这意味着它们是0维、1维或2维的。
# A scalar variable.
a = cp.Variable()
# Vector variable with shape (5,).
x = cp.Variable(5)
# Column vector variable with shape (5, 1).
x = cp.Variable((5, 1))
# Matrix variable with shape (4, 7).
A = cp.Variable((4, 7))
您可以使用您选择的数值库来构造矩阵和向量常量。例如,如果x是表达式A @ x + b中的CVXPY变量,A和b可以是Numpy的ndarrays、SciPy的稀疏矩阵等。A和b甚至可以是不同的类型。
目前以下类型可以用作常量:
NumPy 多维数组
SciPy 稀疏矩阵
这是一个使用向量和矩阵的CVXPY问题的示例:
# Solves a bounded least-squares problem.
import cvxpy as cp
import numpy as np
# Problem data.
m = 10
n = 5
numpy.random.seed(1)
A = np.random.randn(m, n)
b = np.random.randn(m)
# Construct the problem.
x = cp.Variable(n)
objective = cp.Minimize(cp.sum_squares(A @ x - b))
constraints = [0 <= x, x <= 1]
prob = cp.Problem(objective, constraints)
print("Optimal objective value", prob.solve())
print("Optimal variable value")
print(x.value) # A numpy ndarray.
Optimal ojbective value 4.14133859146
Optimal variable value
[ -5.11480673e-21 6.30625742e-21 1.34643668e-01 1.24976681e-01
-4.79039542e-21]
约束条件¶
如示例代码所示,您可以使用 ==、<= 和 >= 在 CVXPY 中构建约束。等式和不等式约束是逐元素的,无论它们涉及标量、向量还是矩阵。例如,约束 0 <= x 和 x <= 1 一起意味着 x 的每个条目都在 0 和 1 之间。
如果你想要表示半定锥约束的矩阵不等式,请参阅半定矩阵。该部分解释了如何表达半定锥不等式。
你不能使用<和>来构造不等式。严格不等式在现实世界中没有意义。此外,你不能将约束条件链接在一起,例如0 <= x <= 1或x == y == 2。Python解释器以CVXPY无法捕获的方式处理链接的约束条件。如果你写了一个链接的约束条件,CVXPY将会抛出异常。
参数¶
Parameters 是常量的符号表示。参数的目的是在不重建整个问题的情况下改变问题中常量的值。在许多情况下,多次求解参数化程序可能比重复求解新问题要快得多:在阅读完本节后,请务必阅读关于Disciplined Parametrized Programming (DPP)的教程。
当你创建一个参数时,你可以选择指定一些属性,比如参数项的符号、参数是否对称等。这些属性在Disciplined Convex Programming中使用,除非指定,否则是未知的。参数在创建后可以随时被赋予一个常数值。这个常数值必须具有与创建参数时指定的相同的维度和属性。
# Positive scalar parameter.
m = cp.Parameter(nonneg=True)
# Column vector parameter with unknown sign (by default).
c = cp.Parameter(5)
# Matrix parameter with negative entries.
G = cp.Parameter((4, 7), nonpos=True)
# Assigns a constant value to G.
G.value = -np.ones((4, 7))
你可以用一个值初始化一个参数。以下代码段是等价的:
# Create parameter, then assign value.
rho = cp.Parameter(nonneg=True)
rho.value = 2
# Initialize parameter with a value.
rho = cp.Parameter(nonneg=True, value=2)
计算权衡曲线是参数的一个常见用途。下面的示例计算了一个LASSO问题的权衡曲线。
import cvxpy as cp
import numpy as np
import matplotlib.pyplot as plt
# Problem data.
n = 15
m = 10
np.random.seed(1)
A = np.random.randn(n, m)
b = np.random.randn(n)
# gamma must be nonnegative due to DCP rules.
gamma = cp.Parameter(nonneg=True)
# Construct the problem.
x = cp.Variable(m)
error = cp.sum_squares(A @ x - b)
obj = cp.Minimize(error + gamma*cp.norm(x, 1))
prob = cp.Problem(obj)
# Construct a trade-off curve of ||Ax-b||^2 vs. ||x||_1
sq_penalty = []
l1_penalty = []
x_values = []
gamma_vals = np.logspace(-4, 6)
for val in gamma_vals:
gamma.value = val
prob.solve()
# Use expr.value to get the numerical value of
# an expression in the problem.
sq_penalty.append(error.value)
l1_penalty.append(cp.norm(x, 1).value)
x_values.append(x.value)
plt.rc('text', usetex=True)
plt.rc('font', family='serif')
plt.figure(figsize=(6,10))
# Plot trade-off curve.
plt.subplot(211)
plt.plot(l1_penalty, sq_penalty)
plt.xlabel(r'\|x\|_1', fontsize=16)
plt.ylabel(r'\|Ax-b\|^2', fontsize=16)
plt.title('Trade-Off Curve for LASSO', fontsize=16)
# Plot entries of x vs. gamma.
plt.subplot(212)
for i in range(m):
plt.plot(gamma_vals, [xi[i] for xi in x_values])
plt.xlabel(r'\gamma', fontsize=16)
plt.ylabel(r'x_{i}', fontsize=16)
plt.xscale('log')
plt.title(r'\text{Entries of x vs. }\gamma', fontsize=16)
plt.tight_layout()
plt.show()
权衡曲线可以很容易地并行计算。下面的代码并行计算了上述LASSO问题中每个\(\gamma\)的最优x。
from multiprocessing import Pool
# Assign a value to gamma and find the optimal x.
def get_x(gamma_value):
gamma.value = gamma_value
result = prob.solve()
return x.value
# Parallel computation (set to 1 process here).
pool = Pool(processes = 1)
x_values = pool.map(get_x, gamma_vals)