Python 中的有界优化中的约束条件求和

Question

我正在尝试在某些限制下最小化最小二乘和。在大多数情况下，除了一个约束之外，这个问题看起来相当简单。我很难找到一种对约束进行建模的方法，即大于给定数字的所有值的总和不能超过常量值。

我目前有一个工作解决方案，它迭代可能的组合并使用 scipy 的最小化函数选择最佳解决方案。然而，就时间（或计算能力）而言，该解决方案并不是最有效的。是否有另一种方法可以一次性解决这个问题？

我一直在研究 Gekko，因为它使一般代码更具可读性，但是一旦我添加一个使用求和的方程，效率就比循环差。

壁虎代码：

w_bar = <List of decimal values>

model = GEKKO(remote=False)

w = model.Array(model.Var, len(w_bar), lb=0, ub=0.3)

# Sum of all values = 1
model.Equation(model.sum(w) == 1)

# Sum of all values greater than const1 cannot exceed const2
# This is the part that degrades the speed greatly
model.Equation(model.sum([model.if3(w[i] - const1, 0, w[i]) for i in range(len(w))]) <= const2)

model.Minimize(model.sum([(w_bar[i] - w[i]) ** 2 for i in range(len(w_bar))]))

model.solve()

使用 scipy，我创建了一个循环，在其中明确指示该值是否受到约束。它有效，但速度很慢。对于大量输入（300+），大约需要 10 秒。

在 Gekko 中，我添加了一个方程来约束，但这极大地降低了性能。对于少量输入，它运行得相对较快（不到 1 秒），但一旦输入超过 300 个，求解时间就会增加到大约 90 秒。

我是否正确思考了方程式？我需要以某种方式重塑这个问题吗？

Answer 1

这需要大量猜测，因为您没有提供问题的规模或实际的约束常数。

我建议您放弃最小二乘的想法，这样您就可以将其构建为混合整数线性程序。这有点棘手，但确实有效。设置变量：

• |err|

  ，

  w - w_bar

  的绝对值，要最小化

• w

  ，限制在 (0, 0.3)

• whi

  ，如果高于 const1，则等于

  w

  ，否则为 0

• ishi

  ，一个二进制表示 w > const1

import numpy as np
from scipy.optimize import milp, Bounds, LinearConstraint
import scipy.sparse as sp

rand = np.random.default_rng(seed=0)
n = 7
w_bar = rand.uniform(low=-0.1, high=0.4, size=n)

# Variables: n |err|, n w, n whi, n ishi

# err, w, whi continuous, ishi integral
integrality = np.zeros(shape=4*n, dtype=np.uint8)
integrality[-n:] = 1

# Constrain error to absolute difference of w and wbar
# err >= w - wbar: -err + w <= wbar
# err >= wbar - w:  err + w >= wbar
err_constraint = LinearConstraint(
    A=sp.bmat((
        (-sp.eye(n), sp.eye(n), sp.csc_array((n, 2*n))),
        ( sp.eye(n), sp.eye(n), sp.csc_array((n, 2*n))),
    ), format='csc'),
    lb=np.concatenate((np.full(n, -np.inf), w_bar)),
    ub=np.concatenate((w_bar, np.full(n, np.inf))),
)

# sum w = 1
one_constraint = LinearConstraint(
    A=sp.hstack((np.zeros((1, n)), np.ones(n), np.zeros(2*n)), format='csc'),
    lb=1, ub=1,
)

const1 = 0.2
const2 = 0.58

# Activate ishi based on value of w and const1
# ishi >= w - const1: w - ishi <= const1
# ishi <= w/const1:   w/const1 - ishi >= 0
ishi_constraint = LinearConstraint(
    A=sp.bmat((
        (sp.csc_array((n, n)), sp.eye(n),        sp.csc_array((n, n)), -sp.eye(n)),
        (sp.csc_array((n, n)), sp.eye(n)/const1, sp.csc_array((n, n)), -sp.eye(n)),
    ), format='csc'),
    lb=np.concatenate((
        np.full(n, -np.inf),
        np.zeros(n),
    )),
    ub=np.concatenate((
        np.full(n, const1),
        np.full(n, np.inf),
    )),
)

# Constrain whi based on whether ishi is set
# whi <= w                  w - whi          >= 0
# whi <= ishi*M                -whi + M*ishi >= 0
# whi >= w - M*(1 - ishi)   w - whi + M*ishi <= M
w_max = 0.3
M = 2*w_max
whi_constraint = LinearConstraint(
    A=sp.bmat((
        (sp.csc_array((n, n)), sp.eye(n),            -sp.eye(n), sp.csc_array((n, n))),
        (sp.csc_array((n, n)), sp.csc_array((n, n)), -sp.eye(n), M*sp.eye(n)),
        (sp.csc_array((n, n)), sp.eye(n),            -sp.eye(n), M*sp.eye(n)),
    ), format='csc'),
    lb=np.concatenate((
        np.zeros(n),
        np.zeros(n),
        np.full(n, -np.inf),
    )),
    ub=np.concatenate((
        np.full(n, np.inf),
        np.full(n, np.inf),
        np.full(n, M),
    )),
)

# Constrain whi sum based on value of const2
# sum whi <= const2
whisum_constraint = LinearConstraint(
    A=sp.hstack((
        sp.csc_array((1, n)), sp.csc_array((1, n)), np.ones(n), sp.csc_array((1, n)),
    ), format='csc'),
    lb=-np.inf, ub=const2,
)


result = milp(
    c=np.concatenate((
        np.ones(n),   # minimize |err|
        np.zeros(3*n),  # don't optimize w, whi or ishi
    )),
    integrality=integrality,
    bounds=Bounds(
        #                  |err|,     w,    whi, ishi
        lb=np.repeat((-np.inf,     0,      0, 0), n),
        ub=np.repeat((+np.inf, w_max, np.inf, 1), n),
    ),
    constraints=(
        err_constraint, one_constraint, ishi_constraint, whi_constraint, whisum_constraint,
    ),
)
assert result.success, result.message
err, w, whi, ishi = np.split(result.x, (n, 2*n, 3*n))

print(f'Cost: {result.fun:.3f}')
print(f'Sum: {w.sum():.2f} == 1')
print(f'const2: {whi.sum():.2f} <= {const2}')
print('err, w, whi, ishi:')
np.set_printoptions(formatter={'float': '{:5.2f}'.format})
print(result.x.reshape((4, -1)).T)

Cost: 0.291
Sum: 1.00 == 1
const2: 0.58 <= 0.58
err, w, whi, ishi:
[[ 0.02  0.20  0.00  0.00]
 [-0.00  0.03  0.00  0.00]
 [ 0.08  0.00  0.00  0.00]
 [ 0.09  0.00  0.00  0.00]
 [ 0.01  0.30  0.30  1.00]
 [ 0.08  0.28  0.28  1.00]
 [ 0.02  0.19  0.00  0.00]]