将 scipy.optimize.minimize 限制为整数值

Restrict scipy.optimize.minimize to integer values

我正在使用 scipy.optimize.minimize 来优化一个答案只能是整数的现实问题。我当前的代码如下所示:

from scipy.optimize import minimize

def f(x):
    return (481.79/(5+x[0]))+(412.04/(4+x[1]))+(365.54/(3+x[2]))+(375.88/(3+x[3]))+(379.75/(3+x[4]))+(632.92/(5+x[5]))+(127.89/(1+x[6]))+(835.71/(6+x[7]))+(200.21/(1+x[8]))

def con(x):
    return sum(x)-7

cons = {'type':'eq', 'fun': con}

print scipy.optimize.minimize(f, [1,1,1,1,1,1,1,0,0], constraints=cons, bounds=([0,7],[0,7],[0,7],[0,7],[0,7],[0,7],[0,7],[0,7],[0,7]))

这产生:

x: array([  2.91950510e-16,   2.44504019e-01,   9.97850733e-01,
     1.05398840e+00,   1.07481251e+00,   2.60570253e-01,
     1.36470363e+00,   4.48527831e-02,   1.95871767e+00]

但我希望它使用整数值进行优化(将所有 x 四舍五入到最接近的整数并不总是给出最小值)。

有没有办法只对整数值使用 scipy.optimize.minimize

(我想我可以创建一个数组,其中包含 x 的所有可能排列,并为每个组合计算 f(x),但这似乎不是一个非常优雅或快速的解决方案。)

浆液

经过一番研究,我认为您的 objective 函数不是线性的。我在 Python pulp library but pulp doesn't like that we're dividing by a float and 'LpAffineExpression'. suggests that linear programming "doesn't understand divisions" but that comment is in context of adding constraints, not the objective function. That comment pointed me to "Mixed Integer Linear Fractional Programming (MILFP)" and on Wikipedia.

中重新创建了问题

如果它真的有效,你可以在 pulp 中使用以下方法(也许有人能弄清楚原因):

import pulp

data = [(481.79, 5), (412.04, 4), (365.54, 3)] #, (375.88, 3), (379.75, 3), (632.92, 5), (127.89, 1), (835.71, 6), (200.21, 1)]
x = pulp.LpVariable.dicts('x', range(len(data)), lowBound=0, upBound=7, cat=pulp.LpInteger)

numerator = dict((i,tup[0]) for i,tup in enumerate(data))
denom_int = dict((i,tup[1]) for i,tup in enumerate(data))

problem = pulp.LpProblem('Mixed Integer Linear Programming', sense=pulp.LpMinimize)

# objective function (doesn't work)
# TypeError: unsupported operand type(s) for /: 'float' and 'LpAffineExpression'
problem += sum([numerator[i] / (denom_int[i] + x[i]) for i in range(len(data))])

problem.solve()

for v in problem.variables():
  print(v.name, "=", v.varValue)

暴力破解 scipy.optimize

您可以在函数中为每个 x 使用 bruteslice 的范围。如果您的函数中有 3 个 x,那么您的范围元组中也会有 3 个 slice。所有这一切的关键是将 1step 大小添加到 slice(start, stop,step ) 所以 slice(#, #, 1).

from scipy.optimize import brute
import itertools

def f(x):
  return (481.79/(5+x[0]))+(412.04/(4+x[1]))+(365.54/(3+x[2]))

ranges = (slice(0, 9, 1),) * 3
result = brute(f, ranges, disp=True, finish=None)
print(result)

itertools解决方案

或者您可以使用 itertools 生成所有组合:

combinations = list(itertools.product(*[[0,1,2,3,4,5,6,7,8]]*3))

values = []
for combination in combinations:
  values.append((combination, f(combination)))

best = [c for c,v in values if v == min([v for c,v in values])]
print(best)

注意:这是您原始函数的缩小版本,用于示例目的。

可能对您的问题有帮助的一件事是:

max([x-int(x)])=0

这不会完全解决您的问题,该算法仍会尝试作弊,并且您将获得具有一定程度错误的值 ~±5e-10,它仍会尝试并针对以下错误进行优化scipy的算法,但总比没有好。

cons = ({'type':'eq', 'fun': con},
        {'type':'eq','fun': lambda x : max([x[i]-int(x[i]) for i in range(len(x))])})

在我知道解决方案的一些优化上测试了这个过程,这个过程比无约束搜索对初始值更敏感,它得到相当准确的答案但是解决方案实际上可能找不到真正的值,你是基本上需要优化过程的大跳跃(它用来确保它没有优化到局部最小值)来搜索样本 space,因为较小的增量通常不足以移动到下一个数字。

暴力破解的通用函数。在 scipy 中比 brute 做得更好,因为 scipy 实际上使用浮点数运行函数,而不仅仅是整数,尽管 range 明确这么说,正如 Jarad 所说

def brute(func, arg_ranges, finish=min):
if isinstance(arg_ranges, dict):
    args = {k:np.unique(np.hstack([a for r in rs for a in r]) if isinstance(rs, list) else [a for a in rs]) for k,rs in arg_ranges.items()}
    print(args)
    return finish([(dict(zip(args.keys(), vs)), func(**dict(zip(args.keys(), vs)))) for vs in itertools.product(*args.values())], key=lambda x: x[1])
elif isinstance(arg_ranges, list):
    return finish([(i, func(i)) for r in arg_ranges for i in r], key=lambda x: x[1])
else:
    return finish([(i, func(i)) for i in arg_ranges], key=lambda x: x[1])

print(brute(lambda x,y: x / (y + 2), {'x':[range(1,5,2), range(0,6,1)], 'y':range(2,5,1)}))

这里有一个解决混合整数非线性规划问题的方法Python Gekko(我维护的一个包):

from gekko import GEKKO

m = GEKKO(remote=False)
x = m.Array(m.Var,9,lb=0,ub=7,integer=True)

def f(x):
    return (481.79/(5+x[0]))+(412.04/(4+x[1]))\
           +(365.54/(3+x[2]))+(375.88/(3+x[3]))\
           +(379.75/(3+x[4]))+(632.92/(5+x[5]))\
           +(127.89/(1+x[6]))+(835.71/(6+x[7]))\
           +(200.21/(1+x[8]))

m.Minimize(f(x))
m.Equation(sum(x)==7)
m.options.SOLVER=1
m.solve()
print(x)

这给出了解决方案:

 ---------------------------------------------------
 Solver         :  APOPT (v1.0)
 Solution time  :  0.0529 sec
 Objective      :  859.5269999999999
 Successful solution
 ---------------------------------------------------


[[0.0] [0.0] [1.0] [1.0] [1.0] [0.0] [1.0] [0.0] [3.0]]