0 引言

现实生活中涌现出越来越多的多目标优化问题(Multi-objective Optimization Problem，MOP)，但是这些目标之间通常是相互制约的，一个目标在得到改善的同时会伴随着其他目标的效果变差，故MOP得到的解不能同时满足全部目标最优，而是一组折衷的解，即Pareto解集或非劣解集^[1]。基于群体搜索的进化算法(Evolutionary Algorithm，EA)是模拟生物进化过程的自组织、自适应的人工智能技术，主要通过选择、交叉、变异、重组等遗传操作求解问题，其最大的优势就是不要求问题满足连续、可微等特定性质。基于EA求解MOP问题的优势，短短数年内涌现出大量优秀的多目标进化算法(Multi-objective Evolutionary Algorithm，MOEA)。

1985年由Schaffer提出的矢量评价遗传算法VEGA^[2]是最早的进化算法。在此后的时间里，学者针对不同的优化问题提出不同的MOEA，其中最具代表性的算法是以Pareto支配为基础的MOEA，如NSGA^[3]、NSGA-2^[4]、强度Pareto进化算法SPEA^[5]、SPEA2^[6]和PESA-Ⅱ^[7]等。上述算法在求解较少数目标优化问题上表现出良好的性能，其中以NSGAⅡ最具代表性。

多目标优化问题在现实生活中的应用比比皆是，例如：电力供应^[8]、海水淡化^[9]、护士排班问题^[10]、汽车控制器优化问题^[11]以及疫情救灾点优化^[12]等。这些实际应用问题的目标数通常是大于等于4的，所以这些问题又被称为高维多目标优化问题(Many-objective Optimization Problem，MaOP)。MaOP在传统的优化过程会面临很多困难，如：①随着目标数的增加，种群中非支配个体所占的比例会急剧增大，严重削弱Pareto选择压力；②在目标空间中可能会出现偏远的解，即支配抵触解(Dominance Resistance Solution，DRS)^[13]；③在高维空间中，因为个体可能相距太远，父代产生的子代有很大的差异，导致后代不稳定且不具有代表性，通常会面临杂交和变异算子失效的情况^[14]；④Pareto前沿面表示困难。假设有m个目标的优化问题，每个目标上有b个解，则共需要mb^m-1个解表示，这在现实应用中十分困难。MaOP算法的提出一般基于上述几个问题，提出针对性的解决方案。

传统的NSGAⅡ算法在高维(目标数大于等于4)目标空间优化时，算法性能急剧下降，且算法的多样性表现较差，不适用于现在的多目标优化问题。受Elarbri等^[15]所提出的RPD-NSGAⅡ启发，本研究采用d₂距离来增强算法在高维空间的多样性，提出算法d₂-NSGAⅡ。新算法拟在保证NSGAⅡ原有精英策略等优点的同时，显著改善算法的分布性。

1 多目标优化问题及相关概念

统一以最小化问题为例，将m个目标函数和n个决策变量，(p+q) 个约束的MOP描述为

