当前在解决目标数超过3个的高维多目标优化问题(Many-objective Optimization Problems,MOPs)中,由于目标函数数量的增加,传统的基于Pareto占优关系生成Pareto最优解时,互不占优解的数目急剧增多,使得Pareto最优解的选择压力降低,难以有效地逼近Pareto前沿^{[1, 2]}.一些经典算法:NSGA-Ⅱ^[3]、SPEA2^[4]和PESA-Ⅱ^[5]对于高维多目标优化问题都不能很好地适应^{[6, 7]}.Saxena和Deb^[8]分别运用PCA和最大变化伸展PCA以及相关熵PCA的方法来减少冗余目标,但是多个目标函数不存在冗余时是无效的.文献[9]提出了Pareto自适应ε占优机制,结合Pareto前沿面的几何形状和前沿面的不同部分的占优强度信息,自适应地调节ε矢量中元素的值.实现了高维多目标的降维,但是这种方法会改变部分解的支配关系,求解有些高维目标优化问题时解的质量较差.文献[10]基于GDDR(Group Divided Dimensional Reduction)算法实现了多目标降维处理,虽然性能上有一定的优势,但是该方法在目标函数进行分组时具有较大的随意性,不同的分组方法将导致不同的优化结果.文献[11]应用文献[12]提出的基于聚集函数的方法,使用加权的最小最大法以及结合对偶优化的权角距离的度量方法,能用于解决MOPs问题,但是处理较低维数时效果欠佳.文献[13]提出的efficiency of order方法能够将Pareto解集中的每个解的目标向量向K维空间进行投影.该方法能够找到更逼近Pareto前沿且分布均匀的非支配解,在解决低维的多目标优化问题时有一定的优势.文献[14]利用超体积、分布度和延展度将高维多目标优化问题转化为具有3个目标函数的多目标优化问题,虽然从理论上证明了可应用于求解高维多目标优化,但是设计高性能的集合进化策略往往是十分困难的,这给此种方法的应用造成了很大的限制.

本文从缩小搜索空间入手,在理论上证明了通过子目标函数值排序进行Pareto最优解求取的可行性,能够有效减少比较次数,加快求解速度.同时,基于SPEA2(Strength Pareto Evolutionary Algorithm2)算法框架,对密度评估策略中的距离和密度计算公式进行重新定义.提出的SOEA(Sub-Object Evolutionary Algorithm)算法同文献[15, 16, 17]中的高维多目标优化算法,在测试函数DTLZ^{[2, 18]}上关于收敛性指标^[17]、覆盖率指标^[19]进行比较,显示出一定的优势.

n个决策变量和m个目标变量的多目标优化问题^[2]可以表述为

定义1 Pareto占优^[3].若x_a,x_b∈X_f是多目标优化问题的两个可行解,则称x_a是比x_b Pareto占优的.

1) 若对于∀i=1,2,…,m;f_i(x_a)≤f_i(x_b)且∃j=1,2,…,m;f_j(x_a)＜f_j(x_b)，则记作x_asx_b,也称为x_a强支配x_b.

2) 若对于∀i=1,2,…,m;f_i(x_a)≤f_i(x_b)记作x_ax_b,也称为x_a弱支配x_b.文中的支配概念,均是该定义下的弱支配.

在证明相关命题之前,需要对R^m空间(m维欧氏空间)上的多目标优化问题,进行如下的假设:

1) 多目标优化问题的各个子目标函数在解空间上是连续的.若存在可去间断点和阶跃间断点,则利用目标函数在解空间连续的方法进行讨论,而不影响问题求解结果.

2) 解空间是Rⁿ上的有界闭子空间.

引理1 多目标优化问题的可行解集是半序集.

证明    设x₁,x₂∈X_f,对应的目标函数为F₁=(f₁(x₁),f₂(x₁),…,f_m(x₁)),F₂=(f₁(x₂),f₂(x₂),…,f_m(x₂)),F₁,F₂∈R^m,“≤”表示实数域条件下的小于等于含义.①对于x₁,有f_i(x₁)≤f_i(x₁),i=1,2,…,m成立,由解支配的定义即得x₁x₁.②对x₁,x₂,若有f_i(x₁)≤f_i(x₂),f_i(x₂)≤f_i(x₁),即f_i(x₂)=f_i(x₁)成立.由①的结论可得x₁x₂(x₂x₁).③对x₁,x₂,x₃∈X_f,若有f_i(x₁)≤f_i(x₂),f_i(x₂)≤f_i(x₃),即有f_i(x₁)≤f_i(x₃)(i=1,2,…,m)成立,则x₁x₃.所以,在可行解集上定义的Pareto占优关系,即Pareto支配“”,是集合X_f上的半序.    证毕

