Mystic差分进化算法约束参数设置问题求助

问题背景 #

尝试用diffevo（差分进化）算法优化带强非线性约束的平滑目标函数，但不清楚如何定义constraints、bounds、init vectors等参数，使用Mystic时遇到约束相关错误，SciPy的多种算法也未成功。

目标函数（预算总和最小化） #

def minimize_me_reach(budgetA, *args):
    # start = time.time()
    # budget = args[0]
    # # trp_total = args[0]
    # freq = args[1]
    # DF_GROUP = args[2]
    # MEMBER = args[3]
    # DF_GROUP['bud_part'] = 0
    # DF_GROUP.loc[DF_GROUP['Prime'] == 1, 'bud_part'] = budgetA
    # DF_GROUP.loc[DF_GROUP['Prime'] == 2, 'bud_part'] = budgetA
    # DF_GROUP['budget_total'] = DF_GROUP['bud_part'] * budget
    # execution of above program is :",(end - start), "sec")
    return np.nan_to_num(budgetA).sum()

约束函数（当前1个，后续计划添加更多） #

def constr_reach(x_all):
    #budget = x_all[-1]
    x = x_all[:-1]
    print(DF_GROUP.info())
    DF_GROUP['bud_part'] = 0
    print(MEMBER.info())
    DF_GROUP.loc[DF_GROUP['Prime'] == 1, 'bud_part'] = x
    DF_GROUP.loc[DF_GROUP['Prime'] == 2, 'bud_part'] = x
    DF_GROUP['budget_total'] = DF_GROUP['bud_part']
    DF_GROUP['trp'] = DF_GROUP['budget_total'] / DF_GROUP['tcpp'] * DF_GROUP['trp_mult']
    DF_GROUP['Q'] = pd.DataFrame([DF_GROUP['trp'] / DF_GROUP['Tvr_Origin'], DF_GROUP['COUNT'] - 1]).min()
    test1 = tv_planet.calc_reach(MEMBER, DF_GROUP)
    sOPT = test1.loc[freq - 1, '%']
    return sOPT

尝试的两种调用方式 #

方式一：使用DifferentialEvolutionSolver2 #

stepmon = VerboseMonitor(50)

solver = DifferentialEvolutionSolver2(dim=53,NP=len(init_bud))
solver.SetConstraints(constr_reach)
#solver.SetInitialPoints(x0=x0, radius=0.1)
solver.SetGenerationMonitor(stepmon)
solver.enable_signal_handler()
solver.Solve(minimize_me_reach, termination=VTR(0.01), strategy=Best1Exp, \
                CrossProbability=1.0, ScalingFactor=0.9, \
                disp=True, ExtraArgs=(reach, freq, DF_GROUP, MEMBER, days))
result = solver.Solution()
iterations = len(stepmon)
cost = stepmon.y[-1]
print("Generation %d has best Chi-Squared: %f" % (iterations, cost))
print(result)

方式二：使用diffev函数 #

result = diffev(minimize_me_reach, x0=x0, npop=len(init_bud), args=(reach, freq, DF_GROUP, MEMBER, days),
                bounds=my_bounds, ftol=0.05, gtol=100,
                maxiter=100, maxfun=None, cross=0.9, scale=0.8, full_output=False,
                disp=1, retall=0, callback=None, strategy=Best1Exp, constraints=constr_reach)

遇到的错误信息 #

