遗传算法

一、介绍

? 遗传算法是用于解决最优化问题的一种搜索算法。遗传算法借用了生物学里达尔文的进化理论：“适者生存，优胜劣汰”，将该理论以算法的形式表现出来就是遗传算法的过程。

二、基本原理

1）程序流程

? 遗传算法是通过大量备选解的变换、迭代和变异，在解空间中并行动态地进行全局搜索的最优化方法，是模拟生物基因遗传的做法，算法的基本原理通过编码组成的群体组成初始群体后，遗传操作的任务就算对群体的个体按照他们对环境适应度（适应度）评估，施加了一定的操作之后从而实现优胜劣汰的进化过程，算法的基本程序框图：

enter

2）物竞天择

? 我们可以假设有一个种群，这个种群存在的目标是帮我寻找一个函数F（x）的极值（可能是max or min），首先我们随机生成一组种群，假设我对变量x的取值范围为【0，8】，这是我们的变量域，因为我们要玩一些不一样的操作，先跟着我，后面你会恍然大悟！我们会有一个二进制域【0，2^8】,我们的二进制域对应的就是我们的变量域，我们将数字8看成8位的二进制数，可以得到2的八次方为：256，但因为我们在二进制当中我们是从0开始的，0-255，有256个数,我们的八位二进制数最多只能取到255,取不到256。

"二进制" 0 0 0 0   0 0 0 0
"二进制数 x-> " 1 1 1 1  1 1 1 1
"x转化为十进制" 1*2^0 + 1*2^1 + 1*2^2 + 1*2^3 + 1*2^4 + 1*2^5 + 1*2^6 + 1*2^7 = 255

? 因为我们实数x的范围只有【0，8】，我们将二进制转化为十进制之后，还有求解二进制域到变量域值，二进制域到变量域之间的映射关系为：

"变量域x_"  x_ = x*8/2^8

? 如此，我们将我们的8看作是染色体的长度，这样我们就可以做一些基因重组和基因变异的操作了！！！，我们先来一组基因重组图：

在这里插入图片描述

? 是不是很熟悉，这个不就是高中生物的杂交实验吗？如果我们把DD和dd看成我们染色体的0和1呢？我们再形象一点

# 定义我们的染色体
x1 = [0 1 0 0 1 0 0 1]
x2 = [1 1 0 1 0 0 0 1]
#OK 我们现在有两个种群了（抽象一点）
# 第一个种群 x1 的染色体为 0 1 0 0  1 0 0 1
# 第一个种群 x2 的染色体为 1 1 0 1  0 0 0 1

#我们尝试做一下杂交实验？学一学孟德尔的八年抗战
"第一代:" 0 1 0 0    1 0 0 1      *      1 1 0 1    0 0 0 1
                   |                              |
                   v                              v
        
"配子:" [0 1 0 0]    [1 0 0 1 ]        [1 1 0 1 ]   [0 0 0 1]  
														 #用x1的第一个配子和x2的第二个配子进行基因重组
														 #用x1的第二个配子和x2的第一个配子进行基因重组

"F1代:" [0 1 0 0][0 0 0 1]              [1 1 0 1 ][1 0 0 1] 

new_x1 = [0 1 0 0 0 0 0 1]
new_x2 = [1 1 0 1 1 0 0 1]

# 不知道你有没有恍然大悟

? 有了二进制的操作，我们做变异的时候也是同样的道理，我们只需要把0变异为1，1变异为0，就能完成我们的变异算法。

3）算法步骤

? 求解函数:f(x) = 9sin(5x) + 7cos(4*x)的最大值，x的取值范围为【0，8】

3.1初始化种群

? 我们先初始化我们的种群，先初始化一些随机解，因为我们要在搜索空间内寻找我们的最优解，可参考粒子群算法的思想。

def initpop(popsize, chromlength):
    """
    Parameters
    ----------
    popsize :  TYPE
        种群的数目.
    chromlength : TYPE
        表示染色体的长度（二值数的长度）， 长度的大小取决于变量的二进制编码的长度     

    Returns 
    -------
    pop : TYPE
        随机种群
    """
    pop=np.round(np.random.rand(popsize,chromlength))
    return pop

3.2二进制编码转化为十进制数

? 设计我们的算法，将二进制的矩阵转化为十进制数，为了得到我们的二进制域。

def decodebinary(pop):
    """
    二进制数转十进制数函数
    Parameters
    ----------
    pop : TYPE
          初始化种群.

    Returns
    -------
    pop2 种群的染色体二进制数转十进制数

    """
    px,py = pop.shape
    pop1 = np.zeros((px,py))
    for i in range(py):
        pop1[:,i] = (2**(py-i-1))*pop[:,i]
    pop2=np.sum(pop1, 1)#对pop1向量的每行求和 标识位0代表每列求和，标识位1代表每行求和
    
    return pop2

