|
待优化函数
本文待优化函数选取自《MATLAB智能算法30个案例分析(第2版)》中的第一个案例
利用遗传算法计算以下函数的最小值:
f(x) = \frac{\sin(10 \pi x)}{x},x\in[1,2]
遗传算法流程
关于遗传算法的原理,书籍和文章均比较多,这里就不再赘述,这里给出简单遗传算法的流程
编码
这里着重说明一下编码方式,本文算法采用二进制编码。假设某一参数的取值范围是 [U_{min},U_{max}] ,我们用长度为 l 的二进制编码符号串来表示该参数,则它总共能够产生 2^l 种不同的编码,若使参数编码时的对应关系如下:
则二进制编码的编码精度为: \delta = \frac{U_{max} - U_{min}}{2^l - 1}
假设某一个体的编码是:X: b_lb_{l-1}b_{l-2}...b_2b_1
则对应的解码公式为:x=U_{min}+(\sum_{i=1}^l b_i \cdot 2^{i-1})\cdot \frac{U_{max}-U_{min}}{2^l - 1}
而 \sum_{i=1}^l b_i \cdot 2^{i-1} 正是二进制对应的十进制数。
程序及运行结果
关于遗传算法各阶段运算,包括选择(复制)运算、交叉运算、变异运算均有不同的实现,本文代码参考了《遗传算法原理及应用》附录中C语言实现的简单遗传算法,有兴趣的读者可以对以上各阶段运算尝试其他的实现方式。 代码如下:
import math
import random
import copy
import matplotlib.pyplot as plt
PI = 3.1415926
# 个体长度
CHROM_LEN = 20
# 种群大小
POP_SIZE = 40
CMIN = 0
# 最大遗传代数
MAX_GENERATION = 40
# 交叉概率
PC = 0.7
# 变异概率
PM = 0.01
# 优化函数
def F(x):
return math.sin(10 * PI * x) / x
# 解码器
def decode(chrom, lb, ub):
# 二进制对应的十进制数
temp = int(chrom, 2)
# 最终解码值
x = lb + temp * (ub - lb) / (math.pow(2, CHROM_LEN) - 1)
return x
# 个体类
class Individual:
def __init__(self):
temp = []
for _ in range(CHROM_LEN):
temp.append(random.randint(0, 1))
self.chrom = "".join([str(t) for t in temp])
self.fitness = 0
# 计算个体适应度
def get_fitness(self, lb, ub):
x = decode(self.chrom, lb, ub)
value = -F(x) + CMIN
self.fitness = max(0, value)
return self.fitness
def __str__(self):
return "chrom:{}, fitness:{}".format(self.chrom, self.fitness)
# 获得当代最佳和最差个体索引
def best_and_worst(population):
# 最佳个体索引
best_idx = 0
# 最差个体索引
worst_idx = 0
for _idx, p in enumerate(population):
if p.fitness > population[best_idx].fitness:
best_idx = _idx
elif p.fitness < population[worst_idx].fitness:
worst_idx = _idx
return best_idx, worst_idx
# 选择(复制)操作
def select(population):
# 新种群
new_pop = []
# 当代个体适应度总和
fitness_sum = max(sum([i.fitness for i in population]), 0.0001)
# 当代个体累计适应度占比
cfitness = []
# 计算相对适应度占比
for j in range(POP_SIZE):
cfitness.append(population[j].fitness / fitness_sum)
# 计算累计适应度占比
for j in range(POP_SIZE):
if j == 0:
continue
cfitness[j] = cfitness[j-1] + cfitness[j]
# 依据累计适应度占比进行选择复制,随机数大于对应的累计适应度占比,则进行复制
for k in range(POP_SIZE):
index = 0
while random.random() > cfitness[index]:
index += 1
# 若无法找到要复制的其他个体,则沿用当前个体
if index >= POP_SIZE:
index = k
break
new_pop.append(copy.deepcopy(population[index]))
return new_pop
# 交叉操作
def crossover(population):
# 随机产生个体配对索引,类似于洗牌的效果
index = [i for i in range(POP_SIZE)]
for i in range(POP_SIZE):
point = random.randint(0, POP_SIZE - i - 1)
temp = index
index = index[point + i]
index[point + i] = temp
for i in range(0, POP_SIZE, 2):
if random.random() > PC:
# 随机选择交叉开始位置
cross_start = random.randint(0, CHROM_LEN - 2) + 1
# 需要交换的基因
cross_gene1 = population[index].chrom[cross_start:]
cross_gene2 = population[index[i + 1]].chrom[cross_start:]
# 交叉操作
population[index].chrom = population[index].chrom[0: cross_start] + cross_gene2
population[index[i + 1]].chrom = population[index[i + 1]].chrom[0: cross_start] + cross_gene1
# 变异操作
def mutation(population):
for individual in population:
# 初始化新染色体
new_chrom_ch = [c for c in individual.chrom]
for i in range(CHROM_LEN):
# 随机数小于变异概率,则进行变异操作
if random.random() < PM:
new_chrom_ch = &#34;1&#34; if individual.chrom is &#34;0&#34; else &#34;0&#34;
# 更新染色体
individual.chrom = &#34;&#34;.join(new_chrom_ch)
# 绘制结果
def draw_result(best):
import numpy as np
# 绘制优化函数
x = np.linspace(1, 2, 100)
y = [F(_x) for _x in x]
plt.plot(x, y)
# 绘制最优解
best_x = decode(best.chrom, 1, 2)
best_y = F(decode(best.chrom, 1, 2))
plt.scatter(best_x, best_y, s=100, c=&#39;red&#39;, marker=&#39;*&#39;, zorder=2)
plt.show()
# plt.savefig(&#39;sga_result.png&#39;, dpi=800)
# 绘制进化过程
def draw_evolution(evolution):
x = [i for i in range(len(evolution))]
plt.plot(x, evolution)
plt.show()
# plt.savefig(&#39;sga_evolution.png&#39;, dpi=800)
def main():
# 种群
population = []
# 下界
lb = 1
# 上界
ub = 2
# 初始化种群
for _ in range(POP_SIZE):
population.append(Individual())
# 计算初始种群适应度
for individual in population:
individual.get_fitness(lb, ub)
# 初始种群最佳和最差个体
best_idx, worst_idx = best_and_worst(population)
# 历史最佳个体
current_best = population[best_idx]
# 进化过程,每一代的最佳个体的函数值
evolution = []
# 循环直到最大代数
for generation in range(MAX_GENERATION):
# 选择复制
population = select(population)
# 交叉
crossover(population)
# 变异
mutation(population)
# 重新计算适应度
for individual in population:
individual.get_fitness(lb, ub)
# 当代种群最佳和最差个体索引
best_idx, worst_idx = best_and_worst(population)
# 利用精英模型执行进化操作,用历史最佳个体代替当代的最差个体
if population[best_idx].fitness > current_best.fitness:
current_best = population[best_idx]
else:
population[worst_idx] = current_best
# 更新进化过程
evolution.append(round(F(decode(current_best.chrom, 1, 2)), 4))
# 绘制进化过程
# draw_evolution(evolution)
# 绘制结果
draw_result(current_best)
# 打印最佳结果
print(&#34;X = {}&#34;.format(round(decode(current_best.chrom, 1, 2), 4)))
print(&#34;Y = {}&#34;.format(round(F(decode(current_best.chrom, 1, 2)), 4)))
if __name__ == &#34;__main__&#34;:
main()代码输出最优解为:
X = 1.1491
Y = -0.8699待优化函数及最优解如下图所示:
每一代最优解的进化过程如下图所示(由于初始种群的随机性,每一次执行产生的结果可能会不同):
搜索最优解的动态图如下图所示(由于初始种群的随机性,每一次执行产生的结果可能会不同):
在这里插入图片描述
问题
由于初始种群的随机性,每一次得到的最优解可能会稍有差异,本文代码有时会找不到全局最优解,稳定性有待提升,在此作者抛砖引玉,希望有实力的读者能进一步优化并留言
笔者水平有限,若有不对的地方欢迎评论指正 |
本帖子中包含更多资源
您需要 登录 才可以下载或查看,没有账号?立即注册
×
|