告别调参玄学：用Python手把手实现NSGA-II多目标优化（附完整代码）-洪萨配资

告别调参玄学：用Python手把手实现NSGA-II多目标优化（附完整代码）

在工程优化和算法研究中，我们常常面临多个相互冲突的目标需要同时优化。比如在机器学习模型调优中，我们既希望模型准确率尽可能高，又希望推理速度足够快；在资源调度场景下，既要考虑计算资源利用率最大化，又要保证任务响应时间最短。这类多目标优化问题就像走钢丝，需要在多个目标之间找到最佳平衡点。

传统单目标优化方法在这里显得力不从心，而NSGA-II（非支配排序遗传算法II）作为多目标优化领域的经典算法，通过独特的非支配排序和拥挤度机制，能够高效地找到一组最优折衷解（Pareto前沿）。本文将从一个具体的优化案例出发，带你从零实现NSGA-II算法，理解每个步骤背后的设计思想，并可视化展示优化结果。

1. 多目标优化基础与NSGA-II核心思想

1.1 什么是Pareto最优

在多目标优化中，当我们要优化多个目标时，往往不存在一个在所有目标上都是最优的解。这时我们引入Pareto最优的概念：

支配关系：解A支配解B，当且仅当解A在所有目标上都不差于解B，且至少在一个目标上严格优于解B
Pareto最优解：如果一个解不被任何其他解支配，那么它就是Pareto最优解
Pareto前沿：所有Pareto最优解在目标空间形成的曲面

# 示例：判断解x1是否支配解x2 def is_dominated(x1, x2, objectives): """x1是否支配x2""" better_in_any = False for f in objectives: if f(x1) > f(x2): # 假设都是最大化问题 return False if f(x1) < f(x2): better_in_any = True return better_in_any

1.2 NSGA-II算法框架

NSGA-II相比原始NSGA有三项重要改进：

快速非支配排序：将计算复杂度从O(mN³)降到O(mN²)
精英保留策略：合并父代和子代种群，保留优秀个体
拥挤度比较：取代共享半径机制，更好地维持种群多样性

算法主要流程如下：

初始化种群
生成子代种群（选择、交叉、变异）
合并父代和子代种群
快速非支配排序
拥挤度计算
选择新一代种群
重复2-6直到满足终止条件

2. Python实现NSGA-II核心组件

2.1 种群初始化与编码

我们采用实数编码方式初始化种群。对于每个个体，随机生成其染色体（决策变量）：

import numpy as np def initialize_population(pop_size, var_num, bounds): """ 初始化种群 :param pop_size: 种群大小 :param var_num: 变量个数 :param bounds: 每个变量的上下界 [(low1, high1), (low2, high2)...] :return: 初始种群 """ population = np.zeros((pop_size, var_num)) for i in range(var_num): population[:, i] = np.random.uniform(bounds[i][0], bounds[i][1], pop_size) return population

2.2 快速非支配排序实现

快速非支配排序是NSGA-II的核心步骤，其关键在于高效计算每个个体的支配关系：

def fast_non_dominated_sort(population, objectives): """ 快速非支配排序 :param population: 种群 :param objectives: 目标函数列表 :return: 前沿等级列表 """ pop_size = len(population) S = [[] for _ in range(pop_size)] # 被个体p支配的解集 n = [0] * pop_size # 支配个体p的解数量 rank = [0] * pop_size # 前沿等级 # 第一遍遍历计算支配关系 for i in range(pop_size): for j in range(i+1, pop_size): if is_dominated(population[i], population[j], objectives): S[i].append(j) n[j] += 1 elif is_dominated(population[j], population[i], objectives): S[j].append(i) n[i] += 1 # 找出第一前沿 current_rank = 1 fronts = [[]] for i in range(pop_size): if n[i] == 0: rank[i] = current_rank fronts[0].append(i) # 迭代计算后续前沿 while fronts[-1]: next_front = [] for i in fronts[-1]: for j in S[i]: n[j] -= 1 if n[j] == 0: rank[j] = current_rank + 1 next_front.append(j) current_rank += 1 if next_front: fronts.append(next_front) return fronts, rank

2.3 拥挤度计算与比较

拥挤度用于衡量个体在目标空间的分布密度，保证前沿解的多样性：

def crowding_distance(front, population, objectives): """ 计算前沿中个体的拥挤度 :param front: 前沿个体索引列表 :param population: 种群 :param objectives: 目标函数列表 :return: 拥挤度列表 """ distances = [0] * len(front) if not front: return distances # 对每个目标函数分别处理 for obj_idx in range(len(objectives)): # 根据当前目标函数值排序 sorted_front = sorted(front, key=lambda x: objectives[obj_idx](population[x])) # 边界个体距离设为无穷大 distances[sorted_front[0]] = float('inf') distances[sorted_front[-1]] = float('inf') # 归一化目标函数值 min_obj = objectives[obj_idx](population[sorted_front[0]]) max_obj = objectives[obj_idx](population[sorted_front[-1]]) norm = max_obj - min_obj if max_obj != min_obj else 1 # 计算中间个体的拥挤度 for i in range(1, len(front)-1): idx = sorted_front[i] next_obj = objectives[obj_idx](population[sorted_front[i+1]]) prev_obj = objectives[obj_idx](population[sorted_front[i-1]]) distances[idx] += (next_obj - prev_obj) / norm return distances