? 我们的pop2是一个转化为十进制数之后的种群，可以理解为种群的表现型，我们的0和1的操作反应的是种群的基因型，因为我们的变量不可能只有一种，扩展我们的思维，我们先设置的染色体长度，因为我们的变量只有一个，它的值域在【0，8】，如果我们有两个变量，不同变量对应不同值域，我们的染色体的长度是否也要切片计算？仔细思考。

def decodechrom(pop,spoint,length):
    """
    将二进制编码转化为十进制数
    Parameters
    ----------
    pop : TYPE
        DESCRIPTION.
    spoint : TYPE
        染色体的起始位.
    length : TYPE
        DESCRIPTION.
    Returns
    -------
    pop2 : TYPE
        DESCRIPTION.

    """
    #对于多个变量而言，如有两个变量，采用20为表示，每个变量10位，则第一个变量从1开始，另一个变量从11开始。本例为一个变量
    #这句话的意思就是加入目标函数需要两个变量，则我可将染色体的数量拆为两个，然后遗传迭代
    #值得注意的是，我的染色体的长度也要跟着变量的数量改变，呈倍数关系
    pop1=pop[:,spoint:spoint+length]
    pop2=decodebinary(pop1)
    
    return pop2

? 这里的spoint代表的是我们的一个终止位，可以对我们的变量进行一个切片，自行理解。

3.3实现目标函数的计算

? 实现目标函数的计算，就是把x的值带进去，很简单。

def calobjvalue(pop):
    """
    实现目标函数的计算
    Parameters
    ----------
    pop : TYPE
        种群.

    Returns
    -------
    objvalue : TYPE
        返回目标函数值.

    """
    temp = decodechrom(pop, 0, 15)
    x=temp*15/32767
    objvalue = 9*np.sin(5*x)+7*np.cos(4*x)
    objvalue = np.reshape(objvalue,(-1,1))#以行的形式输出
    return objvalue

3.4剔除0以外的值

? 计算个体的适应值，在calobjvalue已经计算好，需要将小于0的个体删除，方便后续的概率计算。

def calfitvalue(objvalue):
    """
    计算个体的适应值，在calobjvalue已经计算好，需要将小于0的个体删除，方便后续的概率计算

    Parameters
    ----------
    objvalue : TYPE
        目标函数值.
    Returns
    -------
    fitvalue : TYPE
        个体适应值.

    """
    global Cmin
    Cmin = 0
    fitvalue = np.zeros((objvalue.shape[0],objvalue.shape[1]))
    
    px, py = objvalue.shape
    for i in range(px):
        if objvalue[i,0] + Cmin > 0:
            temp = objvalue[i,0] + Cmin
        else:
            temp = 0
        
        fitvalue[i,0]=temp
    return fitvalue

3.5选择算法

? 选择算法我们采用轮盘赌，决定哪些个体可以进入下一代，用轮盘赌选择复制

def selection(pop, fitvalue):
    """
    选择函数 选择复制，决定哪些个体可以进入下一代
    采用轮盘赌选择

    Parameters
    ----------
    pop : TYPE
        种群的个体.
    fitvalue : TYPE
        个体适应值.

    Returns
    -------
    newpop : dict
        新的种群.

    """
    totalfit = np.sum(fitvalue)
    fitvalue_pro = fitvalue / (totalfit + 0.000001)
    fitvalue_pro_cumnsum = np.cumsum(fitvalue_pro)
    fitvalue_pro_cumnsum = np.reshape(fitvalue_pro_cumnsum,(-1,1))#以行的形式
    
    px, py = pop.shape 
    #轮盘随机概率 从小到大排序
    ms = np.sort(np.random.rand(px,1),0)
    
    fitin = 1 - 1   #pop种群 第几代个体 因为python的下标从0开始。所以是第一代的索引是0 故用1-1
    newin = 1 - 1  #pop种群 第几代个体
    newpop_dict = {}
    #我愿称为[适者生存，优胜劣汰] while循环
    while newin <= px-1:
        if ms[newin, 0] < fitvalue_pro_cumnsum[fitin, 0]:
            newpop_dict[newin] = pop[fitin,:]
            newin += 1
        else:
            fitin += 1
    
    newpop = np.zeros((newin,py))
    for i in range(newin):
        newpop[i,:] = newpop_dict[i]
    
    return newpop

3.6交叉算法（基因重组）

? 物竞天择是遗传学的核心，相应的，选择和变异也是遗传算法的核心，我们对染色体进行一个基因重组的操作，具体代码如下：

