文章目录

• • • 一.多目标优化简介
    • • 1.多目标优化问题
      • 2.多目标优化求解思路
    • 二.NSGA-II算法解析
    • • 1.快速非支配排序(Fast non-dominated sort)
      • 2.拥挤距离计算(Crowding distance assignment)
      • 3.精英保存策略(Elitism)
    • 三.NSGA-II算法实现
    • • 1.测试函数

一.多目标优化简介

1.多目标优化问题

在很多实际工程问题中，我们的优化目标不止一个，而是对多个目标函数求一个综合最优解。例如在物流配送问题中，不仅要求配送路径最短，还可能需要参与运输车辆最少等。

多目标优化问题通常有如下特点：

• 各目标函数往往相互冲突，无法同时取到最优解。这种冲突使得解的搜索变得更加复杂；
• 各目标函数的单位不同，不能简单比较。因此多目标优化问题的合理解集通常是Pareto最优解。

2.多目标优化求解思路

对于多目标优化问题，传统方法是将原问题通过加权方式变换为单目标优化问题，进而求得最优解。该方法有两大问题：

• 权重的量化设定困难；
• 多目标优化问题的求解结果是包含一组非支配解的解集，用这种方法每次只能求得一个不同解(运气好的话),要求解具有足够多样性的Pareto最优解集非常困难

遗传算法具有多点多方向搜索的特征，在一次搜索中可以得到多个Pareto最优解，因此更适合求解多目标优化问题。

而当前用于求解多目标优化问题的遗传算法一般有两种思路:

• Pareto排序评价方法：其思想是构造基于解的优劣关系排序的适应值函数；在子代中更多保留支配解，从而在迭代一定次数后获得Pareto前沿。最经典的是Glodberg提出的Pareto ranking genetic algorithm，另外我们这里介绍的NSGA-II也是基于Pareto排序评价的方法。
• 多目标函数加权方法：给多目标函数加以各种权重，转化为单目标优化问题；相对于传统方法的固定权重，这类方法通常会在遗传操作前以一定规律生成权重，指导个体的多方向搜索。代表性的算法有Ishibuchi的random weight GA和Gen-Cheng的adaptive weight GA等。

二.NSGA-II算法解析

NSGA-II(nondominated sorting genetic algorithm II)是2002年Deb教授提出的NSGA的改进型，这个算法主要解决了第一版NSGA的三个痛点：

• 非支配排序的高计算复杂度
• 共享参数难以确定
• 缺少保存精英策略

针对这三个问题，在NSGA-II中，Deb提出了快速非支配排序算子，引入了保存精英策略，并用“拥挤距离”(crowding distance)替代了共享(sharing)。

在介绍NSGA-II的整体流程前，我们需要先了解快速非支配排序与拥挤距离的定义。

1.快速非支配排序(Fast non-dominated sort)

解的支配关系与Pareto最优解




快速非支配排序步骤：
快速非支配排序就是将解集分解为不同次序的Pareto前沿的过程。
它可以描述为：

1. 为每个解p分配两个关键量：支配p的解个数n_p以及被p支配的解集S_p；
2. 设置i=1，将n_p=0的个体归入F_i；
3. 对于F_i中的个体，遍历每个解p的S_p，将其中每个解的n_p减1；
4. i+=1，将n_p=0的解归入F_i;
5. 重复3、4，直到解集中所有个体都被归入某一个F_i。
   
   DEAP的实现
   DEAP内置了实现快速非支配排序操作的函数

tools.emo.sortNondominated(individuals, k, first_front_only=False)
参数：
individual:被排序的个体列表
k：需要选择的个体数目，注意这里不一定返回正好k个个体，需要看pareto前沿的个体个数
	假设设置k=100，而selectedPop + len(Front[i-1]) <100,返回的个体数是
	selectedPop+len(Front[i-1]) + len(Front[i]),应该是一个大于100的值。
first_front_only:如果设置为True，则对次序为1的前沿排序之后就返回

返回：
	Pareto前沿的列表

2.拥挤距离计算(Crowding distance assignment)

拥挤距离的定义
在NSGA-II中，为了衡量在同一个前沿中各个解质量的优劣，作者为每个解分配了一个拥挤距离，其背后的思想是让求得的Pareto最优解在objective space中尽量分散，也就有更大可能让解在Pareto最优前沿上均匀分布。

DEAP实现
DEAP中内置了计算拥挤距离的函数

tools.emo.assignCrowdingDist(individuals)

参数：
	individual:用来计算拥挤距离的个体列表
返回：
	没有返回值，拥挤距离保存在传入的individuals中的每个个体的individual.fitness.crowding_dist属性中

3.精英保存策略(Elitism)

比较操作
根据快速非支配排序和拥挤距离计算的结果，对族群中的个体进行排序：