3. 遗传算子实现

3.1 选择操作：锦标赛选择

我们采用基于Pareto等级和拥挤度的锦标赛选择：

def tournament_selection(population, rank, crowding_dist, tournament_size=2): """ 锦标赛选择 :param population: 种群 :param rank: 前沿等级 :param crowding_dist: 拥挤度 :param tournament_size: 锦标赛规模 :return: 被选中的个体索引 """ pop_size = len(population) # 随机选择tournament_size个参赛者 contestants = np.random.choice(pop_size, tournament_size, replace=False) # 选择Pareto等级最高的 min_rank = min(rank[c] for c in contestants) candidates = [c for c in contestants if rank[c] == min_rank] # 如果多个同等级，选择拥挤度大的 if len(candidates) > 1: max_dist = max(crowding_dist[c] for c in candidates) candidates = [c for c in candidates if crowding_dist[c] == max_dist] # 如果仍然多个，随机选择一个 return np.random.choice(candidates)

3.2 交叉操作：模拟二进制交叉(SBX)

SBX能够产生靠近父代个体的子代，同时保持种群多样性：

def sbx_crossover(parent1, parent2, bounds, eta_c=20): """ 模拟二进制交叉 :param parent1: 父代1 :param parent2: 父代2 :param bounds: 变量边界 :param eta_c: 交叉分布指数 :return: 两个子代 """ child1, child2 = parent1.copy(), parent2.copy() for i in range(len(parent1)): if np.random.random() < 0.5: # 50%概率进行交叉 if abs(parent1[i] - parent2[i]) > 1e-14: # 避免除零 x1 = min(parent1[i], parent2[i]) x2 = max(parent1[i], parent2[i]) x_low, x_high = bounds[i] # 计算beta值 rand = np.random.random() beta = 1.0 + (2.0 * (x1 - x_low) / (x2 - x1)) alpha = 2.0 - beta ** -(eta_c + 1) if rand <= 1.0 / alpha: beta_q = (rand * alpha) ** (1.0 / (eta_c + 1)) else: beta_q = (1.0 / (2.0 - rand * alpha)) ** (1.0 / (eta_c + 1)) # 生成子代 c1 = 0.5 * (x1 + x2 - beta_q * (x2 - x1)) c2 = 0.5 * (x1 + x2 + beta_q * (x2 - x1)) # 确保在边界内 c1 = min(max(c1, x_low), x_high) c2 = min(max(c2, x_low), x_high) if np.random.random() < 0.5: # 随机交换 child1[i], child2[i] = c1, c2 else: child1[i], child2[i] = c2, c1 return child1, child2

3.3 变异操作：多项式变异

多项式变异能够在父代附近产生多样性，同时偏向于产生接近父代的解：

def polynomial_mutation(individual, bounds, eta_m=20, mutation_prob=1.0): """ 多项式变异 :param individual: 个体 :param bounds: 变量边界 :param eta_m: 变异分布指数 :param mutation_prob: 变异概率 :return: 变异后的个体 """ mutated = individual.copy() for i in range(len(individual)): if np.random.random() < mutation_prob: x = individual[i] x_low, x_high = bounds[i] delta1 = (x - x_low) / (x_high - x_low) delta2 = (x_high - x) / (x_high - x_low) rand = np.random.random() mut_pow = 1.0 / (eta_m + 1.0) if rand < 0.5: xy = 1.0 - delta1 val = 2.0 * rand + (1.0 - 2.0 * rand) * xy ** (eta_m + 1) delta_q = val ** mut_pow - 1.0 else: xy = 1.0 - delta2 val = 2.0 * (1.0 - rand) + 2.0 * (rand - 0.5) * xy ** (eta_m + 1) delta_q = 1.0 - val ** mut_pow x = x + delta_q * (x_high - x_low) x = min(max(x, x_low), x_high) mutated[i] = x return mutated

4. 完整NSGA-II算法实现与案例应用

4.1 算法主流程

将上述组件组合成完整的NSGA-II算法：

