染色体遗传模型的matlab程序,一个用matlab实现的50行的实数染色体遗传算法程序 - 第 2 页 - 计算模拟 - 小木虫 - 学术 科研 互动社区...

【本文属作者原创，但已发表于科学网(链接地址：http://blog.sciencenet.cn/blog-3102863-1029280.html)，现稍作格式上的修该后转载，并发金币祝大家新年快乐！】

1. 引言

遗传算法 (genetic algorithms) 是一种很有意思最优化方法，常用于解决一些传统方法力所不及的多变量最优化问题。这种方法很通用，即用同样的思想可以解决很多不同的问题。只要你能对问题所有可能的解定义一个反映其好坏程度的量，就有可能用遗传算法找到高质量的解。遗传算法的应用十分广泛，编程实现也很容易。本文介绍遗传算法的基本概念和算法流程，给出一个用matlab实现的50行的遗传算法程序，并用一个简单的例子展示它的应用。本文给出的程序适用于待优化函数的各个变量皆为实数的情形，可用于多变量拟合等问题。读者若弄明白了该程序，应该可以写出函数变量取整数值(比如二进制数值)的遗传算法程序，甚至写出整数与实数混合情形的程序。

2 遗传算法的基本概念

2.1 生物学中的相关概念

与遗传算法有关的生物学概念主要有：

a)  染色体(chromosome)。

所有生物都由细胞组成，每一个细胞中都有一套相同的染色体。一条染色体由若干基因(gene) 组成，每个基因控制一种特定的蛋白质，从而决定生物的某种特征。所有染色体合称为基因组(genome)。基因组完全决定了一个生物个体。该个体在微观(基因)层次的表现称为基因型 (genotype)，在宏观(特征)层次的表现称为显型 (phenotype)。在简单的遗传算法中，将基因组中的若干条染色体看作一整条染色体。

b)  个体复制。

在复制的过程中，父母的染色体通过交叉(crossover)产生子女的染色体。染色体还可以以一定的小概率变异(mutate)。

2.2 遗传算法的基本概念

在一般的遗传算法应用中，我们研究的是这样的优化问题：对一个有 N 个参数 x_i (1 <= i <= N) 的实函数 f(x_i) (该函数可以有任何形式；我们不要求它有任何解析性质)，找到一组合适的参数使得该函数的值尽可能地大或者尽可能地小。该函数叫做适应度函数 (fitness function)，其值叫做适应度。

下面我们用这个最优化问题结合上面的生物学概念说明遗传算法中的基本概念：

a)  个体，问题的一个可能的解，即任何一组参数。这组参数也叫基因型，对应的函数值叫做显型。

b)  解空间，又叫搜索空间，即所有参数构成的 N 维空间。

c)  种群 (population)，即一定数量的个体的集合。该数量叫做种群规模 (population size)。 遗传算法的每一次循环中的种群称为一代 (generation)。每一代中最好的解叫做精英 (elite)。

3. 遗传算法的基本流程

首先，随机产生一个初始种群。

遗传演化开始

计算当前种群中所有个体的适应度

根据适应度的大小选择一部分较好的个体，使其存活下来并繁殖子女(所谓“适者生存”)

用遗传操作(genetic operators) 繁殖子女(子女数目等于死去的个体的数目)，包括

交叉：将两个候选者的染色体以某种方式交叉，产生两个新的个体的染色体

变异：以某种方式随机改变一些个体的染色体的部分基因

在每一代中，精英都会记录下来，且不参与变异

若找到了足够好的解，则可提前终止演化

遗传演化结束

4. 一个50行的matlab函数

下面直接给出一个实现简单的遗传算法的matlab函数。此函数虽短小，但却很通用。它适用于任何实数基因的情形。这里的实数是与整数(包括二进制数)对应的。在很多遗传算法的应用中，部分或者全部变量的值仅能取整数。我给出的函数不适用于这种离散变量的情形，但通过少量的修改可以变成适用于该情形的函数。若读者对这样的应用感兴趣，可以自行修改。我的程序中有大量的注释，可以帮助理解。

有两点要注意：

1)该函数假定我们的目标是求一个函数的最小值。若我们目标是求最大值，则可将函数f(x)换成-f(x)。

2) 该函数假定每个变量的取值范围都是[0, 1]。也就是说，矩阵population的每个元素的值都在此区间。在编写适应度函数时要根据具体问题对变量进行变换。

==================================================================

function [best_fitness, elite, generation] = my_ga(number_of_variables, fitness_function, ...

population_size, parent_number, mutation_rate, maximal_generation, minimal_cost)

