乔俊飞 李霏 杨翠丽

大多数的工程和科学应用,如工业制造、城市运输、污水处理和资本运算等,几乎每个重要的现实生活中的决策都存在多目标优化问题.这些目标往往是不可比较,甚至是相互冲突的.因此多目标优化问题(Multiobjective optimization problem,MOP)一直是近几年主要的研究课题之一.为了解决该问题,多目标进化算法(Multiobjective evolutionary algorithm,MOEAs)已被广泛研究.其中Schaffer[1]提出了向量评估遗传算法(Vector evaluated genetic algorithm,VE-GA),该算法对单目标遗传算法进行了改进.但VEGA只能找到Pareto前端的始末端,产生不了均匀分布的解.为此Fonseca等[2]提出了多目标遗传算法(Multiobjective genetic algorithm,MOGA),这种算法效率较高,而且易于实现.但该算法的不足在于如果小生境数目信息是基于目标函数的,那么两个具有相同目标函数向量的不同个体无法在同一代种群中存在.为了解决该问题,Srinivas等[3]提出了非劣分类遗传算法(Nondominated sorting genetic algorithm,NSGA),该算法先按解个体的非劣关系进行排序,再按照共享机制来保持进化的多样性.但是该算法计算效率较低,计算复杂度大,且共享参数σ需要预先确定.为了减小算法的计算复杂度,Deb等[4]提出了改进型非劣分类遗传算法(Nondominated sorting genetic algorithm II,NSGAII),该算法采用非支配排序进行分级,通过计算拥挤距离选择最优解,并将其做为精英解保存起来.但该算法是一种类随机搜索算法,存在操作次数多收敛速度慢和解分布特性较差的问题[5].

由于种群在迭代过程中种群分布不均会导致种群多样性变差从而易使解陷入局部最优,并且在迭代过程中Pareto前沿部分区域易出现空白,收敛速度变慢.为了解决种群分布不均带来的上述的问题,Goldberg等[6]提出了基于共享机制的小生境技术,该方法考虑了一个个体与群体中其他个体之间的相似程度,但计算开销比较大.朱学军等[7]提出用进化群体熵来刻画进化群体的多样性和分布性.然而这种方法缺乏对群体内部个体之间关系的刻画,因此不便于调控群体进化过程中的多样性和分布性.为了解决这个问题,Corne等[8]提出了网格方法,将网格中聚集密度大的个体删除,但不能删除极点个体.Knowles等[9]提出了自适应网格技术,在每一代进化时,根据当前代的个体分布情况自适应地调整边界.Morse[10]提出了聚类分析的方法来保持种群的多样性.Han等[11]提出了基于种群间距信息和种群分布熵的方法.但是以上方法只能删除种群中拥挤距离较小的个体,无法解决解的过程中Pareto前沿部分区域出现空白的问题.

基于以上问题,本文提出了一种基于均匀分布的NSGAII(NSGAII-UID)多目标优化算法.该算法受文献[12]的启发,将种群映射到目标值对应的超平面上,并在该平面上聚类.但是文献[12]中提出的算法仍然存在种群分布不均的问题,从而影响了种群的多样性,为了解决该问题,本文将映射平面均匀分区.当对应区间的分布性不满足时,分布性加强模块激活.由于种群在迭代的过程中对应区间会出现种群个体不足或缺失的状况,此时需要在该区间内放入一些个体.为了解决该问题,本文将所选聚类子群体中拥挤距离最大的点进行局部搜索,采用极限优化变异[13]的方法产生新的个体.实验结果表明,该方法综合评价指标(Inverted generational distance,IGD)值和分布性评价指标(Spacing,SP)值均高于其他算法.因此表明该方法具有较好的种群多样性和分布性,且收敛速度较快.

1 基本概念

1.1 多目标及相关定义

定义1.多目标优化问题(MOP)

一个具有n个决策变量,m个目标函数的多目标优化问题可以描述为:

式中,x称为决策变量,X是n维的决策空间;F(X)∈Rm为m维目标向量;fj(X)(j=1,2,···,m)为第j个目标函数,li和ui分别为第i个决策变量的上界和下界.

