使用联合边界/符号约束进行优化
Optimization with joint bounds / symbolic constraint
我正在尝试对向量 [a1,a2,a3] 执行最小二乘优化,具有以下约束,其中 k
是常数:
-k < a3/a2 < k
我将 post 我的代码的简化版本,希望这能让我清楚地知道我想要的是什么。
from scipy import optimize
def loss_function(a_candidate):
return MyObject(a_candidate).mean_square_error()
def feasibility(a_candidate):
# Returns a boolean
k = 1.66711
a2 = a_candidate[1]
a3 = a_candidate[2]
return -k < a3/a2 < k
def feasibility_scipy(a_candidate):
# As I understand it, for SciPy the constraint has to be a number
boolean = feasibility(a_candidate)
if boolean:
return 0
else:
return 1
# Equality constraint forcing feasibility_scipy to be 0.
constraint = optimize.NonlinearConstraint(feasibility_scipy, 0, 0)
optimize_results = optimize.minimize(
loss_function,
a_init, # Initial guess
constraints=constraint)
由于初始猜测 a_init
的生成方式,它位于 feasibility
为 False
的区域。 (我们首先需要使用数值方法的原因是早期的封闭形式方法返回了一个不可行的解决方案。有可能提供一个非常差的可行猜测,如 (0,0,0),但这将是离真正的解决方案更远)。
由于constraint
的梯度几乎处处为零,优化例程无法找到离开这个不可行(不可接受)区域的出路,并且没有成功终止。使用 SLSQP
它仅在 1 次迭代后停止,并显示消息 Singular matrix C in LSQ subproblem
。使用 trust-constr
求解器,它达到了函数计算的最大次数,我相信它没有离开不可行区域,因为 constr_violation
是 1.0
.
据我了解,在 SciPy 中不可能在 a2
和 a3
上提供 'joint bound'(对发明的术语表示歉意),意思是我被迫使用 NonlinearConstraint
方法。
正确的做法是什么? (一些搜索建议我可能想尝试带有符号约束的 mystic
包。但在我花时间学习这个新包之前,我想看看 Whosebug 是否有基于 SciPy解决方案。或者,如果您在 mystic
中知道如何执行此操作,一些示例代码将非常有帮助。)
没关系,我认为解决方案可能只是:
def a_ratio(a_candidate):
a2 = a_candidate[1]
a3 = a_candidate[2]
return a3/a2
feasibility_constraint = optimize.NonlinearConstraint(a_ratio,-k,k)
而不是使用 -k < a3/a2 < k
我会做这两个线性约束:
-k*a2 <= a3
a3 <= k*a2
我知道这个问题已经得到解答,但您想知道如何在 mystic
中处理它。我是 mystic
的作者。像这样:
>>> import mystic as my
>>> import numpy as np
>>>
>>> truth = np.array([1,1,1])
>>>
>>>
>>> class MyObject(object):
... def __init__(self, candidate):
... self.candidate = candidate
... self.truth = truth
... def mean_square_error(self):
... return ((self.truth - self.candidate)**2).sum()
...
>>>
>>> def loss_function(a_candidate):
... return MyObject(a_candidate).mean_square_error()
...
>>>
然后是解决方案:
>>> equations = """
... -k < a3/a2
... a3/a2 < k
... """
>>>
>>> eqn = my.symbolic.simplify(equations, variables=['a1','a2','a3'], locals=dict(k=1.66711), all=False, target='a3')
>>> c = my.symbolic.generate_constraint(my.symbolic.generate_solvers(eqn, variables=['a1','a2','a3']), join=my.constraints.or_)
>>> res = my.solvers.fmin(loss_function, x0=[2,3,4], constraints=c, disp=True, xtol=1e-8)
Optimization terminated successfully.
Current function value: 0.000000
Iterations: 121
Function evaluations: 229
>>> res
array([1., 1., 1.])
mystic
的不同之处在于它构建了一个“运算符”c
,不允许优化器 select 不可行的解决方案。并不是说它用布尔值拒绝不好的解决方案,而是它将“坏”候选者转换为“好”候选者,如下所示:
>>> c([2,10,4])
[2, 10, 4]
>>> c([2,1,4])
[2, 1, 1.6671099999999974]
我正在尝试对向量 [a1,a2,a3] 执行最小二乘优化,具有以下约束,其中 k
是常数:
-k < a3/a2 < k
我将 post 我的代码的简化版本,希望这能让我清楚地知道我想要的是什么。
from scipy import optimize
def loss_function(a_candidate):
return MyObject(a_candidate).mean_square_error()
def feasibility(a_candidate):
# Returns a boolean
k = 1.66711
a2 = a_candidate[1]
a3 = a_candidate[2]
return -k < a3/a2 < k
def feasibility_scipy(a_candidate):
# As I understand it, for SciPy the constraint has to be a number
boolean = feasibility(a_candidate)
if boolean:
return 0
else:
return 1
# Equality constraint forcing feasibility_scipy to be 0.
constraint = optimize.NonlinearConstraint(feasibility_scipy, 0, 0)
optimize_results = optimize.minimize(
loss_function,
a_init, # Initial guess
constraints=constraint)
由于初始猜测 a_init
的生成方式,它位于 feasibility
为 False
的区域。 (我们首先需要使用数值方法的原因是早期的封闭形式方法返回了一个不可行的解决方案。有可能提供一个非常差的可行猜测,如 (0,0,0),但这将是离真正的解决方案更远)。
由于constraint
的梯度几乎处处为零,优化例程无法找到离开这个不可行(不可接受)区域的出路,并且没有成功终止。使用 SLSQP
它仅在 1 次迭代后停止,并显示消息 Singular matrix C in LSQ subproblem
。使用 trust-constr
求解器,它达到了函数计算的最大次数,我相信它没有离开不可行区域,因为 constr_violation
是 1.0
.
据我了解,在 SciPy 中不可能在 a2
和 a3
上提供 'joint bound'(对发明的术语表示歉意),意思是我被迫使用 NonlinearConstraint
方法。
正确的做法是什么? (一些搜索建议我可能想尝试带有符号约束的 mystic
包。但在我花时间学习这个新包之前,我想看看 Whosebug 是否有基于 SciPy解决方案。或者,如果您在 mystic
中知道如何执行此操作,一些示例代码将非常有帮助。)
没关系,我认为解决方案可能只是:
def a_ratio(a_candidate):
a2 = a_candidate[1]
a3 = a_candidate[2]
return a3/a2
feasibility_constraint = optimize.NonlinearConstraint(a_ratio,-k,k)
而不是使用 -k < a3/a2 < k
我会做这两个线性约束:
-k*a2 <= a3
a3 <= k*a2
我知道这个问题已经得到解答,但您想知道如何在 mystic
中处理它。我是 mystic
的作者。像这样:
>>> import mystic as my
>>> import numpy as np
>>>
>>> truth = np.array([1,1,1])
>>>
>>>
>>> class MyObject(object):
... def __init__(self, candidate):
... self.candidate = candidate
... self.truth = truth
... def mean_square_error(self):
... return ((self.truth - self.candidate)**2).sum()
...
>>>
>>> def loss_function(a_candidate):
... return MyObject(a_candidate).mean_square_error()
...
>>>
然后是解决方案:
>>> equations = """
... -k < a3/a2
... a3/a2 < k
... """
>>>
>>> eqn = my.symbolic.simplify(equations, variables=['a1','a2','a3'], locals=dict(k=1.66711), all=False, target='a3')
>>> c = my.symbolic.generate_constraint(my.symbolic.generate_solvers(eqn, variables=['a1','a2','a3']), join=my.constraints.or_)
>>> res = my.solvers.fmin(loss_function, x0=[2,3,4], constraints=c, disp=True, xtol=1e-8)
Optimization terminated successfully.
Current function value: 0.000000
Iterations: 121
Function evaluations: 229
>>> res
array([1., 1., 1.])
mystic
的不同之处在于它构建了一个“运算符”c
,不允许优化器 select 不可行的解决方案。并不是说它用布尔值拒绝不好的解决方案,而是它将“坏”候选者转换为“好”候选者,如下所示:
>>> c([2,10,4])
[2, 10, 4]
>>> c([2,1,4])
[2, 1, 1.6671099999999974]