在每个迭代步的最后，将父代与子代合为一个族群，依照比较操作对合并后族群中的个体进行排序，然后从中选取数量等同于父代规模的优秀子代，这就是NSGA-II算法中的精英保存策略。
DEAP实现
DEAP内置了实现NSGA-II中基于拥挤度的选择函数用来实现精英保存策略

tools.selNSGA2(individuals, k, nd='standard')

参数：
	individuals:被选择的个体列表，在NSGA-II中是父代与子代的并集
	k：需要选择的个体个数，通常要比被选择的个体数目小；如果与被选择的个体数量相同，那么效果
		等同于基于pareto前沿的排序
	nd：选择使用的非支配(non-dominated)算法，可以设置为'standard'或'log'
返回：
	被选择的个体列表

三.NSGA-II算法实现

1.测试函数

#!usr/bin/env python
# -*- coding:utf-8 _*-
"""
@author: liujie
@software: PyCharm
@file: NSGA-II实现.py
@time: 2020/11/27 16:17
"""
'''
NSGA-II算法实现
'''
import random
import numpy as np
import matplotlib.pyplot as plt
from deap import creator,tools,base,algorithms


# 定义问题
creator.create('MultiObjMin',base.Fitness,weights=(-1.0,-1.0))
creator.create('Individual',list,fitness=creator.MultiObjMin)

# 个体编码
def uniform(low,up):
    # 均匀分布生成个体
    return [random.uniform(a,b) for a,b in zip(low,up)]

pop_size = 100
NDim = 2
# 变量下界
low = [0] * NDim
# 变量上界
up = [1] * NDim
# 生成个体
toolbox = base.Toolbox()
toolbox.register('Attr_float',uniform,low,up)
toolbox.register('Individual',tools.initIterate,creator.Individual,toolbox.Attr_float)
# 生成种群
toolbox.register('Population',tools.initRepeat,list,toolbox.Individual)
pop = toolbox.Population(n = pop_size)

# ZDT3评价函数,ind长度为2
def ZDT3(ind):
    n = len(ind)
    f1 = ind[0]
    g = 1 + 9 * np.sum(ind[1:]) / (n-1)
    f2 = g * (1 - np.sqrt(ind[0]/g) - ind[0]/g * np.sin(10*np.pi*ind[0]))
    return f1,f2

# 注册工具
toolbox.register('evaluate',ZDT3)
# 锦标赛选择
toolbox.register('selectGen1',tools.selTournament,tournsize=2)
# selTournamentDCD(individuals, k)
toolbox.register('select',tools.selTournamentDCD)
# tools.cxSimulatedBinaryBounded(ind1, ind2, eta, low, up)
toolbox.register('mate',tools.cxSimulatedBinaryBounded,eta=20.0,low=low,up=up)
# 变异 - 多项式变异
toolbox.register('mutate',tools.mutPolynomialBounded,eta=20.0,low=low,up=up,indpb=1.0/NDim)

# 用数据记录算法迭代过程
# 创建统计信息对象
stats = tools.Statistics(key= lambda ind : ind.fitness.values)
stats.register('avg',np.mean)
stats.register('std',np.std)
stats.register('min',np.min)
stats.register('max',np.max)

# 创建日志对象
logbook = tools.Logbook()


# 遗传算法主程序
# 参数设置
maxGen = 50
cxProb = 0.7
mutateProb = 0.2

# 遗传算法迭代
# 第一代
fitnesses = map(toolbox.evaluate,pop)
for ind,fit in zip(pop,fitnesses):
    ind.fitness.values = fit
record = stats.compile(pop)
logbook.record(gen=0,**record)

# 快速非支配排序操作
fronts = tools.emo.sortNondominated(pop,k=pop_size)
# 将每个个体的适应度设置为pareto前沿的次序
for idx,front in enumerate(fronts):
    for ind in front:
        ind.fitness.values = (idx + 1),

# 创建子代
offspring = toolbox.selectGen1(pop,pop_size)    # 锦标赛选择
# algorithms.varAnd进化算法的一部分，仅应用变化部分（交叉和变异）,克隆了个体，因此返回的种群独立于输入种群
offspring = algorithms.varAnd(offspring,toolbox,cxProb,mutateProb)

