手搓代码
简单的遗传算法实现(详细注释版):
import random# 定义适应度函数
def fitness_function(individual):# 这里我们使用一个简单的函数作为示例,实际应用中需要根据具体问题定义"""计算个体的适应度值参数:individual (list): 一个包含整数的列表,表示一个个体的染色体返回:int: 个体的适应度值,即染色体中所有整数的和示例:>>> fitness_function([1, 2, 3, 4, 5])15"""return sum(individual)# 定义遗传算法参数
population_size = 50 # 种群大小
chromosome_length = 10 # 染色体长度
mutation_rate = 0.01 # 变异率
generations = 100 # 种群进化的代数# 初始化种群
population = [] # 储存所有种群的一个列表
for _ in range(population_size): # for 循环初始化每个种群individual = [random.randint(0, 1) for _ in range(chromosome_length)] # 每个个体是一个列表,长度是染色体长度,数值是0/1population.append(individual)# 遗传算法主循环
for generation in range(generations): # 对每一代进行for循环# 遍历所有的个体,计算每个个体的适应度,储存为一个列表fitness_values = [fitness_function(individual) for individual in population] # 储存所有个体的自适应度# 选择操作, 根据个体的适应度值从当前种群中选择出一定数量的个体作为下一代的父母# random.choices 可以从指定的序列中随机选择多个元素,并根据指定的权重(这里是 fitness_values)来决定每个元素被选中的概率。k 参数指定了要选择的元素数量# 这里k设置为种群大小,则每个个体都可以成为父母。如果k小于种群大小,那么只有部分个体有可能成为父母parents = random.choices(population, weights=fitness_values, k=population_size)# 交叉操作offspring = [] # 储存所有的子代个体for i in range(0, population_size, 2): # 前后相邻的两个个体作为父母parent1 = parents[i]parent2 = parents[i + 1]# crossover_point: 交叉的点位crossover_point = random.randint(1, chromosome_length - 1) # 产生一个在 [1, 染色体长度-1] 之间的整数child1 = parent1[:crossover_point] + parent2[crossover_point:] # “+”运算符用于连接两个列表child2 = parent2[:crossover_point] + parent1[crossover_point:]offspring.extend([child1, child2]) # extend 可以增加多个元素到列表的末尾# 子代中的每个个体有几率进行变异操作for i in range(population_size):if random.random() < mutation_rate:# mutation_point: 变异的点位mutation_point = random.randint(0, chromosome_length - 1)offspring[i][mutation_point] = 1 - offspring[i][mutation_point]# 更新种群population = offspring# 输出最终种群中的自适应度最优个体
best_individual = max(population, key=fitness_function) # key 关键字可选择,是一个可以迭代应用到population每个元素上的函数,之后再对结果返回最大值
print("Best individual:", best_individual)
print("Fitness value:", fitness_function(best_individual))
使用Geatpy库实现
安装 Geatpy2
pip install geatpy
- 解决多目标问题
"""MyProblem.py"""import numpy as np
import geatpy as eaclass MyProblem(ea.Problem): # Inherited from Problem class."""定义自己的优化问题这个自定义问题类用于定义一个具有特定名称、目标数量、决策变量维度和类型、以及边界条件的优化问题。在实际应用中,你可以根据需要扩展这个类,例如添加目标函数、约束条件等。"""def __init__(self, M): # M is the number of objects. M 表示目标的数量name = "DTLZ1" # Problem's name. 设置问题的名字,自定义maxormins = [1] * M # All objects are need to be minimized. 目标是最小化,所以是1,最大化是-1,这里是最小化,所以是1Dim = (M + 4) # Set the dimension of decision variables. 决策变量的维数,M+4 是决策变量的维数,M 是目标的数量varTypes = [0] * Dim # Set the types of decision variables. 0 means continuous while 1 means discrete. 0 表示连续,1 表示离散lb = [0] * Dim # The lower bound of each decision variable. 设置每个决策变量的下界为0。ub = [1] * Dim # The upper bound of each decision variable. 设置每个决策变量的上界为1。lbin = [1] * Dim # Whether the lower boundary is included. 设置每个决策变量的下界是否包含在内,这里都是包含的。ubin = [1] * Dim # Whether the upper boundary is included. 设置每个决策变量的上界是否包含在内,这里都是包含的。# Call the superclass's constructor to complete the instantiation 调用父类ea.Problem的构造函数,完成实例化。ea.Problem.__init__(self, name, M, maxormins, Dim, varTypes, lb, ub, lbin, ubin)def aimFunc(self, pop): # Write the aim function here, pop is an object of Population class."""这段代码定义了一个名为 aimFunc 的方法** 它是一个目标函数 **用于计算给定种群的目标值并将结果存储在 pop.ObjV 中参数:- pop 是一个种群对象,包含了一组个体(解)。- pop.Phen 是种群的决策变量矩阵"""# 获取决策变量Vars = pop.Phen # Get the decision variables# 计算辅助变量 g,g是一个与决策变量相关的惩罚项XM = Vars[:, (self.M - 1) :] # 取出矩阵的 M-1 列到最后一列g = np.array([100 * (self.Dim - self.M + 1 + np.sum(((XM - 0.5) ** 2 - np.cos(20 * np.pi * (XM - 0.5))), 1))]).Tones_metrix = np.ones((Vars.shape[0], 1)) # 定义辅助矩阵:创建一个与 Vars 行数相同,列数为1 的全1列向量。# 计算目标函数值,赋值给 pop.ObjVpop.ObjV = (0.5 * np.fliplr(np.cumprod(np.hstack([ones_metrix, Vars[:, : self.M - 1]]), 1)) * np.hstack([ones_metrix, 1 - Vars[:, range(self.M - 2, -1, -1)]]) * np.tile(1 + g, (1, self.M)))def calReferObjV(self): # Calculate the theoretic global optimal solution here. 该方法用于计算理论上的全局最优解# self.M 表示目标的数量uniformPoint, ans = ea.crtup(self.M, 10000) # create 10000 uniform points. 创建10000个均匀分布的点。realBestObjV = uniformPoint / 2 # 计算每个目标的理论最优解return realBestObjV # Return the global optimal objective value. 返回全局最优的目标函数值。M = 3 # Set the number of objects.
problem = MyProblem(M) # Instantiate MyProblem class
# Instantiate a algorithm class. 设置算法类,非支配排序遗传算法NSGA-III (多目标进化算法)
algorithm = ea.moea_NSGA3_templet(problem,ea.Population(Encoding="RI", NIND=100), # Set 100 individuals. 设置种群大小为100。MAXGEN=500, # Set the max iteration number. 设置最大迭代次数。logTras=1, # Set the frequency of logging. If it is zero, it would not log. 设置日志记录的频率,如果为0,则不记录日志。
) # 实例化一个算法模板对象。# Do the optimization problem. 运行算法模板,得到帕累托最优解集
res = ea.optimize(algorithm, verbose=False, drawing=1, outputMsg=True, drawLog=True, saveFlag=True
) # 执行算法模板,得到帕累托最优解集。
- 单目标实现
import numpy as np
import geatpy as ea
class Ackley(ea.Problem): # Inherited from Problem class.def __init__(self, D = 30):# 初始化的变量和多目标的一样,变量的值根据需要修改name = 'Ackley' # Problem's name.M = 1 # Set the number of objects.maxormins = [1] * M # All objects are need to be minimized.Dim = D # Set the dimension of decision variables. 根据需要修改,决策变量的维度varTypes = [0] * Dim # Set the types of decision variables. 0 means continuous while 1 means discrete.lb = [-32.768] * Dim # The lower bound of each decision variable. 每个决策变量的下界ub = [32.768] * Dim # The upper bound of each decision variable. 每个决策变量的上界lbin = [1] * Dim # Whether the lower boundary is included. 每个决策变量的下界是否包含ubin = [1] * Dim # Whether the upper boundary is included. 每个决策变量的上界是否包含# Call the superclass's constructor to complete the instantiationea.Problem.__init__(self, name, M, maxormins, Dim, varTypes, lb, ub, lbin, ubin)def aimFunc(self, pop): # Write the aim function here, pop is an object of Population class. 定义目标函数x = pop.Phen # Get the decision variables 获取决策变量n = self.Dim # 决策变量的维度f = np.array([-20 * np.exp(-0.2*np.sqrt(1/n*np.sum(x**2, 1))) - np.exp(1/n * np.sum(np.cos(2 * np.pi * x), 1)) + np.e + 20]).T# CV用于表示约束条件。 这里没有约束条件,因此CV的值为0, 详见Geatpy数据结构。实际中,需要定义一个函数来计算每个个体的违反约束程度,并将其作为CV返回# 假设我们有一个约束条件:所有决策变量的和不能超过某个值constraint = np.sum(x, axis=1) - 10 # 假设约束值为10CV = np.where(constraint > 0, constraint, 0) # 如果违反约束,CV为正值,否则为0return f, CVdef calReferObjV(self): # Calculate the global optimal solution here. 计算全局最优解realBestObjV = np.array([[0]])return realBestObjVproblem = Ackley(30)
# Instantiate a algorithm class. 差分进化(DE)算法
algorithm = ea.soea_DE_rand_1_bin_templet(problem,ea.Population(Encoding='RI', NIND=20),MAXGEN=1000, # Set the max times of iteration.logTras=1) # Set the frequency of logging. If it is zero, it would not log.
algorithm.mutOper.F = 0.5 # Set the F of DE
algorithm.recOper.XOVR = 0.2 # Set the Cr of DE (Here it is marked as XOVR)
# Do the optimization
res = ea.optimize(algorithm, verbose=False, drawing=1, outputMsg=True, drawLog=True, saveFlag=True, dirName='result')