遗传算法求解二维函数最大值（动态展示）

提示：基于前者代码的改进。原代码链接

根据前者提供的代码
在复现的过程中发现了一些改进的点（交叉和变异部分）
并且对每次迭代的结果进行了动态展示。

---

文章目录

• 遗传算法求解二维函数最大值（动态展示）
• 前言
• • 1.导入库
  • 2.定义变量
  • 3.完整代码（含解释）

---

前言

代码运行可得到每一次迭代结果的图形，即为动态寻找最大值的过程。

---

1.导入库

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

2.定义变量

数值可以自行调整。

DNA_SIZE = 24
POP_SIZE = 200
CROSSOVER_RATE = 0.8
MUTATION_RATE = 0.02
N_GENERATIONS = 50
X_BOUND = [-3, 3]
Y_BOUND = [-3, 3]

定义函数F(x, y)

def F(x, y):
    return 3*(1-x)**2*np.exp(-(x**2)-(y+1)**2) - 10*(x/5 - x**3 - y**5)*np.exp(-x**2-y**2) - 1/3**np.exp(-(x+1)**2 - y**2)

创建初始种群pop，POP_SIZE行DNA_SIZE*2列元素为0或1的二维矩阵。

pop = np.random.randint(0, 2, size=(POP_SIZE, DNA_SIZE*2))

构建网格坐标矩阵，画3D展示图是需要用，一维元组需要展开成二维矩阵。

X = np.linspace(X_BOUND[0], X_BOUND[1], 250)
Y = np.linspace(Y_BOUND[0], Y_BOUND[1], 250)
X, Y = np.meshgrid(X, Y)
Z = F(X, Y)

最后生成的三维点坐标即为(X, Y, Z)
函数np.meshgrid(X, Y)的详细说明和用法可参照这篇文章。

3.完整代码（含解释）

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

DNA_SIZE = 20
POP_SIZE = 200
CROSSOVER_RATE = 0.8
MUTATION_RATE = 0.05
N_GENERATIONS = 100
X_BOUND = [-3, 3]
Y_BOUND = [-3, 3]


def F(x, y):
    return 3*(1-x)**2*np.exp(-(x**2)-(y+1)**2) - 10*(x/5 - x**3 - y**5)*np.exp(-x**2-y**2) - 1/3**np.exp(-(x+1)**2 - y**2)


X = np.linspace(X_BOUND[0], X_BOUND[1], 250)
Y = np.linspace(Y_BOUND[0], Y_BOUND[1], 250)
# 生成网格点坐标矩阵
X, Y = np.meshgrid(X, Y)
# 将网格矩阵坐标带入方程，求出方程的值。
Z = F(X, Y)
# 构建3维网格点坐标(X, Y, Z) X是行坐标矩阵，Y是列坐标矩阵，Z是高坐标矩阵。


# 将二进制的DNA转换成实数值，再求解函数值（适应度值）
def translateDNA(pop):
    x_pop = pop[:, 1::2]  # 奇数列表示X
    y_pop = pop[:, ::2]  # 偶数列表示y
    # 矩阵A.dot(矩阵B) 矩阵A点乘矩阵B
    # 将x y 对应到定义的区间内
    trans_x = x_pop.dot(2**np.arange(DNA_SIZE)[::-1]) / \ float(2**DNA_SIZE-1)*(X_BOUND[1]-X_BOUND[0])+X_BOUND[0]
    trans_y = y_pop.dot(2**np.arange(DNA_SIZE)[::-1]) / \ float(2**DNA_SIZE-1)*(Y_BOUND[1]-Y_BOUND[0])+Y_BOUND[0]
    return trans_x, trans_y


# 求解种群中每一个个体的适应度
def get_fitness(pop):
    x, y = translateDNA(pop)
    pred = F(x, y)
    return (pred - np.min(pred)) + 1e-3
    # 将适应度的值>0


def select(pop, fitness):    # nature selection wrt pop's fitness
    pop_index = np.random.choice(np.arange(POP_SIZE), size=POP_SIZE, replace=True, p=(fitness)/(fitness.sum()))
    return pop[pop_index]