$ \left\{ \begin{array} { l } { \operatorname { min } y = F ( x ) = ( f _ { 1 } ( x ) , f _ { 2 } ( x ) , \cdots , f _ { m } ( x ) ) } \\ { x = ( x _ { 1 } , x _ { 2 } , \cdots , x _ { n } , ) } \\ { y = ( y _ { 1 } , y _ { 2 } , \cdots , y _ { n } , ) } \\ { s . t . g _ { i } ( x ) \geq 0 , i = 1, 2 , \cdots , p } \\ { h _ { j } ( x ) = 0 , j = 1, 2 , \cdots , q } \end{array} \right.，$

(1)

式中x∈X⊂Rⁿ，y∈Y⊂R^m，n是决策变量的个数，m是目标向量的个数，x表示决策向量，X是由决策向量x形成的决策空间；y表示目标向量，Y是由目标向量形成的目标空间；F(x)定义了由决策空间X向目标空间Y映射的函数，g(x)定义了p个不等式约束，h(x)定义了q个不等式约束。

定义  1     (可行解集)满足式(1)中约束函数的决策向量的集合，用X_l表示，X_l={ x∈X g(x)≥0∧h(x)=0}。

定义  2     (Pareto支配)存在x₁, x₂∈X_l是问题的解，当满足∀i∈(1, …M): f_i(x₁)≤f_i(x₂)∧∃j∈(1, …M): f_i(x₁) < f_i(x₂)时候，称x₁ Pareto支配x₂(记为x₁$ \prec $x₂)。

定义  3     (Pareto最优解)当不存在其他解x∈X_l，使得x$ \prec $x_l，称x_l为最优解，最优解也叫做非支配解。

定义  4     (Pareto最优解集)满足定义3的解的集合称为Pareto最优解集(Pareto optimal set，PS)，表示为PS={x_l}={x∈X_l|-∃x′∈X_l: x′$ \prec $x}。

定义  5     (Pareto前沿) PS在目标空间的投影叫做Pareto前沿(Pareto optimal front，PF)，表示为PF={F(x)∈R^m|x∈PS}。

多目标进化算法按照不同优化标准可以分为3类：基于支配的MOEA、基于指标的MOEA以及基于分解的MOEA。其中，基于支配的MOEA算法代表有NSGAⅡ和SPEA2，基于分解的MOEA算法代表有MOEA/D^[16]、MOEA/D-URAW^[17]以及MOEA/D-DRA^[18]，基于指标的MOEA算法代表有超体积评价指标(Hypervolume，HV)^[19]、MOEA-HV^[20]以及基于R2指标的R2-HVC^[21]。在基于分解的MOEA算法中，Zhang与Li^[16]于2007年提出了一种优秀的PBI分解方法，该方法因为可以较好地控制种群的收敛性和分布性而广受关注。PBI在应用时需要度量两个距离d₁和d₂，d₁代表解在参考向量方向投影距离原点的欧式距离，d₂代表解距离参考向量的欧式距离，具体度量如图 1所示。通过图示不难看出，一个较小的d₁代表解距离PF近，个体的收敛性更好；同样较小的d₂代表解具有更好的多样性，当d₂=0时，意味着解在参考向量的方向上。

具体数学表示如公式(2)：

$ \left\{ \begin{array}{l} {d_1}\left( x \right) = \left\| {f{{\left( x \right)}^{\rm{T}}}\lambda } \right\|/\left\| \lambda \right\|\\ {d_2}\left( x \right) = \left\| {f\left( x \right){\rm{ - }}{{\rm{d}}_1}\left( x \right)\left( {\frac{\lambda }{{\left\| \lambda \right\|}}} \right)} \right\| \end{array} \right., $

(2)

其中f(x)= f₁(x), f₂(x), …, f_m(x)^T是解x的标准目标向量，PF是真实前沿面，λ是m维的方向向量。

2 基于d₂距离的改进型NSGAⅡ

传统的NSGAⅡ算法在保持种群多样性方面采用的是度量拥挤距离，该方法存在一定缺陷，即欧氏距离的度量在高维空间中不再适用^[22]。个体的拥挤距离存在不能准确度量个体拥挤度的问题, 在种群迭代的时候，该问题会导致分布性差的个体反而被保留，进而导致种群多样性较差。如图 2所示，图中个体为最后非支配层的个体，个体B、C的拥挤距离较大且相等，在迭代的时候这两个个体会被保留下来。但是，由图 2可以看出这两个个体实际距离较为接近，不应该全部被保留；此时根据拥挤距离不能确定个体，可以度量两个个体距离相应权重向量的d₂距离，选择距离小的个体保留(保留个体B)。

设置对照实验，增加统一数据集的维数，计算出目标空间中最近与最远个体之间的差值d_max-d_min，实验参数设置如下：目标数M=30，决策变量数目D=30，种群数N=100，评价次数Evaluation=10 000，p分别取值2.0，1.0，0.5，0.1。通过实验数据(图 3)可以看出，当p=2.0(欧氏距离)以及p=1.0(曼哈顿距离)时，在高维空间d_max与d_min差值不是很明显; 当p=0.5以及p=0.1时，随维数的增加，d_max与d_min差值明显变大，有利于在高维空间中选择分布性较好的个体保留。因此，欧式距离降低了在高维空间中多样性的选择能力，其在高维空间中意义不大。相反地，其他p范式距离(p < 1) 在高维空间上表现更好。

Pareto支配能更好地提高算法的收敛性，而PBI效用函数的d₂距离可以提高种群的多样性分布。为将二者结合并求解高维目标优化问题，本研究将d₂距离的度量应用到NSGAⅡ框架中。同样采用NSGAⅡ非支配排序的方法选择优秀个体，保留F₁, F₂…F_l-1非支配层个体，删除F_l层以后的个体，在第F_l非支配层利用个体距离相应参考向量的d₂距离，选择出剩余优秀个体，使得种群规模保持不变。这样在保证种群收敛性的同时，促进种群的多样性。由此，产生一种基于d₂距离的新算法d₂_NSGAⅡ，算法的示意如图 4所示。

算法中使用的Pareto支配满足反自反、反对称和传递性质。

① 反自反性：种群中任一个体x，满足反自反性，即x-$ \prec $x。

证明：如果x$ \prec $y，则至少存在一个分量i=1, 2, …, m: 满足f_i(x) < f_i(y)；但是Pareto定义中不存在f_i(x) < f_i(y)，所以x不能Pareto支配其本身，即x满足反自反性。

② 反对称性：种群中任意两个体x，y如果满足x$ \prec $y，则y-$ \prec $x。

证明：假设x$ \prec $y，根据定义3，则任意f_i(x), f_i(y), i=1, 2, …, m: 存在f_i(x)≤f_i(y)。所以不存在f_i(x)>f_i(y)，故y不能Pareto支配x，即y-$ \prec $x，满足反对称性质。

③ 传递性：x, y, z是种群中3个个体，如果x$ \prec $y且y$ \prec $z, 则x$ \prec $z。

证明：因为x$ \prec $y，则满足∀i=1, 2, …, m: f_i(x)≤f_i(y)∧∃j=1, 2, …, m: f_j(x) < f_j(y)；同理y$ \prec $z，满足∀i=1, 2, …, m: f_i(y)≤f_i(z)∧∃j=1, 2, …, m: f_j(y) < f_j(z)，因此对于个体x，z，至少存在一个分量i=1, 2, …, m: f_i(x)≤f_i(z)∧∃j=1, 2, …, m: f_j(x) < f_j(z)，即x$ \prec $z成立，个体满足传递性。

除非支配排序外，文章还涉及参考点的产生。因文章算法针对高维多目标优化问题，因此采用Deb和Jain的双层参考点产生方法，产生两层参考点，每层参考点数量计算公式相同，具体如下所示：

$ C = \left( \begin{array}{l} H + M - 1\\ M - 1 \end{array} \right), $

(3)

其中，C代表参考点数量，M代表目标数量，H代表超平面上沿着每个目标考虑的划分数量。

假设在一个3目标优化问题上，在外层参考点上，取H₁=2，则外层参考点数目C₁= $ \left( \begin{array}{l} {H_1} + M - 1\\ M - 1 \end{array} \right)$=6；在内层参考点上，取H₂=1，则内层参考点数目C₂= $ \left( \begin{array}{l} {H_2} + M - 1\\ M - 1 \end{array} \right)$=3；总的参考点数目C=C₁+C₂，产生参考点的过程如图 5所示。在目标数量M=3，H₁=2，H₂=1，两层参考点共(6+3=9)个，双层参考点的具体分布如图 5c所示。不同目标数量M、不同划分数量H的情况下，参考数量多少如表 1所示。

d₂_NSGAⅡ具有以下特点：①利用NSGAⅡ经典框架，保留精英策略，同时采用非支配排序，算法时间复杂度低；②算法最后的非支配层采用d₂距离度量，可以在收敛性的基础之上最大化地提高多样性。以下给出算法的具体流程。

算法：d₂_NSGAⅡ算法

输入：种群规模N，决策变量数目D，MaOP问题的目标数目M，最大迭代次数Tmaxgen，沿目标方向划分H份

输出：最末代种群P_maxgen

1.初始化

  1.1初始化迭代计数器t=0；

  1.2在MaOP问题的可行决策空间内随机产生N个初始点，形成初始化种群P₀；

  1.3计算P₀中所有解的目标值向量F₁(0), …, F_N(0)；

  1.4产生参考点(M, H)→Set；

2.WHILE(t < T_maxgen)

3.构建交配池：P_t^mating=Mating_selection(P_t)；

4.重组运算：P_t^{recombination}=Recombination(P_t^mating)；

5.变异运算：P_t^offspring=Mutation(P_t^{recombination})；

6.合并子种群和父种群：R_t=P_t∪Q_t；

7.计算最小向量值→F_min；

8.计算最大向量值→F_max；

9.归一化处理(R_t, F_min, F_max)→R_t；

10.保留R_t中边界解；

11.对R_t集合中的个体进行非支配排序；

12.确定最小的k值，使其满足|F₁∪F₂∪…∪F_k|≥N; //F_i表示第i层非Pareto支配层；

13.计算非支配层(F_k)个体的d₂距离d₂(F_k, Set)；

14.关联F_k中个体与Set集合中的点；

15.IF |F₁∪F₂∪…∪F_k|>N THEN

16.根据距离从非Pareto支配层F_k中删除(|F₁∪F₂∪…∪F_k|-N)个具有较大d₂距离的个体；

17.END IF

18.P_t+1=|F₁∪F₂∪…∪F_k|；//获得规模为N的下一代种群;

19.更新迭代计数器t=t+1；

20.END WHILE.

21.输出最末代种群P_maxgen。

算法的第1步为初始化，随机产生个体，并计算其目标值大小，同时产生参考点；第2步进入循环，进行重组变异，产生新的后代；第3步Mating_selection利用二元竞标赛法选取优秀的N/2个体；第4步对交配池中的个体进行仿二进制交叉(SBX)，得到规模为N的种群P_t^{recombination}；第5步执行变异操作，得到后代Q_t，规模为N。然后进行父种群以及子代种群的合并，生成规模为2N的新种群。为了便于目标值的比较，对目标值进行归一化处理。

归一化的目的是为了消除个体极端目标值，将数据标准化，便于比较，使实验数据更加可靠真实。首先，确定种群中个体每个目标分量f_i的最小值f_i^min，构造向量F^min= f₁^min, f₂^min, …, f_M^min T，同样构造目标最大值的向量F^max= f₁^max, f₂^max, …, f_M^max T，其中f_i^max是种群R_t中f_i的最大值。然后将R_t中所有个体目标值归一化，具体归一化公式如下所示：

$ {{\bar f}_i}(x) = \frac{{{f_i}(x) - f_i^{{\rm{ min }}}}}{{f_i^{{\rm{ max }}} - f_i^{{\rm{ min }}}}}{\rm{ }}。$

(4)

算法第13步通过公式(2)计算F_k层个体距离参考向量的d₂距离，找到距离个体最近的参考向量进行关联。如图 6所示，度量个体a、b、c、e所有权重向量的d₂距离，不难发现，个体b距离参考向量λ₁最近，个体c距离参考向量λ₂最近，将F_k层个体与相应权重向量关联。算法中，先进性非支配排序获得若干非支配层后，再对F_k非支配个体进行关联参考点，这样可以避免计算所有个体的d₂距离，减少计算复杂度。

同时，保留目标空间中的边界值，将在每个目标上的最大值与最小值的个体保留下来(边界值)。将边界个体的d₂值设置为0，在迭代的时候将d₂值为0的个体保留下来。与NSAGⅡ同样地进行非支配排序，确定最后的非支配层使得种群规模保持N，并在最后非支配层根据d₂值由小到大筛选优秀的个体，使种群的规模保持为N，进入下一代迭代。

d₂_NSGAⅡ的时间复杂度主要由两部分组成：非支配排序、d₂距离度量。假设MOP有M个目标，种群规模N。非支配排序中，种群中所有个体之间相互比较，确定所有个体计数S_q(被个体q支配的个体数目) 和个体q支配的个体数目n_q，需要两两个体比较，时间复杂度为O(MN²)。在找到第一非支配层个体后，剩余个体的n_q计数最多为N-1，故剩余个体最多需要遍历(N-1)次才会找到所在非支配层。因为除去第一非支配层之外最多剩余个体为(N-1)，故此过程时间复杂度为O(N²)。因此，非支配排序的时间复杂度为O(MN²)。在最差情况下，即第一非支配层即为最后非支配层，需度量所有个体距离参考向量的距离，时间复杂度为O(N²)。综上两部分时间复杂度比较，取最高次幂，故d₂_NSGAⅡ的时间复杂度为O(MN²)。

3 测试实验与数据分析

为验证所提出的d₂_NSGAⅡ算法的有效性，采用3组对比算法在DTLZ系列测试函数上进行实验数据的对照以及分析，所采用的算法依次是NSGAⅡ、NSGA、VaEA^[23]、MOEAD以及MOEADFFRMAB^[24]。

3.1 测试问题

研究采用的是DTLZ系列测试函数，具体为DTLZ1-DTLZ7。DTLZ测试函数的前沿面具有多样性，具备不同的特征，可以从多个角度检验算法在高维多目标情况下优化问题的能力，而且采用的测试函数在目标空间和决策空间上均具有可扩展性。DTLZ系列测试函数的具体特征描述如表 2所示。

3.2 性能指标

采用IGD指标测试新算法的性能。该指标反映了真实的Pareto前沿到近似Pareto前沿之间的距离，可以综合度量算法的收敛性与多样性。IGD指标具体工作机制如下：在真实的Pareto前沿上采用一组分布较为均匀的点，并计算他们到近似Pareto解之间的距离。IGD指标的值越小表明算法收敛性与多样性越好。具体的计算公式如下：

$ {\rm{IGD}} = \frac{1}{{|P|}}\sum\limits_{i = 1}^{|P|} { Dist{ _i}} , $

(5)

$ {\rm{Dis}}{{\rm{t}}_i} = \mathop {min}\limits_{j = 1}^{\left| A \right|} \sqrt {{{\sum\limits_{m = 1}^M {(\frac{{{f_m}({p_i}) - {f_m}({a_j})}}{{f_m^{ max } - f_m^{ min }}})} }^2}} , $

(6)

上述公式中，P是MOP问题的解集；Dist_i是目标函数归一化处理之后的最小欧氏距离；f_m^max和f_m^min分别表示P在第m个目标上获得的最大值与最小值。A为近似Pareto前沿，p_i∈P, i=1, 2, …, |P|；a_j∈A, j=1, 2, …, |A|。在对比实验中，实验进行30次，取结果平均值作为最后的实验结果；同时利用显著水平为0.05的Wilcoxon秩和检验分析所有对比算法结果的差异。

3.3 实验结果与分析

实验对比了6种算法在DTLZ6测试函数上，目标数为8的IGD值变化情况(图 7)。从图 7可以看出，随着评估次数的增大，6种算法的IGD值整体在变小。评估次数在3 000之前，d₂_NSGAⅡ的下降速度远远大于NSGAⅡ、VaEA、NSGAⅢ，但是稍微小于MOEAD和MOEADFRRMAB；评估次数在3 000之后，d₂_NSGAⅡ的下降速度仅仅略逊色于MOEADFRRMAB，但是远远超过其他4种算法。图 7算法的IGD变化情况印证了d₂_NSGAⅡ算法在算法多样性与收敛性方面的效果。

表 3给出d₂_NSGAⅡ与NSGAⅡ、NSGAⅢ、VaEA、MOEAD、MOEADFRRMAB 6种算法在测试函数上的IGD值，以及他们的优劣对比，表现优秀的IGD值用加粗字体标注。将目标数分别设置为3, 5, 8, 10 (表 3)，在28组测试函数中，d₂_NSGAⅡ共取得7个最佳IGD值，NSGAⅡ共获得1个最佳IGD值，VaEA取得1个最佳值，NSGAⅢ共取得4个最佳IGD值，MOEAD共取得11个最佳IGD值，MOEADFRRMAB取得4个最佳IGD值。从表中的Wilcoxon秩和检验结果看，相比于MOEAD、MOEADFRRMAB、VaEA、NSGAⅢ、NSGAⅡ，d₂_NSGAⅡ净胜得分(数据优于对比算法的个数减去差于对比算法的个数)为-2, 0, 6, 1, 16。由此可见，除算法性能稍微劣于MOEAD、MOEADFRRMAB，d₂_NSGAⅡ均优于其他几种算法，故认为d₂_NSGAⅡ算法在求解DTLZ系列测试函数具有较为优秀的性能。

4 结论

为解决NSGAⅡ在高维空间中效果不佳的问题，研究提出基于d₂距离的NSGAⅡ算法，在保证算法收敛性的同时增加其多样性。通过实验数据验证，该算法在MaOP上取得了不错的效果，但是随着目标数目的增加，算法性能下降。究其原因，随着目标数目的增加，个体之间变得相互非支配，这是传统非支配方法面临的统一问题。未来可从支配关系以及距离度量方面进一步改进与探索，同时可以设计一些更加复杂的MaOP问题检验算法的有效性，使其效果更加优秀。