File "d:\tvplan\tv-planner-backend\TV_PLANER2\top_to_bottom.py", line 1179, in optimizeReach_DE
    solver.Solve(minimize_me_reach, termination=VTR(0.01), strategy=Best1Exp, \
  File "C:\Users\Mi\AppData\Local\Programs\Python\Python310\lib\site-packages\mystic\differential_evolution.py", line 661, in Solve
    super(DifferentialEvolutionSolver2, self).Solve(cost, termination,\
  File "C:\Users\Mi\AppData\Local\Programs\Python\Python310\lib\site-packages\mystic\abstract_solver.py", line 1180, in Solve
    self._Solve(cost, ExtraArgs, **settings)
  File "C:\Users\Mi\AppData\Local\Programs\Python\Python310\lib\site-packages\mystic\abstract_solver.py", line 1123, in _Solve
    stop = self.Step(**settings) #XXX: remove need to pass settings?
  File "C:\Users\Mi\AppData\Local\Programs\Python\Python310\lib\site-packages\mystic\abstract_solver.py", line 1096, in Step
    self._Step(**kwds) #FIXME: not all kwds are given in __doc__
  File "C:\Users\Mi\AppData\Local\Programs\Python\Python310\lib\site-packages\mystic\differential_evolution.py", line 552, in _Step
    self.trialSolution[candidate][:] = constraints(self.trialSolution[candidate])
TypeError: can only assign an iterable

解决建议与示例 #

核心问题分析 #

报错的原因是Mystic的约束函数要求返回与输入维度一致的可迭代对象（用来修正候选解），但你的constr_reach返回的是单个数值，不符合要求。需要将约束定义为等式/不等式约束，或者构造能修正候选解的函数。

示例：带非线性约束的差分进化优化 #

import numpy as np
import mystic as my
from mystic.monitors import VerboseMonitor
from mystic.solvers import DifferentialEvolutionSolver2
from mystic.strategy import Best1Exp
from mystic.constraints import as_constraint, inequality

# 目标函数：最小化变量求和
def objective(x):
    return np.sum(x)

# 非线性约束：x[0]^2 + x[1]^2 >= 4
def constraint_func(x):
    return x[0]**2 + x[1]**2 - 4  # 返回值>=0表示满足约束

# 构造约束对象
@as_constraint(inequality(constraint_func))
def constrained_solution(x):
    return x

# 定义变量边界：每个变量范围[0,5]
bounds = [(0, 5)] * 2

# 初始化求解器
stepmon = VerboseMonitor(10)
solver = DifferentialEvolutionSolver2(dim=2, NP=20)
solver.SetConstraints(constrained_solution)
solver.SetStrictRanges(bounds)
# 初始化种群：随机生成20个候选解
solver.SetInitialPoints(np.random.uniform(0, 5, size=(20, 2)))

# 执行求解
solver.Solve(objective, termination=my.termination.VTR(0.001),
             strategy=Best1Exp, CrossProbability=0.9, ScalingFactor=0.8,
             disp=True)

# 输出结果
print("最优解:", solver.Solution())
print("目标函数值:", objective(solver.Solution()))

针对你的问题的调整建议 #

1. 修正约束函数定义
   假设你需要sOPT >= 目标reach值，修改约束函数并包装成Mystic可识别的约束对象：

   from mystic.constraints import as_constraint, inequality
   
   def constraint_reach(x, target_reach, freq, DF_GROUP, MEMBER):
       # 操作DF_GROUP的副本，避免修改全局数据污染种群迭代
       df_copy = DF_GROUP.copy()
       df_copy['bud_part'] = 0
       df_copy.loc[df_copy['Prime'] == 1, 'bud_part'] = x
       df_copy.loc[df_copy['Prime'] == 2, 'bud_part'] = x
       df_copy['budget_total'] = df_copy['bud_part']
       df_copy['trp'] = df_copy['budget_total'] / df_copy['tcpp'] * df_copy['trp_mult']
       df_copy['Q'] = pd.DataFrame([df_copy['trp'] / df_copy['Tvr_Origin'], df_copy['COUNT'] - 1]).min()
       test1 = tv_planet.calc_reach(MEMBER, df_copy)
       sOPT = test1.loc[freq - 1, '%']
       return sOPT - target_reach  # 返回值>=0表示满足约束
   
   # 构造带参数的约束对象
   def make_constraint(target_reach, freq, DF_GROUP, MEMBER):
       def _constraint(x):
           return constraint_reach(x, target_reach, freq, DF_GROUP, MEMBER)
       return as_constraint(inequality(_constraint))
   
   # 使用时传入参数生成约束
   constrained = make_constraint(reach, freq, DF_GROUP, MEMBER)
   solver.SetConstraints(constrained)
   

2. 设置变量边界
   定义每个变量的取值范围，例如：

   # 假设53个变量的范围都是[0, 1000]
   my_bounds = [(0, 1000)] * 53
   # 方式一中设置边界
   solver.SetStrictRanges(my_bounds)
   # 方式二中直接传入bounds参数即可
   

3. 初始化种群
   如果有初始候选解x0（形状为(npop, 53)），可以通过solver.SetInitialPoints(x0)设置；如果没有，求解器会自动随机生成符合边界的初始种群。

内容的提问来源于stack exchange，提问作者Pavel_K