引理2 子目标的最优解集是全序集.

证明    设{x_p}是子目标f_i(x)可行解的集合,则{x_p}⊂X_f.由于X_f≠∅.对∀x₁,x₂∈{x_p},x₁≠x₂,有f_i(x₁)≤f_i(x₂)或f_i(x₂)≤f_i(x₁),则{x_p}是全序集.由于f_i(x)是解空间上的连续函数,解空间是Rⁿ上的有界闭子空间,故{x_p}存在上界.由Zorn引理^[20]可知X_f至少存在一个极大元.    证毕

引理3 极大元是非支配解.

证明    对于多目标优化问题的可行解集X_f,X_f已被证明为半序集.若x^*是X_f的极大元,则由极大元的定义^[20],对∀x∈X_f,若xx^*,则x^*=x.故对∀x∈X_f,有x^*x成立.也即:¬∃x∈X_f使得xx^*成立,则x^*是多目标优化问题的非支配解.    证毕

定理1 设{f₁,f₂,…,f_m}是优化问题的目标函数集,X_f为其解空间,x^**为其Pareto非支配集;x^*₁,x^*₂,…,x^*_m,是使得各目标函数取得极值的解,则集合Y={x^*₁,x^*₂,…,x^*_m}是非支配解集,且Y是x^**的子集.

证明    对∀x_i∈X_f,取集合Y中的任一元素x^*_j,由定义,x^*_j使得f_j(·)取得极值,则由支配的定义,对∀x_i∈X_f,使得f_j(x_i)≤f_j(x^*_j)不成立,故x_ix^*_j不成立,即x^*_j是极大元.由引理3知x^*_j是非支配解,故集合Y={x^*₁,x^*₂,…,x^*_m}是非支配解集.

往证,Y是x^**的子集.对∀x^*_j∈Y,由于Y是非支配解集,故x^*_j∈x^**.若有x^*′_j∈x^**,则对∀x_i∈X_f,∃l,l∈{1,2,…,m},使得f_l(x_i)≥f_l(x^*′_j)≥f_l(x^*_j)成立,则f_l(x^*′_j)不是f_l(·)上的极值,故x^*′_j∉Y,则Y⊂x^**.    证毕

定理2 多目标优化Pareto非支配解x^**的求取,可分解为f₁,f₂,…,f_m个子目标进行,且此种求取非支配解的方法是完备的.

证明    对∀x_i,x_j∈X_f,若对∀l=1,2,…,m,有fⁱ_l≤f^j_l成立,则x_ix_j.若fⁱ_l≤f^j_l成立,则记F(fⁱ_l,f^j_l)=1,否则,记R(fⁱ_l,f^j_l)=0.(fⁱ_l,f^j_l)=1⇔x_ix_j,可见关系可分解为m个子目标函数值比较结果的逻辑与关系.若(fⁱ_l,f^j_l)=0,则x_ix_j不成立.

对于x^*_i(i=1,2,…,m)是使得各目标函数取得极值的解,则集合Y={x^*₁,x^*₂,…,x^*_m}是非支配解集,且Y是x^**的子集,已由定理1得证.若x^*′_i是在目标函数f_i(·)上取值仅次于的极值点x^*_i的解.由于x^*_i是Pareto非支配解,x^*′_i在目标函数f_i(·)上的取值仅次于x^*_i,故只需比较j≠i上的f_j(x^*_i)与f_j(x^*′_i)的取值大小,若f_j(x^*′_i)≤f_j(x^*_i)成立,则x^*_ix^*′_i不成立,由于x^*_i是非支配解,故x^*′_i也是非支配解.同理,若x^*″_i是在目标函数f_i(·)上的取值次于x^*′_i的解.若有f_k(x^*″_i)≤f_k(x^*′_i)、 f_j(x^*″_i)≤f_j(x^*_i),j≠i,k≠j成立,则x^*′_ix^*″_i,x^*_ix^*″_i不成立,则x^*″_i也是非支配解.继续此过程求取更多Pareto非支配解x^*′″_i,…,x^*(n′)_i,记这些非支配解构成的集合为{x^*(n′)_i}.