# 原作者是随机确定一点，从这点到最后一段即为交叉的DNA片段，
# 但每次最后一段DNA都会进行交叉，丧失了部分随机性。莫烦Python的处理很高超，
# 使用了布尔值找出随机点的位置对DNA进行交叉操作，但每次都有一半的DNA片段进行交叉操作。
# 笔者的思路是对DNA进行剪接，找到两处随机点位置，对中间的片段进行交叉操作。
def crossover(pop):
    cro_pop_index = []
    for cro_index, father in zip(range(POP_SIZE), pop):  # 遍历种群中的每一个个体，将该个体作为父亲
        child = father  # 孩子先得到父亲的全部基因
        # 判断是否对这个个体进行交叉操作
        if np.random.rand() < CROSSOVER_RATE:  # 产生子代时不是必然发生交叉，而是以一定的概率发生交叉
            mother = pop[np.random.randint(POP_SIZE)]  # 在种群中选择另一个个体，并将该个体作为母亲
            cross_point1 = np.random.randint(low=0, high=DNA_SIZE*2)  # 随机产生的交叉点1
            cross_point2 = np.random.randint(low=0, high=DNA_SIZE*2)  # 随机产生的交叉点2
            cross_start = min(cross_point1, cross_point2)
            cross_end = max(cross_point1, cross_point2)
            child[cross_start:cross_end] = mother[cross_start:cross_end]
            # 记录进行交叉操作的个体位置
            cro_pop_index.append(cro_index)
            pop[cro_index] = child
    return pop, cro_pop_index


def mutation(child):
    for point_mut in range(DNA_SIZE*2):
        if np.random.rand() < MUTATION_RATE:  # 以MUTATION_RATE的概率进行变异
            # 方案1
            child[point_mut] = 1 if child[point_mut] == 0 else 0
            # 方案2
            # child[point_mut] = child[point_mut] ^ 1
            # 将变异点的二进制为反转 child[] ^ 1    0变 1，1变 0。
    return child


# plt.ico() plt.ioff()将画图模式改为交互模式，程序遇到plt.show不会暂停，而是继续执行
plt.ion()
# 设定初始sca对象，以免报错；无实际作用
sca = plt.scatter(1, 1)
# 构建初始解
pop = np.random.randint(0, 2, size=(POP_SIZE, DNA_SIZE*2))
# 遗传算法迭代N代
for N in range(N_GENERATIONS):
    x, y = translateDNA(pop)
    # 画图展示每一次结果
    if 'sca' in globals():
        sca.remove()
    fig = plt.figure(dpi=300)
    # sca = ax.scatter(x, y, F(x, y), c='red', marker='o', s=15)
    ax = Axes3D(fig)
    sca = ax.scatter(x, y, F(x, y), c='red', marker='o', s=10)
    # ax2=ax1.twin()
    ax.plot_wireframe(X, Y, Z, rstride=10, cstride=10, color='gray')
    # plt.gca().view_init(elev, azim) 设置网格线密度
    plt.gca().view_init(10, 5)
    plt.show()
    plt.pause(0.01)
    pop, cro_pop_index = crossover(pop)  # 交叉
    for mut_index in cro_pop_index:  # 变异
        child = pop[mut_index]
        child = mutation(child)
        pop[mut_index] = child
    fitness = get_fitness(pop)
    pop = select(pop, fitness)  # 选择生成新的种群
plt.ioff()


max_fitness_index = np.argmax(fitness)
# print("max_fitness:", fitness[max_fitness_index])
print("最优的基因型：", pop[max_fitness_index])
print("x={:.2f}, y={:.2f}".format(x[max_fitness_index], y[max_fitness_index]))
print("F(x, y)= {:.3f}".format(F(x[max_fitness_index], y[max_fitness_index])))
test = str(F(x[max_fitness_index], y[max_fitness_index]))
with open(r'H:\遗传算法.txt', 'w') as f:
    f.write(test)

---