# 第二代之后的迭代
for gen in range(1,maxGen+1):
    # 合并父代与子代
    combinedPop = pop + offspring
    # 评价族群-更新新族群的适应度
    fitnesses = map(toolbox.evaluate,combinedPop)
    for ind,fit in zip(combinedPop,fitnesses):
        ind.fitness.values = fit

    # 快速非支配排序
    fronts = tools.emo.sortNondominated(combinedPop,k=pop_size,first_front_only=False)
    
    # 拥挤距离计算
    for front in fronts:
        tools.emo.assignCrowdingDist(front)

    # 环境选择--精英保留
    pop = []
    for front in fronts:
        pop += front

    # 复制
    pop = toolbox.clone(pop)
    # 基于拥挤度的选择函数用来实现精英保存策略
    pop = tools.selNSGA2(pop,k=pop_size,nd='standard')

    # 创建子代
    offspring = toolbox.select(pop,pop_size)
    offspring = toolbox.clone(offspring)
    offspring = algorithms.varAnd(offspring,toolbox,cxProb,mutateProb)

    # 记录数据-将stats的注册功能应用于pop，并作为字典返回
    record = stats.compile(pop)
    logbook.record(gen = gen,**record)

# 输出计算过程
logbook.header = 'gen','avg','std','min','max'
print(logbook)

# 输出最优解
bestInd = tools.selBest(pop,1)[0]
bestFit = bestInd.fitness.values
print('当前最优解:',bestInd)
print('对应的函数最小值为:',bestFit)

front = tools.emo.sortNondominated(pop,len(pop))[0]
gridPop
for ind in front:
    plt.plot(ind.fitness.values[0],ind.fitness.values[1],'r.',ms=2)
plt.xlabel('f_1')
plt.ylabel('f_2')
plt.tight_layout()
plt.show()

运行结果

gen	avg     	std     	min      	max     
0  	2.21829 	2.4389  	0.0341521	9.1206  
1  	1.34439 	2.04583 	0.0341521	9.12091 
2  	0.902985	1.61735 	0.0229718	9.41342 
3  	0.682112	1.227   	-0.356872	9.41342 
4  	0.728242	1.54161 	-0.356872	9.41342 
5  	0.597391	1.24451 	-0.617481	9.33351 
6  	0.511788	1.05542 	-0.688099	9.63604 
7  	0.390928	0.565392	-0.688099	2.21917 
8  	0.367941	0.490075	-0.688686	2.11735 
9  	0.300265	0.450415	-0.70712 	1.08968 
10 	0.298629	0.459386	-0.708115	1.08968 
11 	0.304469	0.459981	-0.708115	1.36692 
12 	0.309114	0.467134	-0.734622	1.3744  
13 	0.325322	0.470062	-0.745101	1.3744  
14 	0.274998	0.445207	-0.763896	0.98771 
15 	0.269556	0.447157	-0.76485 	0.98771 
16 	0.258597	0.447632	-0.771994	0.98771 
17 	0.276563	0.435062	-0.771994	0.98771 
18 	0.271015	0.422953	-0.773203	0.9877  
19 	0.269312	0.44875 	-0.773225	0.992029
20 	0.263029	0.434753	-0.773321	0.992029
21 	0.261087	0.431257	-0.773321	0.992029
22 	0.253956	0.450363	-0.773321	0.992029
23 	0.262906	0.433012	-0.773321	0.997527
24 	0.26855 	0.431864	-0.773321	0.99865 
25 	0.277768	0.43518 	-0.773321	0.99865 
26 	0.268571	0.456225	-0.773349	0.998648
27 	0.262892	0.447244	-0.773349	0.998644
28 	0.267743	0.449131	-0.773349	0.998644
29 	0.263141	0.455621	-0.773349	0.999438
30 	0.258008	0.459385	-0.773365	0.999438
31 	0.272148	0.441205	-0.773365	0.999438
32 	0.265067	0.435414	-0.773368	0.999438
33 	0.263523	0.436773	-0.773368	0.999438
34 	0.267321	0.427655	-0.773368	0.998816
35 	0.258097	0.438233	-0.773368	0.998401
36 	0.258395	0.442298	-0.773368	0.998401
37 	0.264326	0.439118	-0.773368	0.998378
38 	0.269527	0.428407	-0.773368	0.998378
39 	0.260322	0.446813	-0.773368	0.998378
40 	0.272001	0.425237	-0.773368	0.998378
41 	0.272534	0.42689 	-0.773368	0.998372
42 	0.275443	0.433953	-0.773368	1.53123 
43 	0.272329	0.450508	-0.773368	1.53123 
44 	0.266652	0.437602	-0.773368	1.00659 
45 	0.263218	0.455956	-0.773368	1.00659 
46 	0.256601	0.458183	-0.773368	0.9984  
47 	0.258821	0.449391	-0.773368	0.9984  
48 	0.254844	0.451403	-0.773368	0.9984  
49 	0.266876	0.442811	-0.773368	0.9984  
50 	0.261986	0.447793	-0.773368	0.998367
当前最优解: [2.6721132552406934e-06, 2.0752896061476544e-07]
对应的函数最小值为: (2.6721132552406934e-06, 0.998367206028216)

可视化

可以看出NSGA-II算法得到的pareto最优前沿质量很高：最优解均匀分布在不连续前沿的各个线段上；同时在最优前沿以外没有个体存在。