基于遗传算法的多目标优化算法

一、概述

1、引言

  优化目标可以理解为目标函数，在多目标优化问题中优化目标个数在2个及以上。因此，多目标优化问题和单目标优化相比，最大的区别在于多目标优化问题是一个向量优化（目标函数向量）的问题。而向量之间仅仅存在偏序关系，难以直接比较向量之间的大小，这就导致该类优化问题的性质非常不好。现实问题中，多个优化目标之间或多或少都会存在矛盾。
  如果我们想要选购一款性价比高的手机，自然就要在价格和性能这两个优化目标之间进行权衡。选取手机内存(64g、128g、256g)、尺寸(普通、pro、promax)作为该问题的决策变量，不同的决策变量值就可以构成许多个决策向量=[内存，尺寸]，这些决策向量构成决策空间。显而易见的是，内存越大-价格越高-性能越好；尺寸越高-性能越好-价格越贵。因此，我们无法找到一个确切的解作为该优化问题的最终解，只能获得可行域范围内的解集供用户挑选（易看出，此处的解就是一个决策向量)。

2、多目标优化数学模型

  决策变量、目标函数、约束条件是多目标优化问题的三要素。标准的多目标优化问题的数学模型描述如下：
$$\begin{gathered} \min [f_1(x),f_2(x),\cdot \cdot \cdot ,f_n(x)] \\ s.t.\{\begin{array}{l}lb\le x\le ub\\ g_i(x)\le 0\\ h_j(x)=0\end{array}\text{\hspace{1em}(1.1)} \end{gathered}$$
符号解释：

• $f_i(x)$为待优化的目标函数
• $x$为待优化的变量
• $lb和ub$为变量$x$的上界和下界
• $g_i(x)$为变量$x$的不等式约束
• $h_j(x)$为变量$x$的等式约束

  在下图所示的优化问题中，目标函数$f_1$和$f_2$是相互矛盾的。因为因为$A_1<B_1$且$A_2>B_2$，也就是说，某一个目标函数的提高需要以另一个目标函数的降低作为代价，称这样的解A和解B是非劣解(noninferiority solutions)，或者说是Pareto最优解(Pareto optima)。多目标优化算法的目的就是要寻找这些Pareto最优解。

二、基本概念

2.1 支配(dominate)与非劣(non- inferior)

  在多目标优化问题中,如果个体$p$至少有一个目标比个体$q$的好，而且个体$p$的所有目标都不比个体$q$的差，那么称个体$p$支配个体$q$($p$ dominates $q$)，或者称个体$q$受个体$p$支配($p$ is dominated by $q$) ，也可以说，个体$p$非于个体$q$($p$ is non- inferior to $q$)。

2.2 序值(rank)和前端(front)

  如果$p$支配$q$，那么$p$的序值比$q$的低。如果$p$和$q$互不支配，或者说，$p$和$q$相互非劣，那么$p$和$q$有相同的序值。序值为1的个体属于第一前端，序值为2的个体属于第二前端，依次类推。显然，在当前种群中，第一前端是完全不受支配的，第二前端受第一前端中个体的支配。这样，通过排序，可以将种群中的个体分到不同的前端。

2.3 拥挤距离(crowding distance)

  拥挤距离用来计算某前端中的某个体与该前端中其他个体之间的距离，用以表征个体间的拥挤程度。显然，拥挤距离的值越大，个体间就越不拥挤，种群的多样性就越好。需要指出的是，只有处于同一前端的个体间才需要计算拥挤距离，不同前端之间的个体计算拥挤距离是没有意义的。

2.4 最优前端个体系数(ParetoFraction)

  最优前端个体系数定义为最优前端中的个体在种群中所占的比例，即最优前端个体数=$\min (ParetoFraction\times 种群大小，前端中现存的个体数目)$，其取值范围为0~1。需要指出的是，ParetoFraction的概念是函数gamultiobj所特有的，在NSGA-II中是没有的，这也是为什么称函数gamultiobj(MATLAB自带的基于遗传算法的多目标优化算法，是基于NSGA -II改进的一种多目标优化算法，也是本文下面求解问题所用的算法)是一种多目标优化算法的原因。

三、MATLAB实现

3.1 问题描述

  待优化的多目标问题表述如下：
$$\begin{array}{rl}\min f_1(x_1,x_2)=& x_1^4-10x_1^2+x_1{x}_2+x_2^4+x_1^2{x}_2^2\\ \min f_2(x_1,x_2)=& x_2^4-x_1^2{x}_2^2+x_1^4+x_1{x}_2\\ s.t.\{\begin{array}{l}-5\le x_1\le 5\\ -5\le x_2\le 5\end{array}& \end{array}\text{\hspace{1em}(3.1)}$$

3.2 gamultiobj函数介绍

  gamultiobj函数在：MATLAB安装目录\toolbox\globaloptim中。其使用方法如下：
$$[x,fval]=gamultiobj(fun,nvars,A,b,Aeq,beq,lb,ub,nonlcon,intcon,options)$$
符号解释：

• $x$：返回找到一个局部Pareto集合$x$
• $fval$：返回$x$中所有解定义的所有适应度函数的值。$fval$有nf列，其中nf是目标的数量，与$x$具有相同的行数。
• $fun$：目标函数
  例如：

function f = my_first_multi(x)
 %目标函数
f(1) = x(1)^4 - 10*x(1)^2+x(1)*x(2) + x(2)^4 - (x(1)^2)*(x(2)^2);
f(2) = x(2)^4 - (x(1)^2)*(x(2)^2) + x(1)^4 + x(1)*x(2);

fitnessfcn = @my_first_multi;   % 适应度函数的函数句柄

• $nvars$：变量$x$的个数
• $A,b$：$x$的线性不等式约束，$A\ast x\le b$，没有就用[]
  例如：A = [1,2;3,4;5,6] and b = [10;20;30] 代表
  $x_1+2x_2\le 10$
  $3x_1+4x_2\le 20$
  $5x_1+6x_2\le 30$
• $Aeq,beq$：$x$的线性等式约束，$Aeq\ast x=beq$，没有就用[]
  例如：
  Aeq = [1,2,3;2,4,1] and beq = [10;20]代表
  $x_1+2x_2+3x_3=10$
  $2x_1+4x_2+x_3=20$
• $lb,ub$：$x$的上下界，$lb\le x\le ub$，没有就用[]
• $nonlcon$：$x$的非线性约束，函数$nonlcon$接受$x$并返回向量$c$和$ceq$，分别表示非线性不等式和等式。使c(x)≤0和ceq(x) = 0。
  例如：
  $x_1^2-x_2\le 0$
  $x_1^2+2x_1{x}_2-5\le 0$
  $x_1-x_2^2=0$
  代码如下：

function [c,ceq] = nonlcon(x)
c(1) = x(1)^2-x2;
c(2) = x(1)^2+2*x1*x2-5;
ceq = x(1)-x(2)^2;
end

• $intcon$：$x$的整数约束，要求$intcon$中列出的变量取整数值。当$intcon$为非空时，对于$nonlcon$，ceq必须返回空值。
• $options$：设置优化参数

这里只是做简单介绍，详细介绍参考MATLABA官方

3.3 代码实现

  (1)使用函数gamultiobj求解多目标优化问题的第一步是编写目标函数的M文件。对
于以上问题，函数名为my_first_multi，目标函数代码如下:

function f = my_first_multi(x)
 %目标函数
f(1) = x(1)^4 - 10*x(1)^2+x(1)*x(2) + x(2)^4 - (x(1)^2)*(x(2)^2);
f(2) = x(2)^4 - (x(1)^2)*(x(2)^2) + x(1)^4 + x(1)*x(2);

  (2)使用命令行方式调用gamultiobj函数，代码如下:

clear
clc
fitnessfcn = @my_first_multi;   % Function handle to the fitness function
nvars = 2;                      % Number of decision variables
lb = [-5,-5];                   % Lower bound
ub = [5,5];                     % Upper bound
A = []; b = [];                 % No linear inequality constraints
Aeq = []; beq = [];             % No linear equality constraints
options = gaoptimset('ParetoFraction',0.3,'PopulationSize',100,'Generations',200,'StallGenLimit',200,'TolFun',1e-100,'PlotFcns',@gaplotpareto);

[x,fval] = gamultiobj(fitnessfcn,nvars, A,b,Aeq,beq,lb,ub,options);

  其中，fitnessfcn即(1)中定义的目标函数M文件，设置的最优前端个体系数ParetoFrac-tion为0.3，种群大小PopulationSize为100，最大进化代数Generations为200，停止代数StallGenLimit 也为200，适应度函数值偏差TolFun为$e^{-100}$，绘制Pareto前端。

3.4 结果分析

  可以看到，在基于遗传算法的多目标优化算法的运行过程中，自动绘制了第一前端中个体的分布情况，且分布随着算法进化一代而更新一次。当选代停止后，得到如下图所示的第一前端个体分布图。同时，Worksapce中返回了函数gamultiobj得到的Pareto解集x及与x对应的目标函数值，如下表所列。需要说明的是，由于算法的初始种群是随机产生的，因此每次运行的结果不一样。

序号	$x_1$	$x_2$	$f_1$	$f_2$
1	-0.707106331	0.707073	-5.24999	-0.25
2	0.832482315	-0.83055	-7.14362	-0.21335
3	-0.707106331	0.707073	-5.24999	-0.25
4	2.671846723	-1.97666	-38.3334	33.05425
5	1.796546904	-1.46041	-26.8172	5.458635
6	2.409171646	-1.8517	-36.9589	21.08214
7	2.368949823	-1.72953	-36.5618	19.55743
8	1.702078369	-1.53457	-24.4664	4.504309
9	0.989973221	-0.73077	-9.80161	-0.00114
10	2.220327613	-1.5795	-34.5771	14.72148
11	2.340450713	-1.83376	-36.1758	18.6013
12	1.643693127	-1.48438	-23.2559	3.76141
13	2.096188757	-1.6257	-32.6685	11.27152
14	1.078856237	-0.8303	-11.5075	0.131821
15	2.173833905	-1.67333	-33.9537	13.30179
16	2.638950538	-1.96228	-38.3096	31.33104
17	2.012820536	-1.63771	-31.0694	9.44507
18	1.609579116	-1.07152	-22.5765	3.330951
19	2.053805095	-1.35124	-31.5318	10.64938
20	1.319236996	-1.3451	-16.0248	1.37911
21	2.259777888	-1.46052	-34.6319	16.4341
22	1.889380772	-1.56903	-28.6464	7.051187
23	1.813599357	-1.40256	-27.2172	5.674235
24	1.744700635	-1.37822	-25.7525	4.687286
25	2.539551817	-1.95918	-37.8967	26.59656
26	1.961798584	-1.5357	-30.2018	8.284762
27	1.42453827	-1.09142	-18.7281	1.564972
28	2.460348783	-1.78242	-37.414	23.11921
29	2.601732392	-1.95616	-38.2195	29.47064
30	1.200014386	-1.02306	-13.9661	0.434279

  从上图可以看到，第一前端的Pareto最优解分布均匀。从上表可以看到，返回的Pareto最优解个数为30个，而种群大小为100，可见，ParetoFraction为0.3的设置发挥了作用。另外，个体被限制在了[-5,5]的上下限范围内。