往证,完备性.对∀x∈x^**,若x∈Y,则∃i∈{1,2,…,m}使得f_i(x)=minf_i(·),则通过子目标排序的方法一定能找到minf_i(·),即能找到非支配解x.若x∈x^**,x∉Y,不妨设在某子目标i∈{1,2,…,m}上有f_i(x)≥f_i(x^*″)≥f_i(x^*′),由于f_i(x^*′)=minf_i(·),故x^*′必定被找到,且x^*′x^*″不成立,由于x^*″在f_i(·)上的取值仅次于x^*′,故x^*″必定被找到,依此法x必定被找到,则对∀x∈x^**都可被找到.    证毕

SOEA算法是以进化种群中子目标函数值排序为基础.根据求取的是目标函数的最大(最小)值,对子目标函数值进行排序.通过子目标值的排序,快速找出非支配解.以子目标函数f₁为例,要求取目标函数的最小值.进化过程如图 1所示,通过排序,得到使f₁取最小值的解空间中的点A,A为通过目标函数f₁获得的第1个非支配解.继续该过程得到C、C₁在f₁上的函数值仅比A大,但是C、C₁在子目标f₂上的函数值比A小,故C、C₁相对于A为非支配解.但是C支配C₁,故C为非支配解.由此可见,子目标进化算法的方向性很强.

图 1 进化过程示意图Fig. 1 Schematic of evolution process

图选项

SOEA算法利用了SPEA2算法的框架,与SPEA2算法不同之处有:利用子目标函数值排序,快速获得非支配解集,而不需要通过进化种群中的解之间两两比较来获得非支配解;利用归一化函数差值的Minkowski距离进行解之间的距离计算,通过“k近邻”距离来评判密度值的大小.而不是通过计算进化种群中解之间的欧氏距离.利用从外部种群中获取下一代进化的解,而不是在外部集合和进化种群中选取下一代的进化种群.

输入值:进化种群的大小M,外部集合的势,最大进化代数P.

输出值:非支配解集U_non.

第1步 初始化,产生初始种群U₀,并且产生一个空的外部集合=∅,设置进化代数p=0.

第2步 利用子目标优化方法求取进化种群U_p中的Pareto非支配解.

第3步 将产生的非支配解加入外部集合,如果非支配解的个数大于外部集合的大小,则利用截断方法进行操作;如果产生的解的个数小于外部集合的大小,则利用适应度的值,选择满足条件的解进入外部集合.

第4步 如果达到最大进化代数P,则停止计算.将外部集合中的非支配解输出到U_non中作为最终的非支配解集.

第5步 计算外部集合中解的适应度,对外部集合利用锦标赛的方法进行选择,当被选择的解的个数等于进化种群的大小M时停止选择.

第6步 利用交叉、变异操作对进化种群中的解进行操作,产生新一代种群U_p+1,转入第2步.

在求取外部集合的解时,设置一个临时集合用于存放由进化种群中选出的解,这些解在满足一定的条件时将被选入外部集合,临时集合作为求取过程中的辅助存储空间.根据外部集合的势=||,对于有m个子目标的多目标优化问题,首先由每个子目标对外部集合提供M_i=/m个元素.利用子目标快速选择方法生成非支配解,生成的非支配解的数目为|U_pi|,当|U_pi|=/m时,暂时中止在此子目标上生成Pareto非支配解.将以上生成的非支配解加入临时集合,同时将这些非支配解从进化种群中删去.依据子目标f₁,f₂,…,f_m依次生成非支配解集U_p1,U_p2,…,U_pm.按最后一个子目标生成非支配解时,如果产生的非支配解的个数已经等于外部集合的大小.则依次选择进化种群中的解,按照该解中子目标的值最小(或者最大)的子目标,与临时集合中相应的非支配解进行比较,直至进化种群中所有的解比较完毕.此时,进化种群中的全部非支配解被选入临时集合.在选择非支配解进入临时集合的过程中,如果进化种群中的解选择完毕,也没有达到外部集合的个数,则利用适应度值选择支配解进入外部集合.

