Python 遗传算法 Genetic Algorithm

粒子群算法常常用于在连续解空间内搜索，而在不连续、离散的空间内常常会出现搜索越界的问题

例如旅行商问题，寻找可以遍历 15 个地点的最短路径（当然可以用二进制状态压缩 + 动态规划解决），以 {0, 1, ..., 14} 表示这些地点，并以 {0, 1, ..., 14} 的一种排列方式为一个解

当这个问题的解集映射在 15 维空间中时，这个空间中的可行解将非常的稀疏，从而阻碍粒子群的搜索

遗传算法有几个关键词：

• 保优：当种群更新时，不改变最优的几个个体，为交叉提供优质基因并与新个体进行比较
• 天择：根据每个个体的适应度，使用轮盘赌法进行挑选
• 交叉：对天择生成的每个个体，以一定概率与原来的个体进行交叉
• 变异：对天择生成的每个个体，以一定概率进行基因突变

下面是我编写的遗传算法模板，在使用时需要重写 new_unit（群体初始化方法）、cross（两个体交叉方法）、variation（个体变异方法）、fitness（群体适应度计算方法） 函数（后面我将会以旅行商问题进行举例）

其中的 fit 方法为主函数，记群体规模为 n，for 循环体的内容为：

• 对当前的群体进行重叠检测，去除重复的个体
• 计算每个个体的适应度，排序后对适应度前 5% 的个体进行“保优”，得到规模 0.05n 的新群体
• 对原有群体进行天择，选出 0.95n 的个体（可重复），根据概率对每个个体进行交叉、变异操作，加入新群体得到规模 n 的群体

import random
 
import numpy as np
from tqdm import trange
 
DTYPE = np.float16
 
 
class GeneticOpt:
    ''' 遗传算法
        :param n_unit: 染色体群体规模
        :param n_gene: 染色体的基因数
        :param well_radio: 最优个体比例
        :param cross_proba: 交叉概率
        :param var_proba: 变异概率'''
 
    def __init__(self,
                 n_unit: int,
                 n_gene: int,
                 well_radio: float = 0.05,
                 cross_proba: float = 0.4,
                 var_proba: float = 0.3):
        self.n_unit = n_unit
        self.n_gene = n_gene
        self._well_radio = well_radio
        self._cross_proba = cross_proba
        self._var_proba = var_proba
 
    def new_unit(self, size: int) -> np.ndarray:
        ''' 初始化染色体群体
            :return: [size, n_gene]'''
        raise NotImplementedError
 
    def cross(self, unit: np.ndarray, other: np.ndarray) -> np.ndarray:
        ''' 交叉遗传
            :return: [n_gene, ]'''
        raise NotImplementedError
 
    def variation(self, unit: np.ndarray) -> np.ndarray:
        ''' 基因突变
            :return: [n_gene, ]'''
        gene_idx = np.arange(self.n_gene)
        l = random.choice(gene_idx)
        r = random.choice(gene_idx[l:])
        np.random.shuffle(unit[l: r + 1])
        return unit
 
    def fitness(self, group: np.ndarray) -> np.ndarray:
        ''' 适应度函数 (max -> best)'''
        raise NotImplementedError
 
    def fit(self, epochs: int,
            patience: int = np.inf,
            prefix: str = 'GA_fit') -> np.ndarray:
        ''' :param epochs: 训练轮次
            :param patience: 允许搜索无进展的次数'''
        cur_group = self.new_unit(self.n_unit)
        pbar = trange(epochs)
        last_fitness, angry = - np.inf, 0
        # 最优个体数, 随机选取数
        n_well = max(2, round(self.n_unit * self._well_radio))
        n_choose = self.n_unit - n_well
        for _ in pbar:
            cur_group = np.unique(cur_group, axis=0)
            # 计算每个个体的适应度并排序
            fitness = self.fitness(cur_group)
            order = np.argsort(fitness)[::-1]
            cur_group, fitness = cur_group[order], fitness[order]
            # 收敛检测
            angry = 0 if fitness[0] > last_fitness else angry + 1
            last_fitness = fitness[0]
            if angry == patience: break
            # 保留一定数量的个体
            new_group = cur_group[:n_well]
            pbar.set_description((f'%-10s' + '%-10.4g') % (prefix, fitness[0]))
            # 使用轮盘赌法进行筛选
            proba = fitness - fitness.min()
            proba = proba / proba.sum()
            for pc, pv in np.random.random([n_choose, 2]):
                unit = random.choices(cur_group, weights=proba)[0]
                # 交叉遗传 / 基因突变
                if pc <= self._cross_proba:
                    unit = self.cross(unit.copy(), random.choices(cur_group, weights=proba)[0].copy())
                if pv <= self._var_proba:
                    unit = self.variation(unit.copy())
                # 拼接新个体
                new_group = np.concatenate([new_group, unit[None]])
            cur_group = new_group
        return cur_group[0]

求解示例

对于 15 个地点的旅行商问题，重写的函数思路如下：

• init_adj：使用实例属性 pos 记录 15 个地点的位置，实例属性 adj 记录这 15 个地点的邻接矩阵
• new_unit：生成 n 个 [0, 1, ..., 14]，使用 np.random.shuffle 进行打乱
• fitness：对每个个体，依次遍历个体中的地点叠加距离（越小表示该解越优），并取负值（越大表示该解越优，符合 fit 函数的设计）
• cross：因为旅行商问题中的解在进行交叉时（交换片段），容易出现“重复经过一地点”的情况，故此处不使用交叉
• variation：随机选取区间的左右边界，使用 np.random.shuffle 对该区间的基因进行打乱（已编写在模板中）

if __name__ == '__main__':
    import matplotlib.pyplot as plt
 
 
    class ShortestPath(GeneticOpt):
 
        def init_adj(self):
            # 初始化邻接矩阵
            self.pos = np.random.random([self.n_gene, 2]) * 10
            self.adj = np.zeros([self.n_gene] * 2, dtype=DTYPE)
            for i in range(self.n_gene):
                for j in range(i + 1, self.n_gene):
                    self.adj[i][j] = self.adj[j][i] = \
                        np.sqrt(((self.pos[i] - self.pos[j]) ** 2).sum())
 
        def new_unit(self, size):
            ''' 初始化染色体群体'''
            group = []
            for _ in range(size):
                unit = list(range(self.n_gene))
                np.random.shuffle(unit)
                group += [unit]
            return np.array(group, dtype=np.int32)
 
        def fitness(self, group):
            ''' 适应度函数 (max -> best)'''
            group = np.concatenate([group, group[:, :1]], axis=-1)
            return - self.adj[group[:, :-1], group[:, 1:]].sum(axis=-1)
 
 
    ga = ShortestPath(80, 15, cross_proba=0, var_proba=0.6)
    ga.init_adj()
    unit = ga.fit(500)
 
    # 绘制最优路径
    fig = plt.subplot()
    for key in 'right', 'top':
        fig.spines[key].set_color('None')
    plt.plot(*ga.pos[unit].T, c='deepskyblue')
    plt.scatter(*ga.pos.T, marker='p', c='orange')
    plt.show()