手把手Python 实现 NSGA-III 算法

2025-2-7

给大家介绍一个多目标优化算法界的 “超级巨星” ——NSGA-III（Non-dominated Sorting Genetic Algorithm III）！这个算法真的非常厉害，它不仅能在多个目标之间找到最优解，还能优雅地处理复杂的多目标优化问题。要是你还不知道它，那可就太 “out” 了！今天，我就带大家深入了解 NSGA-III 算法，并用 Python 实现它，保证让你一看就懂！

初识 NSGA-III：多目标优化的“神器”

首先，咱们先来说说 NSGA-III 为啥这么牛。在现实生活中，很多问题都涉及到多目标优化，比如设计一辆汽车，既要追求高速度，又要保证低油耗，还要考虑成本。这时候，我们可能需要权衡这些目标，找到一个折中的最优解。NSGA-III 就是专门用来解决这种多目标优化问题的。

传统的单目标优化算法，比如梯度下降法，只能处理一个目标。而 NSGA-III 可以同时处理多个目标，通过模拟自然选择和遗传进化的过程，逐步逼近最优解。而且，NSGA-III 还引入了参考点和参考方向的概念，使得优化过程更加高效和稳定。

深入 NSGA-III 算法原理

算法核心流程

那 NSGA-III 算法到底是怎么运行的呢？咱们来一步步分解：

初始化种群：首先，我们需要随机生成一批解作为初始种群。这就像一群种子，我们将要通过算法的培育，让这些种子成长为最优解。
非支配排序：在多目标优化中，没有一个解一定比其他解好，因此，我们需要对解进行非支配排序。简单来说，就是把解按照它们的好坏（即是否被其他解支配）分为不同的层级，这样可以更好地筛选出优质的解。
拥挤度计算：为了保持种群的多样性，我们需要计算每个解的拥挤度。拥挤度大的解表示周围有较多的解，这样的解可能会被剔除；而拥挤度小的解则表示周围解较少，更有可能被保留下来。
选择、交叉和变异：接下来，我们通过选择操作，从当前种群中挑选出优质的解，然后进行交叉和变异操作，生成新的子代种群。这就像生物进化中的遗传和变异过程，通过这种方式，我们可以在解空间中不断探索新的可能。
精英保留：为了保留优秀的解，我们将当前种群和子代种群合并，然后再次进行非支配排序和拥挤度计算，最终选出下一代种群。这样可以确保每一代种群都比上一代更优秀。

通过这样的循环迭代，NSGA-III 算法就能逐渐逼近最优解。

参考点与参考方向

NSGA-III 算法的核心之一就是参考点和参考方向。在优化过程中，参考点可以看作是我们期望的最优解的目标值，而参考方向则指导算法在解空间中向哪个方向搜索。通过这种方式，NSGA-III 算法可以更有效地找到 Pareto 前沿（即所有非支配解组成的集合）。

NSGA-III 算法简介

NSGA-III（Non-dominated Sorting Genetic Algorithm III）是一种用于多目标优化的进化算法，由 Deb 和 Jain 在 2014 年提出。它是 NSGA-II 的扩展，旨在解决高维多目标优化问题（目标数大于 3）。NSGA-III 的核心思想是通过引入参考点（Reference Points）和参考方向（Reference Directions）来引导种群的进化，从而在 Pareto 前沿上均匀分布解。

Python 实现

以下是 NSGA-III 算法的 Python 实现，我们将从头开始，不依赖外部库（除了基础库）。

1. 导入基础库

import numpy as npimport matplotlib.pyplot as pltfrom sklearn.cluster import KMeans

2. 个体类

定义一个个体类，用于存储个体的基因、目标值、非支配层级和拥挤度。

class Individual:    def __init__(self, genes):        self.genes = genes  # 解的基因        self.objectives = None  # 解的目标值        self.rank = None  # 非支配层级        self.crowding_distance = 0  # 拥挤度        self.ref_dir = None  # 关联的参考方向

3. 初始化种群

生成随机种群。

def initialize_population(pop_size, num_vars, lower_bound, upper_bound):    population = []    for _ in range(pop_size):        genes = np.random.uniform(lower_bound, upper_bound, size=num_vars)        individual = Individual(genes)        population.append(individual)    return population

4. 非支配排序

实现非支配排序。

def non_dominated_sort(population):    for individual in population:        individual.rank = None        individual.domination_count = 0        individual.dominated_solutions = []
    for i in range(len(population)):        for j in range(len(population)):            if i == j:                continue            if dominates(population[i], population[j]):                population[i].dominated_solutions.append(population[j])            elif dominates(population[j], population[i]):                population[i].domination_count += 1
    front = [individual for individual in population if individual.domination_count == 0]    current_rank = 1    F = [[]]
    while front:        for individual in front:            individual.rank = current_rank        F.append(front)        next_front = []        for individual in front:            for dominated in individual.dominated_solutions:                dominated.domination_count -= 1                if dominated.domination_count == 0:                    next_front.append(dominated)        front = [individual for individual in next_front if individual not in front]        current_rank += 1
    return F[:-1]  # 去掉空的最后一层