def crossover(pop, pc):
    """
    交叉算法 实现基因重组
    Parameters
    ----------
    pop : TYPE
        种群.
    pc : TYPE
        交叉概率.
    Returns
    -------
    newpop : TYPE
        DESCRIPTION.

    """
    px, py = pop.shape 
    newpop = np.zeros((px, py))
    # seletion_litst = [i for i in range(px)]
    # sl = seletion_litst[0:px:2]
    for i in range(px-1):
        #是否能够进行基因重组
        if pc > np.random.rand(1)[0]:
            cpoint = int(np.round(np.random.rand(1)[0] * py) - 1)
            if cpoint == 0:
                cpoint = 1
                
            newpop[i, :][0:cpoint] = pop[i, 0:cpoint]
            newpop[i, :][cpoint+1:py] = pop[i+1, cpoint+1:py]
            
            newpop[i+1, :][0:cpoint] = pop[i+1, 0:cpoint]
            newpop[i+1, :][cpoint+1:py] = pop[i, cpoint+1:py]
            
        else:
            newpop[i,:] = pop[i,:]
            newpop[i+1,:] = pop[i+1,:]
    
    return newpop

3.7变异算法（基因突变）

? 0变1，1变0 你上你也行

def mutation(pop, pm):
    """
    变异算法 实现基因突变

    Parameters
    ----------
    pop : TYPE
        种群.
    pm : TYPE
        变异概率.

    Returns
    -------
    newpop : TYPE
        新的变异种群.

    """
    px, py = pop.shape 
    newpop = np.zeros((px, py))
    
    for i in range(px):
        if pm > np.random.rand(1)[0]:
            mpoint = int(np.round(np.random.rand(1)[0] * py) - 1)
            if mpoint == 0:
                mpoint = 1
            newpop[i,:] = pop[i,:]
            if newpop[i, mpoint] == 0:
                newpop[i, mpoint] = 1

        newpop[i,:] = pop[i,:]
    return newpop

3.7最优计算（基因突变）

? 简单的排序取最大（不要忘记我们粒子群算法里面的全局极值和个体极值）

def best(pop, fitvalue):
    """
    最优的个体及其适应值
    Parameters
    ----------
    pop : TYPE
        种群.
    fitvalue : TYPE
        适应值.
    Returns
    -------
    bestindividual : 最大适应值
        DESCRIPTION.
    bestfit : TYPE
        最大适应值的个体.

    """
    px, py = pop.shape 
    bestindividual = pop[0,:]
    bestfit = fitvalue[0]
    
    for i in range(1,px):
        if fitvalue[i] > bestfit:
            bestindividual = pop[i,:]
            bestfit = fitvalue[i]
    
    return bestindividual,bestfit

3.8主函数

? 不要忘记导入我们的numpy 和matplotlib库，将上面的函数和导入的库，放在主函数的上面，然后运行代码，不过我这里的取值是【0，15】，上面说【0，8】便于读者理解。

import numpy as np
import matplotlib.pyplot as plt

if __name__ == "__main__":
    popsize=30; #群体大小
    chromlength=15; #字符串长度（染色体的长度）
    pc=0.7; #交叉概率 
    pm=0.005 #变异概率
    
    #初始化种群
    pop=initpop(popsize, chromlength)
    poptest=pop
    #开始迭代 
    epoch = 500
    
    x_ = []
    y_ = []
    
    for i in range(epoch):
        objvalue = calobjvalue(pop)     #计算目标函数值
        fitvalue = calfitvalue(objvalue)    #计算群体中每个个体的适应度
        
        newpop = selection(pop, fitvalue)
        newpop1 = crossover(newpop, pc)
        newpop2 = mutation(newpop1, pm)
        
        objvalue = calobjvalue(newpop2)     #计算目标函数值
        fitvalue = calfitvalue(objvalue)    #计算群体中每个个体的适应度
        
        bestindividual,bestfit = best(newpop2, fitvalue) #求出群体中适应值最大的个体及其适应值
        x_.append(bestfit)
        # y_.append(bestfit)
        pop = newpop2
        
        plt.plot(x_, lw=3, label='funcotio_max_value_x')
        plt.show()
        
    reshape_best_infrivedual = np.reshape(bestindividual,(1,-1))#以行的形式
    best_x_value = encode(reshape_best_infrivedual)
    plt.plot(x_, lw=3, label='funcotio_max_value_x')
    plt.legend()
    plt.show()
    
    print("最好的个体是 ",bestindividual)
    print("函数：f(x) = 9sin(5x)+7cos(4x)的极值为 ",bestfit[0]) #最后拟合之后随便取一个值即可
    print("x的取值应为： ",best_x_value[0]) #最后拟合之后随便取一个值即可
    
    # plt.plot(fv2, lw=3, label='fv2')
    # plt.legend()
    # plt.show()