智能优化算法:差分进化+蚁群

Posted by : on

Category : 优化


一、差分进化算法(Differential Evolution,DE)

差分进化算法种群中的每个个体对应一个解向量。它的进化流程则与遗传算法非常类似,都包括变异、杂交和选择操作,但这些操作的具体定义与遗传算法有所不同。差分进化算法被称为一种非常高效的全局优化器,因此在同等条件下,差分进化算法的收敛速度、精度能够比遗传算法、粒子群算法好。

1. 算法步骤

步骤1:初始化参数

在差分进化算法中,需要初始化以下参数:

  • F:变异参数,大多数的情况都是取0到1,主要可调超参数
  • prob_mut:交叉概率,大多数的情况都是取0到1,主要可调超参数
  • n_dim:个体维度
  • size_pop:种群数,通常是Np=5n_dim~10n_dim
  • max_iter:迭代轮数
  • lb, ub:决策空间的上下界,用于初始化种群

初始化种群:

  • 利用lb,ub作为上下界,构建一般符合均匀概率分布,尺寸为size_pop x n_dim的矩阵。

例子: 对于三元优化问题min f(x1, x2, x3) = x1^2 + x2^2 + x3^2 约束为x1*x2 >= 1, x1*x2 <= 5, x2 + x3 = 1, 0 <= x1, x2, x3 <= 5 由于三个变量,n_dim为3,可以选择size_pop为50。根据约束条件,lb为0 ub为5

步骤2:种群变异

对其进行变异产生变异向量,为后期的产生子代种群建立基础。一般的变异公式如下,把种群中两个个体之间的加权差向量加到第三个个体上来产生变异向量 \(V_{i}\left(g+1\right)=X_{r1}\left(g\right)+F\left(X_{r2}\left(g\right)-X_{r3}\left(g\right)\right)\) 根据上下界约束,将更新后的种群V中超过上下界限制的数替换为随机的范围内的样本

步骤3:种群交叉

将变异向量的参数与另外预先决定的目标向量的参数按照一定的规则混合起来产生子种群