5. 拥挤度计算

实现拥挤度计算，用于保持种群多样性。

def calculate_crowding_distance(population):    distance = np.zeros(len(population))    num_objs = len(population[0].objectives)
    for i, obj in enumerate(zip(*[ind.objectives for ind in population])):        obj = np.array(obj)        sorted_idx = np.argsort(obj)        distance[sorted_idx[0]] = np.inf        distance[sorted_idx[-1]] = np.inf        normalized_obj = (obj - np.min(obj)) / (np.max(obj) - np.min(obj) + 1e-10)        for j in range(1, len(sorted_idx) - 1):            distance[sorted_idx[j]] += normalized_obj[sorted_idx[j+1]] - normalized_obj[sorted_idx[j-1]]
    return distance

6. 参考方向生成

使用 Das-Dennis 方法生成参考方向。

def generate_reference_directions(num_objs, num_partitions):    if num_objs == 1:        return np.array([[1.0]])    ref_dirs = []    step = 1.0 / num_partitions    for i in range(num_partitions + 1):        xi = i * step        remaining = 1.0 - xi        if remaining < 0:            continue        if num_objs > 2:            sub_dirs = generate_reference_directions(num_objs - 1, num_partitions)            for sub_dir in sub_dirs:                ref_dir = np.array([xi] + list(sub_dir))                ref_dirs.append(ref_dir)        else:            ref_dirs.append(np.array([xi, 1.0 - xi]))    return np.array(ref_dirs)

7. 关联参考方向

为每个解找到最近的参考方向。

def associate_to_reference_directions(population, ref_dirs):    for individual in population:        min_distance = np.inf        best_dir = None        for ref_dir in ref_dirs:            distance = np.linalg.norm(individual.objectives - ref_dir)            if distance < min_distance:                min_distance = distance                best_dir = ref_dir        individual.ref_dir = best_dir

8. 选择、交叉、变异

实现遗传操作。

def select_tournament(population, tournament_size):    selected = []    for _ in range(tournament_size):        candidate = np.random.choice(population)        selected.append(candidate)    return min(selected, key=lambda x: (x.rank, -x.crowding_distance))
def crossover(parent1, parent2, eta_c):    child1_genes = np.zeros_like(parent1.genes)    child2_genes = np.zeros_like(parent1.genes)    for i in range(len(parent1.genes)):        if parent1.genes[i] > parent2.genes[i]:            uj = np.random.random()            if uj < 0.5:                beta = (2 * uj) ** (1 / (eta_c + 1))            else:                beta = (1 / (2 - 2 * uj)) ** (1 / (eta_c + 1))            child1_genes[i] = 0.5 * ((1 + beta) * parent1.genes[i] + (1 - beta) * parent2.genes[i])            child2_genes[i] = 0.5 * ((1 - beta) * parent1.genes[i] + (1 + beta) * parent2.genes[i])    return child1_genes, child2_genes
def mutate(child, eta_m, mutation_rate):    for i in range(len(child.genes)):        if np.random.random() < mutation_rate:            delta = np.random.random() ** (1 / (eta_m + 1))            child.genes[i] += delta * (child.genes[i] - 0.5)    return child

9. 主算法

整合所有步骤，实现 NSGA-III。

def nsga3(pop_size, num_vars, num_objs, lower_bound, upper_bound, problem, num_generations, ref_dirs):    population = initialize_population(pop_size, num_vars, lower_bound, upper_bound)
    for generation in range(num_generations):        # 计算目标值        for individual in population:            objectives = problem(individual.genes)            individual.objectives = objectives
        # 非支配排序和拥挤度        fronts = non_dominated_sort(population)        for front in fronts:            crowding_distance = calculate_crowding_distance(front)            for i, individual in enumerate(front):                individual.crowding_distance = crowding_distance[i]
        # 关联参考方向        associate_to_reference_directions(population, ref_dirs)
        # 选择、交叉、变异        offspring = []        while len(offspring) < pop_size:            parent1 = select_tournament(population, tournament_size=2)            parent2 = select_tournament(population, tournament_size=2)            child1_genes, child2_genes = crossover(parent1, parent2, eta_c=2)            child1 = Individual(child1_genes)            child2 = Individual(child2_genes)            child1 = mutate(child1, eta_m=20, mutation_rate=0.1)            child2 = mutate(child2, eta_m=20, mutation_rate=0.1)            offspring.extend([child1, child2])
        # 合并种群        combined_population = population + offspring
        # 淘汰        fronts = non_dominated_sort(combined_population)        new_population = []        for front in fronts:            if len(new_population) + len(front) > pop_size:                sorted_front = sorted(front, key=lambda x: (x.rank, -x.crowding_distance))                new_population.extend(sorted_front[:pop_size - len(new_population)])                break            new_population.extend(front)        population = new_population
    return population