def nsga2(pop_size, var_num, bounds, objectives, max_gen=100, crossover_prob=0.9, mutation_prob=1.0): """ NSGA-II主算法 :param pop_size: 种群大小 :param var_num: 变量个数 :param bounds: 变量边界 :param objectives: 目标函数列表 :param max_gen: 最大迭代次数 :param crossover_prob: 交叉概率 :param mutation_prob: 变异概率 :return: 最终种群，前沿历史 """ # 初始化 population = initialize_population(pop_size, var_num, bounds) history = [] for gen in range(max_gen): # 评估目标函数 obj_values = [[f(ind) for f in objectives] for ind in population] # 非支配排序 fronts, rank = fast_non_dominated_sort(population, objectives) # 计算拥挤度 crowding_dist = [0] * pop_size for front in fronts: front_dist = crowding_distance(front, population, objectives) for i, dist in zip(front, front_dist): crowding_dist[i] = dist # 记录前沿 history.append([population[front[0]] for front in fronts if front]) # 生成子代 offspring = [] while len(offspring) < pop_size: # 选择 parent1_idx = tournament_selection(population, rank, crowding_dist) parent2_idx = tournament_selection(population, rank, crowding_dist) parent1 = population[parent1_idx] parent2 = population[parent2_idx] # 交叉 if np.random.random() < crossover_prob: child1, child2 = sbx_crossover(parent1, parent2, bounds) else: child1, child2 = parent1.copy(), parent2.copy() # 变异 child1 = polynomial_mutation(child1, bounds, mutation_prob=mutation_prob) child2 = polynomial_mutation(child2, bounds, mutation_prob=mutation_prob) offspring.extend([child1, child2]) # 合并种群 combined_pop = np.vstack((population, offspring[:pop_size])) # 环境选择 new_population = [] remaining = len(combined_pop) fronts, _ = fast_non_dominated_sort(combined_pop, objectives) for front in fronts: if len(new_population) + len(front) <= pop_size: new_population.extend(combined_pop[front]) else: # 需要按拥挤度选择部分个体 front_dist = crowding_distance(front, combined_pop, objectives) sorted_front = sorted(zip(front, front_dist), key=lambda x: x[1], reverse=True) needed = pop_size - len(new_population) new_population.extend(combined_pop[f[0]] for f in sorted_front[:needed]) break population = np.array(new_population) return population, history

4.2 案例：优化机器学习模型超参数

假设我们要优化一个随机森林模型的两个目标：

最大化模型准确率
最小化预测延迟

定义目标函数：

from sklearn.ensemble import RandomForestClassifier from sklearn.datasets import load_iris from sklearn.model_selection import cross_val_score import time # 加载数据 iris = load_iris() X, y = iris.data, iris.target def evaluate_model(n_estimators, max_depth): """ 评估随机森林模型 :param n_estimators: 树的数量 :param max_depth: 最大深度 :return: (准确率, 延迟) """ model = RandomForestClassifier(n_estimators=int(n_estimators), max_depth=int(max_depth) if max_depth > 1 else None, random_state=42) # 计算准确率 accuracy = cross_val_score(model, X, y, cv=5).mean() # 计算预测延迟 model.fit(X, y) start = time.time() for _ in range(100): model.predict(X[:1]) latency = (time.time() - start) / 100 return accuracy, -latency # 第二个目标取负，转为最大化问题 # 定义NSGA-II的目标函数 def objective1(individual): acc, _ = evaluate_model(individual[0], individual[1]) return acc def objective2(individual): _, lat = evaluate_model(individual[0], individual[1]) return lat # 参数边界 bounds = [(10, 100), # n_estimators (1, 20)] # max_depth

运行优化：

# 运行NSGA-II final_pop, history = nsga2(pop_size=20, var_num=2, bounds=bounds, objectives=[objective1, objective2], max_gen=50) # 提取Pareto前沿 fronts, _ = fast_non_dominated_sort(final_pop, [objective1, objective2]) pareto_front = final_pop[fronts[0]]

4.3 结果可视化

使用Matplotlib可视化Pareto前沿：

import matplotlib.pyplot as plt # 计算目标函数值 obj1_values = [objective1(ind) for ind in final_pop] obj2_values = [objective2(ind) for ind in final_pop] pf_obj1 = [objective1(ind) for ind in pareto_front] pf_obj2 = [objective2(ind) for ind in pareto_front] # 绘制结果 plt.figure(figsize=(10, 6)) plt.scatter(obj1_values, obj2_values, c='blue', alpha=0.5, label='Population') plt.scatter(pf_obj1, pf_obj2, c='red', s=100, label='Pareto Front') plt.xlabel('Model Accuracy') plt.ylabel('Negative Latency (higher is better)') plt.title('NSGA-II Optimization Results for Random Forest') plt.legend() plt.grid(True) plt.show()

5. 工程实践建议与调优技巧

5.1 参数设置经验法则

根据实际项目经验，NSGA-II的关键参数设置可以参考以下范围：

参数	推荐范围	说明
种群大小	50-200	复杂问题需要更大种群
交叉概率	0.8-0.95	保持较高交叉概率
变异概率	1/n	n为变量个数，确保每个变量都有变异机会
交叉分布指数(η_c)	15-30	值越大子代越接近父代
变异分布指数(η_m)	15-30	同上