Matlab多目标函数折中求解需求及GA实现代码

问题背景 #

已通过Matlab的gamultiobj函数完成双目标遗传算法优化，得到帕累托最优前沿。现需找到折中解决方案，使两个目标函数各自达到其最优值的50%。

现有优化代码 #

目标函数定义 #

function [output]= multi_objective_function_2(input)
    x1 = input(1);     
    x2 = input(2);    
    x3 = input(3);    
    x4 = input(4);    
    x5 = input(5);    
    x6 = input(6);    
    x7 = input(7);   

    % 原目标为最大化以下两个线性组合，此处取负转为最小化问题适配gamultiobj
    f1 = -(1.0125*x1 + 14.8125*x2 + 1.3113*x3 + 0.783*x4 + 10.365*x5 + 2.435*x6 + 1.78*x7);   
    f2 = -(3303*x1 + 7510*x2 + 41370*x3 + 1628*x4 + 23046*x5 + 5412*x6 + 58113*x7);           
    output = [f1 f2];
end

主优化脚本 #

clear;
clc;
tic
numberOfVariables = 7;
% 不等式约束 A*x <= b
A = [-1 0 0 0 0 0 0; 
     0 -1 0 0 0 0 0; 
     0 0 -1 0 0 0 0; 
     0 0 0 -1 0 0 0; 
     0 0 0 0 -1 0 0; 
     0 0 0 0 0 -1 0; 
     0 0 0 0 0 0 -1; 
     1715.78 5586.32 1654.93 1963.37 4414.64 6973.88 4416.31; 
     -1.0125 -14.8125 -1.3113 -0.783 0 0 0];
b = [-3; -45; -2445; -71; -10402.5; -5201; -10402.5; 208000000; -3931];
% 等式约束 Aeq*x = beq
Aeq = [1 1 1 1 0 0 0; 
       0 0 0 0 1 1 1];
beq = [26006.25; 26006.25];
% 变量上下界
LB = [0 0 0 0 0 0 0];
UB = [26006.25 26006.25 26006.25 26006.25 26006.25 26006.25 26006.25];
% 配置参数，绘制帕累托前沿
options = gaoptimset('PlotFcns',@gaplotpareto);
% 执行多目标优化
[x,fval,exitflag,output] = gamultiobj(@multi_objective_function_2,numberOfVariables,A,b,Aeq,beq,LB,UB,options);
% 输出结果
disp('     x1        x2        x3         x4          x5         x6         x7  ');
disp(x)
disp('         F1      F2');
disp(fval)
disp(output)
fprintf('The number of points on the Pareto front was: %d\n', size(x,1));
toc

50%折中解决方案实现 #

核心思路 #

要实现两个目标各达50%最优值，需先获取每个目标的单独最优值，再通过加权求和法或帕累托前沿筛选找到满足条件的解：

1. 单独求解每个目标的最优值（最大化原目标，即最小化f1/f2）；
2. 以50:50的权重组合两个目标函数，重新优化得到折中解；
3. 验证解是否满足两个目标均达到各自最优值的50%。

完整实现代码 #

clear;
clc;
tic
numberOfVariables = 7;
A = [-1 0 0 0 0 0 0; 
     0 -1 0 0 0 0 0; 
     0 0 -1 0 0 0 0; 
     0 0 0 -1 0 0 0; 
     0 0 0 0 -1 0 0; 
     0 0 0 0 0 -1 0; 
     0 0 0 0 0 0 -1; 
     1715.78 5586.32 1654.93 1963.37 4414.64 6973.88 4416.31; 
     -1.0125 -14.8125 -1.3113 -0.783 0 0 0];
b = [-3; -45; -2445; -71; -10402.5; -5201; -10402.5; 208000000; -3931];
Aeq = [1 1 1 1 0 0 0; 
       0 0 0 0 1 1 1];
beq = [26006.25; 26006.25];
LB = [0 0 0 0 0 0 0];
UB = [26006.25 26006.25 26006.25 26006.25 26006.25 26006.25 26006.25];

% --------------------------
% 步骤1：求解单个目标的最优值
% --------------------------
options_single = gaoptimset;
% 求解目标1的最大值（对应最小化f1）
[x1_opt, f1_opt] = ga(@(x) multi_objective_function_2(x)(1), numberOfVariables, A, b, Aeq, beq, LB, UB, options_single);
obj1_max = -f1_opt;
% 求解目标2的最大值（对应最小化f2）
[x2_opt, f2_opt] = ga(@(x) multi_objective_function_2(x)(2), numberOfVariables, A, b, Aeq, beq, LB, UB, options_single);
obj2_max = -f2_opt;

fprintf('目标1最优值：%.4f\n目标2最优值：%.4f\n', obj1_max, obj2_max);
fprintf('目标1 50%阈值：%.4f\n目标2 50%阈值：%.4f\n', 0.5*obj1_max, 0.5*obj2_max);

% --------------------------
% 步骤2：求解50%权重的折中解
% --------------------------
% 定义加权目标函数：0.5*f1 + 0.5*f2，最小化该组合
折中目标函数 = @(x) 0.5*multi_objective_function_2(x)(1) + 0.5*multi_objective_function_2(x)(2);
options_compromise = gaoptimset;
[x_compromise, f_compromise] = ga(折中目标函数, numberOfVariables, A, b, Aeq, beq, LB, UB, options_compromise);

% 计算折中解对应的原目标值
obj1_compromise = -f_compromise(1);
obj2_compromise = -f_compromise(2);

% --------------------------
% 步骤3：输出并验证结果
% --------------------------
disp('==================== 折中解决方案 ====================');
disp('     x1        x2        x3         x4          x5         x6         x7  ');
disp(x_compromise);
disp('         原目标1值      原目标2值');
disp([obj1_compromise, obj2_compromise]);
disp('         目标1达成率(%)  目标2达成率(%)');
disp([obj1_compromise/obj1_max*100, obj2_compromise/obj2_max*100]);

toc

帕累托前沿筛选补充 #

若需直接从已得到的帕累托前沿中筛选最接近50%达成率的点，可使用以下代码（需先运行原多目标优化脚本）：

% 基于已得到的帕累托前沿x和fval筛选
obj1_values = -fval(:,1);
obj2_values = -fval(:,2);
% 计算每个点的目标达成率
rate1 = obj1_values / obj1_max;
rate2 = obj2_values / obj2_max;
% 寻找最接近(0.5,0.5)的点
distance = sqrt((rate1-0.5).^2 + (rate2-0.5).^2);
[~, idx] = min(distance);
x_compromise_pareto = x(idx,:);
obj1_compromise_pareto = obj1_values(idx);
obj2_compromise_pareto = obj2_values(idx);

disp('==================== 帕累托前沿筛选的折中解 ====================');
disp(x_compromise_pareto);
disp([obj1_compromise_pareto, obj2_compromise_pareto]);
disp([rate1(idx)*100, rate2(idx)*100]);

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