% number_of_variables = 求解问题的参数个数

% fitness_function = 自定义适应度函数名

% population_size = 种群规模(每一代个体数目)

% parent_number = 每一代中保持不变的数目(除了变异)

% mutation_rate = 变异概率

% maximal_generation = 最大演化代数

% minimal_cost = 最小目标值(函数值越小，则适应度越高)

cumulative_probabilities = cumsum((parent_number:-1:1) / sum(parent_number:-1:1)); %一个从0到1增长得越来越慢的函数

best_fitness = ones(maximal_generation, 1); % 用于记录每一代的最佳适应度

elite = zeros(maximal_generation, number_of_variables); % 用于记录每一代的最优解

child_number = population_size - parent_number; % 每一代子女的数目

population = rand(population_size, number_of_variables); % 初始化种群

for generation = 1 : maximal_generation % 演化循环开始

cost = feval(fitness_function, population); % 计算所有个体的适应度

[cost, index] = sort(cost); % 将适应度函数值从小到大排序

population = population(index(1:parent_number), ; % 先保留一部分较优的个体

best_fitness(generation) = cost(1); % 记录本代的最佳适应度

elite(generation,  = population(1, ; % 记录本代的最优解(精英)

if best_fitness(generation) < minimal_cost; break; end % 若最优解已足够好，则停止演化

for child = 1:2:child_number % 染色体交叉开始

mother = min(find(cumulative_probabilities > rand)); % 选择一个较优秀的母亲

father = min(find(cumulative_probabilities > rand)); % 选择一个较优秀的父亲

crossover_point = ceil(rand*number_of_variables); % 随机地确定一个染色体交叉点

mask1 = [ones(1, crossover_point), zeros(1, number_of_variables - crossover_point)];

mask2 = not(mask1);

mother_1 = mask1 .* population(mother, ; % 母亲染色体的前部分

mother_2 = mask2 .* population(mother, ; % 母亲染色体的后部分

father_1 = mask1 .* population(father, ; % 父亲染色体的前部分

father_2 = mask2 .* population(father, ; % 父亲染色体的后部分

population(parent_number + child,  = mother_1 + father_2; % 一个孩子

population(parent_number+child+1,  = mother_2 + father_1; % 另一个孩子

end % 染色体交叉结束

% 染色体变异开始

mutation_population = population(2:population_size, ; % 精英不参与变异

number_of_elements = (population_size - 1) * number_of_variables; % 全部基因数目

number_of_mutations = ceil(number_of_elements * mutation_rate); % 将要变异的基因数目

mutation_points = ceil(number_of_elements * rand(1, number_of_mutations)); % 确定要变异的基因

mutation_population(mutation_points) = rand(1, number_of_mutations); % 对选中的基因进行变异操作

population(2:population_size,  = mutation_population; % 发生变异之后的种群

% 染色体变异结束

end % 演化循环结束

==================================================================

5. 一个简单的例子

考虑一个浅显的例子：计算函数 f(x) = \sum_{i=1}^{10}  x_i^2 的最小值，其中每个变量的取值区间都是[-1, +1]。我们知道该问题的最优解：每个x_i都等于0。如果用遗传算法研究这个问题，我们可以写出如下适应度函数(注意matlab矢量化编程的技巧)：

%请将该函数放在名为my_fitness的文件内

function y = my_fitness(population)

population = 2 * population - 1;

y = sum(population.^2, 2);

写好了适应度函数，再写一个脚本调用遗传算法函数并画图即可完成任务。脚本内容如下：

clear; close all;

% 调用 my_ga 进行计算

[best_fitness, elite, generation] = my_ga(10, 'my_fitness', 100, 50, 0.1, 10000, 1.0e-6);

% 最佳适应度的演化情况

figure

loglog(1 : generation, best_fitness(1 : generation), 'linewidth',2)

xlabel('Generation','fontsize',15);

ylabel('Best Fitness','fontsize',15);

set(gca,'fontsize',15,'ticklength',get(gca,'ticklength')*2);

% 最优解的演化情况

figure

semilogx(1 : generation, 2 * elite(1 : generation,  - 1)

xlabel('Generation','fontsize',15);

ylabel('Best Solution','fontsize',15);

set(gca,'fontsize',15,'ticklength',get(gca,'ticklength')*2);

此脚本的计算只需一秒钟左右。计算结束后，会出现两幅图。第一幅图展示最佳适应度(对于最小化问题，函数值越小，适应度越高)随演化代数增加而变化的情况。第二幅图展示最优解变量随演化代数增加而变化的情况。