定义 2.Pareto-占优对于给定的两点x,x∗∈Xf,x∗是Pareto-占优(非支配)的,当且仅当式(2)成立,记为x∗>x.

定义 3.Pareto-最优解

若对于任意解x,不存在x’∈Ω使得F(x’)=(f1(x’),f2(x’),···,fm(x’))占优于F(x)=(f1(x),f2(x),···,fm(x)),则称x为Pareto最优解或者非劣解.

定义4.Pareto-最优解集

所有Pareto-最优解组成的集合Ps称为Pareto-最优解集.

定义 5.Pareto-前沿

Pareto-最优解集合Ps中的解对应的目标函数值组成的集合PF称为Pareto-前沿,即:

1.2 模糊 C-均值聚类算法

设n个数据样本为P={p1,p2,···,pn},c(2≤c≤n)是要将数据样本分成不同类型的数目,{A1,A2,···,Ac}表示相应的C个类别,U是其相似分类矩阵,各类别的聚类中心为{v1,v2,···,vn},µk(pi)是样本pi对于类Ak的隶属度(简写为µk).则目标函数Jb可以表达为:

其中,dik是欧几里得距离,用来度量第i个样本pi与第k类中心点之间的距离;m是样本的特征数;b是加权参数,取值范围是1≤b≤∞.模糊C-均值聚类方法就是寻找一种最佳的分类,以使该分类能产生最小的函数值Jb.它要求一个样本对于各个聚类的隶属度值和为1,即满足:

设Ik={i|2≤c

2 基于均匀分布的NSGAII算法

为了解决种群分布不均的问题,本文提出了将目标空间均匀分区,选取相同数量个体的方法.由于目标空间一般为曲线或者曲面,想要均匀划分较难操作.为了解决这个问题,本文受文献[12]的启发,在文献[12]中,种群中的个体映射在一个超平面上,并在目标空间中聚类.为了验证种群的多样性将集群质量度量标准.但是该方法在迭代过程中存在种群分布不均的问题,从而影响种群的多样性,进而造成收敛速度慢,Pareto前沿有些区域出现空白的现象.针对以上问题本文提出了将种群映射到目标函数值对应的超平面H,在H上进行均匀分区以增加种群的多样性的方法.该方法将种群在映射平面上进行聚类,并记录聚类中心,计算每个分区中的聚类中心个数,从而判断是否满足均匀分布条件.当条件不满足时,分布性加强模块激活.由于在计算过程中区间内有时会出现个体数量不足或者为空的情况,为了在该区间获取缺失的个体,本文采用极限优化变异产生新的个体.当该区间个体不为空时,将聚类子群体按照拥挤距离从大到小排序,选取第一个个体进行变异.当该区域为空时,将离该区间最近的个体进行变异.该方法不但可以使种群跳出局部最优,并且能够较快地找到Pareto前端.

2.1 映射

设目标向量为:

定义超平面[12]H为k维欧氏空间的k−1维线性子空间,bs=−1,h·i为欧氏内积

目标向量f映射到超平面[12]H为PI,公式如下所示:

2.2 分布性

由于NSGAII在运算过程中种群的分布性较差,为了解决这个问题,本文加入了分布性判断模块.该模块用来判断每个均匀分段区间内聚类子群体的分布情况,同时选取聚类子群体中种类数量最大值作为阈值来判断分布性.当区间内种群类别小于该阈值时,则判断该区间不满足分布性,分布性加强模块激活.否则,该区域满足分布性,且选取每个聚类子群体中拥挤距离最大的个体放入精英档案.

该算法首先将目标区间均匀分区,其中每个区间的种群分布如下:

n为分段区间的数量,Doi为第i个区间所包含的子群体所有类别的集合,i∈[1,n],如下所示:

NRPij为第i个区间内的种群中心编号,其中i∈[1,n],j∈[1,r],r为种群中心数量,ABnrk为第i个区间第r个聚类中心所对应的个体,k∈[1,m],m为该类子群体个体数量.区间内种群中心数量最大值为maxnr.判断种群分布性的方法如下所示:

其中,ir为第i个区间内的种群中心数量,maxDon为分段区间内种群中心数量最大值.

情况1.若该条件成立,则分布性加强模块激活.

情况2.将该类子群体按照拥挤距离从大到小排序,选择每个聚类子群体第一个个体.

2.3 分布性加强

当种群不满足上述分布性条件时,分布性加强模块被激活.在种群选取的过程中,每个区间选择相同数量的个体.同时为了增加种群的多样性,从每个聚类子群体中按照方法1)与2)选取拥挤距离较大的个体放入精英档案.然而在运算的过程中,区间内有时会出现种群数量不足或空缺的状况.在以往文献中多是采用将拥挤距离较小的个体删除的方法,当种群空缺或不足时,并未有相关文献进行介绍.为了解决这个问题,本文提出了对所选精英个体进行极限优化变异的方法,来填补欠缺的个体.由于在计算过程中区间内种群中心的数量不同,采取的策略也不相同.具体方法如下所示:

1)第i个区间种群中心个数量ir为0:

其中icm为在第i个区间内第c个聚类中心对应的子群体中的个体数量,c∈[1,ir]∪0.

情况1.将该类子群体按照拥挤距离从大到小排序,取前maxnr个个体.

情况2.选取所有个体,剩余maxnr-icm个个体由拥挤距离最大的个体变异得到.(详见第2.4节)

情况3.由于该区间无个体,本文选取离该区间最近的两个个体变异得到maxnr个个体(详见第2.4节).若该区间连续20次为空则该区间设定为空.

2)第i个区间种群中心个数量ir不为0:

当该区间种群中心个数不为零时,需要从聚类中心对应的子群体中选择相应的个体.选择方式如下所示:

S-num为每个聚类子群体应当平均选择的个体数目,Renr为需要选择的剩余个体的集合,|Renr|为需要选择的剩余个体的集合的数量,Renri为剩余的个体对应编号,i∈[1,t],t为余数的个数.根据余数的是否为0,选择方式也不相同,具体如下:

a)余数Rer为0:

情况2.选取所有个体,剩nrm-S-num个个体由变异得到.(详见第2.4节)

b)余数Rer不为0:

当剩余个体数量不为0时,选择的方式如下所示:

2.4 局部变异策略

由于在迭代过程中有的时候区间会出现种群数量不足或为空的状况.为了保证每个区间个体选取数量相同.本文采取最优个体极限优化变异策略.首先找到需要选择子群体中拥挤距离最大的个体或者离该区间最近的两个个体,然后对其进行变异,本文采用两种变异策略来产生局部解[13]:

1)第一种变异策略:

对于选定的对象,每次只有一个决策变

量发生变异,这样有效的提高了种群的局部搜索能力,从而提高了解的计算精度.设当前选定个体ABnrm=(ABnrm1,···,ABnrmi,···,ABnrmj),j为决策变量的个数.种群总数为N,则产生局部解的个数为µ,µ的值由上文(详见第2.3节)的需要动态决定.具体如下:

其中,ABnrmi为待变异个体的决策变量;α为 (0,1)区间的随机变量;β∈R,且为形状参数,本文通过多次实验得到β设为11为最佳.θmax(ABnrmi)为决策变量ABnrmi可变动区间的最大值.该变异方法因为每次只将一个决策变量变异,所以具有较强的局部微调能力.但是该方法只能够在小范围内搜索.为了避免种群陷入局部最优,同时提高搜索速度,本文加入了第二种变异策略.该策略将产生20%×µ个局部解[14](如非整数则取整),具体方法如下所示:

2)第二种变异策略:

其中,λ∈(0,1.2)为一随机数.以上两种变异策略共产生λ+[0.2λ]个局部解.

2.5 算法流程

根据具体问题初始化参数:设定目标函数有m个,初始种群数量为N,精英解数量为N/2,函数最大调用次数设为Imax,决策变量维数为j,决策变量的取值上下界为u=(u1,u2,···,ui,···,uj),下界为ll=(l1,l2,···,li,···,lj),形状参数设为g,交叉参数θc,交叉概率Pc,变异参数θm,变异概率Pm,具体的算法流程如下所示:

步骤 1.初始化种群P={p1,p2,···,pm,···,pN},其中pm=(x1,x2,···,xj),m=1,2,···,N,其中xj∈(lm,um).计算种群中每个个体对应的目标函数值.

步骤2.对种群P中的解进行非支配排序,排序后当前种群的所有非支配解记为Pc.

步骤 6.当区间内的子群体不满足均匀分布条件时,需要按照第2.3节所示增加种群的分布性,若该区间个体数量不足时,缺少的数量通过对拥挤距离最大的个体变异得到,详见第2.4节.当前选择出来的用来增加分布性的种群记为Pe,种群数量为n×maxnr.

步骤7.将种群Pc中包含在Pe中的个体去掉,然后按照拥挤距离选取N/2−n×maxnr个个体,记为Ps.

步骤 8.合并Pe与Ps,得到种群Pm.将Pm进行交叉变异,得到新的种群Ph.合并Pm和Ph为下一代种群PT.

步骤9.重复步骤2至步骤8,当达到最大迭代次数或者预设的目标时,进行下一步.

步骤10.当前PT中的非支配解即为得到的最优解.

3 仿真实验

本文采用Matlab R2013b版本,处理器为3.60GHz,8.00GB内存,Microsoft实验环境.为了验证算法的收敛性和多样性采用了以下性能评价指标对算法进行了验证.

1)综合评价指标(Inverted generational distance,IGD):

IGD用来评价算法的性能,它的值越小则解的收敛性和多样性就越好.具体计算公式如下所示:

2)分布性评价指标(Spacing,SP):

SP用来测量已知帕累托前沿相邻解间距离的方差.定义如下所示:

q为非支配解的个数,q=2,3,···,n,为ui的平均值.

3.1 实验1

由于双目标ZDT系列函数和三目标DTLZ系列函数被广泛地用于验证多目标优化算法[13],本实验分别采用(ZDT1∼ZDT4,ZDT6)函数和(DTLZ1、DTLZ2)函数来验证算法的性能.该系列函数的特征、维度和种群规模如表1所示.

表1 测试函数参数Table 1 Paramter setting of the test function