10. 测试问题

定义一个简单的多目标优化问题（例如 DTLZ1）。

def dtlz1(x):    n_vars = len(x)    k = n_vars - 1    g = 100 * (k + sum((x_i - 0.5)**2 - np.cos(20 * np.pi * (x_i - 0.5)) for x_i in x[:-1]))    f1 = 0.5 * x[-1] * (1 + g)    f2 = 0.5 * (1 - x[-1]) * (1 + g)    return np.array([f1, f2])

11. 可视化结果

运行算法并可视化结果。

# 参数设置pop_size = 100num_vars = 2num_objs = 2lower_bound = np.array([0.0, 0.0])upper_bound = np.array([1.0, 1.0])num_generations = 100
# 生成参考方向ref_dirs = generate_reference_directions(num_objs, num_partitions=12)
# 运行 NSGA-IIIpopulation = nsga3(pop_size, num_vars, num_objs, lower_bound, upper_bound, dtlz1, num_generations, ref_dirs)
# 提取目标值objectives = [ind.objectives for ind in population]
# 可视plt.scatter([o[0] for o in objectives], [o[1] for o in objectives])plt.xlabel("Objective 1")plt.ylabel("Objective 2")plt.title("NSGA-III Results")plt.show()

深度解析

1. 关键代码解析

非支配排序

def non_dominated_sort(population):    for individual in population:        individual.rank = None        individual.domination_count = 0        individual.dominated_solutions = []
    for i in range(len(population)):        for j in range(len(population)):            if i == j:                continue            if dominates(population[i], population[j]):                population[i].dominated_solutions.append(population[j])            elif dominates(population[j], population[i]):                population[i].domination_count += 1
    front = [individual for individual in population if individual.domination_count == 0]    current_rank = 1    F = [[]]
    while front:        for individual in front:            individual.rank = current_rank        F.append(front)        next_front = []        for individual in front:            for dominated in individual.dominated_solutions:                dominated.domination_count -= 1                if dominated.domination_count == 0:                    next_front.append(dominated)        front = [individual for individual in next_front if individual not in front]        current_rank += 1
    return F[:-1]  # 去掉空的最后一层

非支配排序：通过比较每个个体和其他个体的关系，确定非支配层级。
支配关系：如果个体 A 的所有目标值都不劣于个体 B，且至少有一个目标值优于个体 B，则认为个体 A 支配个体 B。

拥挤度计算

def calculate_crowding_distance(population):    distance = np.zeros(len(population))    num_objs = len(population[0].objectives)
    for i, obj in enumerate(zip(*[ind.objectives for ind in population])):        obj = np.array(obj)        sorted_idx = np.argsort(obj)        distance[sorted_idx[0]] = np.inf        distance[sorted_idx[-1]] = np.inf        normalized_obj = (obj - np.min(obj)) / (np.max(obj) - np.min(obj) + 1e-10)        for j in range(1, len(sorted_idx) - 1):            distance[sorted_idx[j]] += normalized_obj[sorted_idx[j+1]] - normalized_obj[sorted_idx[j-1]]
    return distance

拥挤度计算：通过目标值的排序和差值计算拥挤度，确保种群的多样性。

参考方向关联

def associate_to_reference_directions(population, ref_dirs):    for individual in population:        min_distance = np.inf        best_dir = None        for ref_dir in ref_dirs:            distance = np.linalg.norm(individual.objectives - ref_dir)            if distance < min_distance:                min_distance = distance                best_dir = ref_dir        individual.ref_dir = best_dir

参考方向：为每个解找到最近的参考方向，用于维持解在 Pareto 前沿上的均匀分布。

2. 算法流程图

3. 算法性能分析

多样性：通过拥挤度和参考方向确保解在 Pareto 前沿上的均匀分布。
收敛性：通过非支配排序引导种群向 Pareto 前沿收敛。
计算复杂度：非支配排序和拥挤度计算的复杂度分别为 O(MN^2) 和 O(MN)，其中 M 为目标数，N 为种群大小。

总结

NSGA-III 是一种强大的多目标优化算法，适用于高维目标优化问题。通过参考方向和非支配排序，NSGA-III 能够在 Pareto 前沿上生成分布均匀的解。

来源：进修编程

THE END

Ollama部署Deepseek-R1大模型实操教程

<<上一篇

用OpenCV和Python构建自己的图像分类标注工具

下一篇>>