\[\left.U_{i,j}\left(g+1\right)=\left\{\begin{matrix}V_{i,j}\left(g+1\right)&ifrand\left(0,1\right)\leq CR\\x_{i,j}\left(g\right)&otherwise\end{matrix}\right.\right.\]

步骤4:最优种群选择

因为要优化函数f的最小值$Y_{min}$,因此选择交叉后的子种群U和原始种群X中更小的更新种群X,这种更新是无约束优化问题,如果有约束,需要考虑约束条件,可以将不满足约束条件的个数乘一个大的权重加到函数上使得最终最小的结果都是满足约束的。

\[\left.X_i\left(g+1\right)=\left\{\begin{matrix}U_i\left(g+1\right)&iff\left(U_i\left(g+1\right)\right)\leq f\left(X_i\left(g\right)\right)\\X_{i}\left(g\right)\end{matrix}\right.\right.\]

本轮迭代输出的最小的Y对应的输入X和Y进行保存

循环上述步骤得到最优的函数解。

主要调整超参数: F:变异参数,大多数的情况都是取0到1 prob_mut:交叉概率,大多数的情况都是取0到1

2. 原理分析

启发式搜索算法是在保留一定随机性的条件下如何更快的搜索到优化函数的最优解,对于差分进化算法的核心公式就是种群变异过程中的变异公式使用了差分的方法。

从几何的角度理解,两项做差代表一个变化向量,即移动的方向,加在一个点上代表将该点移动。因此变异过程实际是一个对原始随机数据点的移动过程。通过这种移动去寻找在约束条件下的最优解。种群交叉即以一定概率保留移动后的点和移动前的点。

差分进化这种方法本质上就是随机寻找点,进行试探目标函数,个人认为理论性不强,但是实验证明这种搜索效果非常好。模型的调参主要针对F:变异参数和prob_mut:交叉概率。

其优势为:

  • 结构简单,容易使用。主要通过差分变异算子来进行遗传操作。
  • 性能优越。具有较好的可靠性,鲁棒性,高效性。
  • 自适应性。差分变异算子可以是一个常数,也可以是一个具有变异步长和搜索方向的自适应能力。
  • 具有内在并行性,可协同搜索。在同一要求下,差分进化算法拥有更快的收敛速度

3. 代码

class DE(GeneticAlgorithmBase):
    def __init__(self, func, n_dim, F=0.5,
                 size_pop=50, max_iter=200, prob_mut=0.3,
                 lb=-1, ub=1,
                 constraint_eq=tuple(), constraint_ueq=tuple(),
                 n_processes=0):
        super().__init__(func, n_dim, size_pop, max_iter, prob_mut,
                         constraint_eq=constraint_eq, constraint_ueq=constraint_ueq, n_processes=n_processes)

        self.F = F
        self.V, self.U = None, None
        self.lb, self.ub = np.array(lb) * np.ones(self.n_dim), np.array(ub) * np.ones(self.n_dim)
        self.crtbp()

    def crtbp(self):
        # 初始化种群
        self.X = np.random.uniform(low=self.lb, high=self.ub, size=(self.size_pop, self.n_dim))
        return self.X

    def chrom2x(self, Chrom):
        pass

    def ranking(self):
        pass

    def mutation(self):
        '''
        V[i]=X[r1]+F(X[r2]-X[r3]),
        where i, r1, r2, r3 are randomly generated
        '''
        X = self.X
        # i is not needed,
        random_idx = np.random.randint(0, self.size_pop, size=(self.size_pop, 3))

        r1, r2, r3 = random_idx[:, 0], random_idx[:, 1], random_idx[:, 2]

        # 这里F用固定值,为了防止早熟,可以换成自适应值
        self.V = X[r1, :] + self.F * (X[r2, :] - X[r3, :])

        # the lower & upper bound still works in mutation
        mask = np.random.uniform(low=self.lb, high=self.ub, size=(self.size_pop, self.n_dim))
        self.V = np.where(self.V < self.lb, mask, self.V)
        self.V = np.where(self.V > self.ub, mask, self.V)
        return self.V

    def crossover(self):
        '''
        if rand < prob_crossover, use V, else use X
        '''
        mask = np.random.rand(self.size_pop, self.n_dim) < self.prob_mut
        self.U = np.where(mask, self.V, self.X)
        return self.U

    def selection(self):
        '''
        greedy selection
        '''
        X = self.X.copy()
        f_X = self.x2y().copy()
        self.X = U = self.U
        f_U = self.x2y()

        self.X = np.where((f_X < f_U).reshape(-1, 1), X, U)
        return self.X

    def run(self, max_iter=None):
        self.max_iter = max_iter or self.max_iter
        for i in range(self.max_iter):
            self.mutation()
            self.crossover()
            self.selection()

            # record the best ones
            generation_best_index = self.Y.argmin()
            self.generation_best_X.append(self.X[generation_best_index, :].copy())
            self.generation_best_Y.append(self.Y[generation_best_index])
            self.all_history_Y.append(self.Y)

        global_best_index = np.array(self.generation_best_Y).argmin()
        self.best_x = self.generation_best_X[global_best_index]
        self.best_y = self.func(np.array([self.best_x]))
        return self.best_x, self.best_y

二、蚁群算法(ACA, Ant Colony Algorithm)

蚁群算法是对自然界蚂蚁的寻径方式进行模拟而得出的一种仿生算法。蚂蚁在运动过程中会留下信息素,蚂蚁在运动时又可以感知这种物质,来确定自己的运动方向。因此,大量蚂蚁组成的蚂蚁群体的集体行为表现出一种信息的正反馈现象。蚁群算法最常用来解决TSP的最短路径优化问题。

  • 时间上来说,路径短释放的信息素累加更快
  • 浓度上来说:路径短释放信息素浓度高

1. 算法步骤

步骤1:初始化参数

  • size_pop:蚂蚁数量
  • max_iter:迭代轮数,建议200
  • n_dim:样本维度(TSP问题中城市)
  • \(\alpha\):信息素因子,其值越大,表示信息素的浓度在转移中起的作用越大,主要可调超参数一般[1, 4]
  • \(\beta\):启发函数因子,其值越大,表示启发函数在转移中起的作用越大,主要可调超参数一般[0, 5]
  • \(\eta\):启发函数矩阵(TSP问题中定义为两个城市间距离的倒数)
  • \(\tau\):信息素浓度矩阵
  • \(\rho\):信息素挥发速度,更新信息素浓度矩阵,主要可调超参数,一般[0.2, 0.5]
  • Q:常数,表示蚂蚁循环一次所释放的信息素总量,主要可调超参数,一般[10, 1000]

步骤2:计算转移概率

根据上一轮迭代更新的信息素计算状态转移概率,根据状态转移概率获取蚂蚁的轨迹。

\[P_{ij}^k=\begin{cases}\frac{[\tau_{ij}(t)]^\alpha\times[\eta_{ij}(t)]^\beta}{\sum_{s\in allow_k}[\tau_{is}(t)]^\alpha\times[\eta_{is}(t)]^\beta},s\in allow_k\\0,s\notin allowk&\end{cases}\]

获取轨迹类似于贪心算法的过程,设置初始节点后根据概率,在可取结点范围内归一化后,取出下一个节点的标号。获得轨迹后计算该轨迹的行驶总距离。在所有蚂蚁中可以获取最优的距离作为该轮迭代的最优结果输出。

步骤3:更新信息素

\[\begin{cases}\tau_{ij}(t+1)=(1-\rho)\tau_{ij}(t)+\Delta\tau_{ij}\\\Delta\tau_{ij}=\sum_{k=1}^m\Delta\tau_{ij}^k\end{cases},0<\rho<1\]

计算信息素变化量\(\Delta\tau_{ij}\),共有三种方法,对应三种不同蚁群模型,三种模型都是只更新该条路径对应位置的信息素

1. Ant-Cycle Mode

\(L_{k}\)代表该路径的总长度,假设每一只蚂蚁的信息素释放速度是相同的,所以总量Q除以长度代表每一段的更新 \(\Delta\tau_{ij}^k=\begin{cases}\frac{Q}{L_k}&蚂蚁k从i到j\\0&其他\end{cases}\)

2. Ant-Quantity Model

\(d_{ij}\)表示城市i到j的距离,假设每两个坐标之间的信息素总量一定,信息素更新与城市之间距离成反比。 \(\Delta\tau_{ij}^k=\begin{cases}\frac{Q}{d_{ij}}&蚂蚁k从i到j\\0&其他\end{cases}\)

3. Ant-Density Model

直接用信息素释放量进行更新,代表所有蚂蚁的信息素释放量一致,与路径无关。

\(\Delta\tau_{ij}^k=\begin{cases}{Q}&蚂蚁k从i到j\\0&其他\end{cases}\) 更新后信息素用于下一轮迭代的状态转移概率计算

2. 原理分析

蚁群算法核心就是状态转移概率矩阵的建立使用了信息素浓度矩阵(不断更新)和启发函数矩阵(固定,描述模型的基线关系),更新信息素矩阵的方法有三种,一般使用前两种进行更新。

启发函数矩阵确定了一个还不错的初始化结果,信息素矩阵不断学习得到一个启发式的好结果,相比于差分优化算法的差分枚举寻优,蚁群算法是需要学习并不断更新模型的,所以效率可能不如差分优化,收敛速度慢,容易陷入局部优化问题。

3. 代码

class ACA_TSP:
    def __init__(self, func, n_dim,
                 size_pop=10, max_iter=20,
                 distance_matrix=None,
                 alpha=1, beta=2, rho=0.1,
                 ):
        self.func = func
        self.n_dim = n_dim  # 城市数量
        self.size_pop = size_pop  # 蚂蚁数量
        self.max_iter = max_iter  # 迭代次数
        self.alpha = alpha  # 信息素重要程度
        self.beta = beta  # 适应度的重要程度
        self.rho = rho  # 信息素挥发速度

        self.prob_matrix_distance = 1 / (distance_matrix + 1e-10 * np.eye(n_dim, n_dim))  # 避免除零错误

        self.Tau = np.ones((n_dim, n_dim))  # 信息素矩阵,每次迭代都会更新
        self.Table = np.zeros((size_pop, n_dim)).astype(int)  # 某一代每个蚂蚁的爬行路径
        self.y = None  # 某一代每个蚂蚁的爬行总距离
        self.generation_best_X, self.generation_best_Y = [], []  # 记录各代的最佳情况
        self.x_best_history, self.y_best_history = self.generation_best_X, self.generation_best_Y  # 历史原因,为了保持统一
        self.best_x, self.best_y = None, None

    def run(self, max_iter=None):
        self.max_iter = max_iter or self.max_iter
        for i in range(self.max_iter):  # 对每次迭代
            prob_matrix = (self.Tau ** self.alpha) * (self.prob_matrix_distance) ** self.beta  # 转移概率,无须归一化。
            for j in range(self.size_pop):  # 对每个蚂蚁
                self.Table[j, 0] = 0  # start point,其实可以随机,但没什么区别
                for k in range(self.n_dim - 1):  # 蚂蚁到达的每个节点
                    taboo_set = set(self.Table[j, :k + 1])  # 已经经过的点和当前点,不能再次经过
                    allow_list = list(set(range(self.n_dim)) - taboo_set)  # 在这些点中做选择
                    prob = prob_matrix[self.Table[j, k], allow_list]
                    prob = prob / prob.sum()  # 概率归一化
                    next_point = np.random.choice(allow_list, size=1, p=prob)[0]
                    self.Table[j, k + 1] = next_point

            # 计算距离
            y = np.array([self.func(i) for i in self.Table])

            # 顺便记录历史最好情况
            index_best = y.argmin()
            x_best, y_best = self.Table[index_best, :].copy(), y[index_best].copy()
            self.generation_best_X.append(x_best)
            self.generation_best_Y.append(y_best)

            # 计算需要新涂抹的信息素
            delta_tau = np.zeros((self.n_dim, self.n_dim))
            for j in range(self.size_pop):  # 每个蚂蚁
                for k in range(self.n_dim - 1):  # 每个节点
                    n1, n2 = self.Table[j, k], self.Table[j, k + 1]  # 蚂蚁从n1节点爬到n2节点
                    delta_tau[n1, n2] += 1 / y[j]  # 涂抹的信息素
                    # 这里使用的事Ant Quantity System模型, Q=1
                n1, n2 = self.Table[j, self.n_dim - 1], self.Table[j, 0]  # 蚂蚁从最后一个节点爬回到第一个节点
                delta_tau[n1, n2] += 1 / y[j]  # 涂抹信息素

            # 信息素飘散+信息素涂抹
            self.Tau = (1 - self.rho) * self.Tau + delta_tau

        best_generation = np.array(self.generation_best_Y).argmin()
        self.best_x = self.generation_best_X[best_generation]
        self.best_y = self.generation_best_Y[best_generation]
        return self.best_x, self.best_y

About Zifu Zhang
Zifu Zhang

I am a master of EE in beihang University, interested in deep learning, image compression and AIGC

Email : zifuzhang@buaa.edu.cn

Website : https://bblgbr.github.io

About Zifu Zhang

Hi, my name is Zifu Zhang.

Star
Categories
Useful Links