为了检验NSGAII-UID算法的性能,该算法分别和一种基于密度的局部搜索NSGAII算法 (NSGAII-DLS[14])、标准 NSGAII[15]、定向搜索算法 NSGAII-els[16]、随机局部搜索算法 HMOEA/D[17]、自适应多目标粒子群算法(AMOPSO[11])、多目标差分进化算法(MODE[16]6种算法进行了对比.

本实验采用模拟二进制交叉和多项式变异方法,交叉参数ηc和变异参数ηm均设为20,交叉和变异概率分别为0.9和1/n,n为决策变量的数量;形状参数q设为11.在进行ZDT实验时,不同的多目标优化算法采用相同的种群数量:最大计算量设为Imax=5000,种群规模设为100,由于通过交叉变异会产生50个子代,所以共计迭代了100次.为了验证算法的有效性,实验结果如图1∼图5所示.

在进行DTLZ系列实验时,参数选择为:最大计算量设为Imax=15000,种群规模设为100,由于通过交叉变异会产生50个子代,所以共计迭代了300次.为了验证算法的有效性,实验结果如图6∼图7所示.

图1 ZDT1优化效果Fig.1 The optimization effect of ZDT1

图2 ZDT2优化效果Fig.2 The optimization effect of ZDT2

由图1∼图7可视,对于凹、凸和非连续的多目标函数,NSGAII-UID可以较好地逼近pareto前端且分布较均匀,NSGAII-UID与其他算法对比结果如表2所示.该表分别列出了当函数调用次数分别为 5000次和 15000次时,实验结果连续30次的最大值、最小值和平均值.从该表可以看出NSGAII-UID在5个测试函数(ZDT1、ZDT2、ZDT4、DTLZ2、DTLZ7)中 IGD的最大值、最小值和平均值均优于其他算法.在ZDT3和ZDT6中由于Pareto前沿在指定区间内存在空白区域,NSGAII-UID的IGD值略小于AMOPSO[11]算法,在后续的工作中将会就此问题进行更深一步的研究.因此从以上结果可以看出和对比算法相比,该算法具有较好的精度和收敛性.

表2 NSGAII-UID与其他局部搜索算法的ZDT和DTLZ系列实验IGD结果Table 2 ZDT and DTLZ series performance IGD comparison of different algorithms

图3 ZDT3优化效果图Fig.3 The optimization effect of ZDT3

图4 ZDT4优化效果Fig.4 The optimization effect of ZDT4

表3 NSGAII-UID与其他局部搜索算法的ZDT和DTLZ系列实验SP结果Table 3 ZDT series performance SP comparison of different algorithms

由图1∼图7可视,对于凹、凸和非连续的多目标函数,NSGAII-UID可以较好地逼近Pareto前端且分布较均匀,NSGAII-UID与其他算法对比结果如表2所示.该表分别列出了当函数调用次数分别为5000次和15000次时,实验结果连续30次的最大值、最小值和平均值.从该表可以看出NSGAII-UID在5个测试函数(ZDT1、ZDT2、ZDT4、DTLZ2、DTLZ7)中IGD的最大值、最小值和平均值均优于其他算法.在ZDT3和ZDT6中由于Pareto前沿在指定区间内存在空白区域,NSGAII-UID的IGD值略小于AMOPSO[11]算法,在后续的工作中将会就此问题进行更深一步的研究.因此从以上结果可以看出和对比算法相比,该算法具有较好的精度和收敛性.

图5 ZDT6优化效果Fig.5 The optimization effect of ZDT6

图6 DTLZ2优化效果图Fig.6 The optimization effect of DTLZ2

图7 DTLZ7优化效果Fig.7 The optimization effect of DTLZ7

3.2 实验2

参照文献[11]的实验,对算法NSGAII-UID的分布性进行检验,不同多目标优化算法采用相同的参数,具体参数设置同实验1所示.该实验共与4个算法 (AMOPSO[11]、cdMOPSO[19]、NSGAII[8]、MODE[17])进行了比较,取连续30次实验平均值作为结果如表3所示.

由表3可见NSGAII-UID算法在6个测试函数(ZDT1、ZDT2、ZDT4、ZDT6、DTLZ2、DTLZ7)的SP最大值、最小值和平均值均小于其他对比算法.在ZDT3中SP值略小于AMOPSO[11]算法,由以上结果显示该算法和其他算法相比算法相比具有更好的分布性.

3.3 实验3

为了验证算法的收敛速度,本文采用记录算法达到指定性能指标时所调用的函数次数的方法.参照文献[9]的实验,设定当IGD值达到0.01时停止优化,记录函数调用次数.实验结果如表4所示.

表4 NSGAII-UID与其他局部搜索算法的函数计算次数结果Table 4 Function calculation comparison of different algorithms

由表4可见:在ZDT1-ZDT4和ZDT6的实验中,当IGD达到设定值时,连续实验10次的平均值,NSGA2-UID所用函数调用次数均小于其他方法,因此NSGA2-UID达到指定标准消耗的计算量更少,即该算法收敛到Pareto前沿速度更快.

4 结论

针对NSGAII在种群进化过程中会出现解分布不均的问题,本文提出了一种基于均匀分布的NSGAII(NSGAII-UID)多目标优化算法.为了使解能够在目标空间均匀分布,而真正的Pareto前沿大都是曲线或者是曲面,本文采用映射的方法,将解映射到目标空间对应的超平面,并在该平面均匀分区.当对应分区的解不满足均匀性时,均匀性加强模块被启用.同时采用聚类分析的方法来维持和增加进化种群的多样性和分布性.为了验证算法的性能,本文采用5个ZDT函数和示,两个DTLZ函数来进行实验,实验结果显示该算法和其他多目标优化算法相比具有更好的收敛性和分布性,同时收敛速度较快.