解的选择是以适应度为准则进行评判的,解的适应度计算由解的强度和解的密度两部分组成.强度由解之间的支配关系确定;密度评估以解的距离计算为前提.密度值小,说明该解与空间中其他解的距离大,从而优先考虑被保留.关于距离的计算也有很多不同的计算方法,NSGAⅡ(Nondominated Sorting Genetic Algorithm Ⅱ)在子目标函数值排序的基础上,通过相邻的两个函数值的差值求取绝对值,归一化处理后,对每个子目标上的这些差值进行求和,从而得出个体的排挤距离.SPEA2算法中,通过计算两个非支配解之间的欧氏距离,然后对与本非支配解相关的距离进行排序,利用“k近邻”的方法,求取在排序中的第k个距离值的大小,通过第k个距离值的大小计算密度值.Pareto存档进化策略(Pareto Archived Evolution Strategy,PAES)算法利用空间超格机制来保持种群的多样性.超格的数目需要预先设定,但是当目标的数目增多时,超格的数目难以确定,并没有很好的方法可循.以上的密度评估策略应用在高维多目标优化问题上时,获得的效果并不理想.

根据以上方法的不足,提出了更为广义的基于归一化Minkowski距离的“k近邻”方法.流程如下.

第1步 对于临时集合中的解,找出每个子目标上的最大值和最小值,计算两个解的同一子目标上函数值的差值,对这一差值进行归一化处理.

第2步 利用Minkowski距离计算公式,对两个解的子目标进行归一化处理后的差值进行距离计算.

第3步 对于临时集合中的每个解,对计算所得的与其他解之间的Minkowski距离按升序进行排列.

第4步 k值取外部集合大小的方根,k=.记与解X_i距离为k的解的距离为σ^k_i.密度评估按D(i)=计算.

子目标函数差值的归一化是对每个子目标函数f_i,在临时集合中分别求取最大值f^max_i和最小值f^min_i,则两个解X_k、X_l在此子目标函数上的差值归一化后的形式为

Minkowski距离计算公式为

当p值较小时,对差值的影响不明显,相当于均衡考虑各个目标;当p值较大时,差值大的目标对距离的影响明显.可以看出,NSGAⅡ和SPEA2中的距离计算公式是式(2)在p取值为1和2时的特例.在高维多目标优化时,非支配解的数目增长很快,均衡的进化方法势必导致收敛到Pareto前沿的速度过慢,因此适当地增大p的取值可以加速收敛速度.本文中采用随机的p值选取方法,在[1,H]的范围内随机选取一个p值,H为某正整数,作为每一代进化的Minkowski距离计算公式中的幂次,在同一代进化中p取相同的值.

1) 非支配解个数小于外部集合时解的选取.

利用子目标快速选择方法中产生的非支配解的数目不足时,通过计算进化种群中的支配解的适应度,选择适应度值小的支配解进入外部集合,直至满足外部集合的大小要求.适应度计算公式为

2) 非支配解个数大于外部集合时解的选取.

由于此时临时集合中所有的解都是非支配解,适应度值计算公式F(i)=R(i)+D(i),退化为F(i)=D(i).将密度评估值进行升序排列,取密度值小的前个解进入外部集合.非支配解x_k1、x_k2到x_i的距离相同时,即σ^k1_i=σ^k2_i.统计非支配解x_k1、x_k2在各个子目标函数上的差值的正负.若在f_j(x_k1)－f_j(x_k2)(j=1,2,…,m)中,取正值的子目标函数占多数,则当优化是求取最小化各目标时,选取x_k2进入外部集合;反之,选取x_k1.

对于含有m个子目标,进化种群中有M个解,外部集合的势为的多目标进化问题.初始化的时间复杂度是O(M),求取非支配解的时间复杂度是O(KM);求取外部集合的时间复杂度是O(L²)+O(LlogL),其中,L是临时存储集合0＜L≤M;交叉、变异及种群更新时间复杂度是O(M).总的时间复杂度为:O(M)+O(KM)+O(L²)+O(LlogL)+O(M),按照时间复杂度的运算规则,最差的时间复杂度应为:max{O(KM),O(L²)}.

对于SPEA2算法,初始化的时间复杂度是O(M).非支配解求取由适应度计算决定,适应度计算包含强度和距离计算.计算强度R(i)的时间复杂度为O((M+)²);计算距离D(i)的时间复杂度为O((M+)²log(M+)).求取外部集合的时间复杂度为O((M+)²log(M+));交叉、变异以及种群更新的时间复杂度为O(M).所以,SPEA2算法的平均时间复杂度为:O((M+)²log(M+)).由时间复杂度的分析可以看出,本文提出的算法在时间复杂度方面与SPEA2算法相比具有明显的优势.

解决高维多目标优化问题的基于子集覆盖实现目标减少的近似k-MOSS算法,时间复杂度为O(m³L²),其中m为子目标的个数,L为临时集合的大小.基于特征选择的目标减少算法,通过选择最拥挤的邻居群,保留该邻居群的中心,然后删除其余的目标,算法的时间复杂度为O(m²L).

1) 测试函数.

测试函数使用高维多目标进化领域被广泛使用的DTLZ1~DTLZ7函数,函数的具体形式可参见文献[2, 18],在此不再赘述.

2) 收敛性指标^[17].

令G^*={g₁,g₂,…,g|_G^*|}是理想的均匀分布的Pareto最优解集，Q={a₁,a₂,…,a_Q}是通过进化算法得到的Pareto最优解集.则对于集合Q中的解a_i,通过式(3)得到该解距离G^*的最小欧氏距离:

IGD评价方法用于衡量所求的Pareto前沿面到真实的Pareto前沿面的收敛程度,该指标的值越小,表明算法得到的解的收敛性越好,越接近理想Pareto前端,定义为

3) 覆盖率指标C^[19].

E、F是优化问题的两个解集,覆盖率指标C定义为

设置进化种群的规模为600,外部种群的大小为200,进化代数为400代,每个算法独立运行20次,评价指标为20次运行的平均值.为了测试SOEA算法的性能,选择类似的解决高维多目标优化问题的算法进行比较,分别是由Jaimes和Coello提出的基于特征选择的算法(记为JC),以及由Brockhoff和Zitzler提出的基于最小子集覆盖的目标减少算法k-MOSS(记为BZ).利用以上3种解决高维多目标优化问题的算法对DTLZ1~DTLZ7测试函数进行计算,比对各算法对不同的测试函数的计算性能.由于以上3种算法对7个测试函数进行求解后,求解结果表现出一定的相似性,选取SOEA算法、BZ算法、JC算法对DTLZ2、DTLZ5、DTLZ7测试函数的测试结果在文中列出.

SOEA算法在进化的过程中,并没有减少目标函数的个数,所以在进化过程中解之间的支配关系不会发生变化.求得的非支配解必定是原来的多目标优化问题的Pareto最优解,从表 1中的数据也能说明这一点.可以看出SOEA算法的IGD值比BZ算法和JC算法得到的值要小,说明SOEA算法获得的Pareto前沿与最优Pareto前沿接近程度更好.JC算法由于删除了过多的目标函数,解的支配情况与BZ算法相比改变程度会大一些,得到的IGD值也比BZ算法得到的要大.

由表 2可见,当目标数为12维时,通过对获得的非支配解进行统计,得出利用SOEA算法求得的非支配解,可以支配BZ算法求得的非支配解,占BZ算法求得解的0.544 3,而BZ算法中可支配SOEA算法中的解,只占到SOEA算法求得解的0.396 2.可以得出利用SOEA算法所求得的非支配解集比BZ算法的质量要高，与JC算法比较结果的统计数据也显示出这一结论.

通过图 2可见,JC算法的运行时间比另两个算法占优,这是因为该算法是直接在目标函数上进行处理的,通过删除冗余的目标函数,来达到降低维数的目的,算法的时间复杂度相当于处理低维的目标优化问题.SOEA算法和BZ算法均需要进行支配关系的比较,时间复杂度与外部集合的大小相关.可以看出SOEA算法相比较于BZ算法,基本上是正比例相关的.SOEA算法在时间上要优于BZ算法.

图 2 求解DTLZ2~DTLZ7时各算法运行时间Fig. 2 DTLZ2-DTLZ7 computation time of every algorithm

图选项

通过对3个算法在测试函数DTLZ5、DTLZ7上相应值的统计结果,显示了与在DTLZ2上类似的结论,如表 3～表 6和图 2所示.

在分析高维多目标优化问题特点和研究现有多目标优化求解算法的基础上,得出如下结论:

1) 高维多目标优化问题Pareto非支配解x^*的求取,可以通过分解为f₁,f₂,…,f_m个子目标,并在子目标函数值排序的基础上获得,且此法获得的Pareto非支配解完备.

2) 解的适应度值求取时,在基于归一化函数差值的Minkowski距离计算方法中采用随机在[1,H]范围内选取p值的方法可以加快算法的收敛速度.

3) 将本文提出的SOEA算法与BZ、JC算法比较,得出SOEA算法在收敛性和覆盖率两个指标上均优于BZ、JC算法.运行的时效性JC算法最优,SOEA算法次之,BZ算法最差.