张麟

张麟的博客

他的个人主页  他的博客

PyGEP说用说明

张麟  2010年12月11日 星期六 21:38 | 3469次浏览 | 3条评论

PyGEP的一个说用说明。及一些例子。

PyGEP说用说明
简介
PYGEP包是一个基因表达式编程（gene express programming）的一个开源python代码包
里面最核心的程序就是population.py和chromosome.py，只要读懂这两个程序，整个GEP的思想也在其中了
还是先看一下该方法的框架吧!


从流程图可以看出，基本上的步骤就是
1、生成初始种群（由染色体组成chromosomes），每个染色体表达了问题的一个解。这本质上就是一个测试空间。我们所有的工作就是在该测试空间下通过一系列的操作以寻找出符合我们要求的解。
2、对每个染色体进行解码（Execute），得到我们的一个评估结果。每个染色体对应于一个适应性函数（fitness），和遗传算法（GA）类似。
3、对当前种群的改造，生成新的种群，重复2-3直到得到我们需要的解。

PyGEP结构介绍

一、生成种群

初始化种群时需要如下参数的设置，这些参数将决定如何对现在种群进行改造或者叫遗传。

  1. - exclusion_level: selection pressure (typically 1.5)

  3. - mutation_rate: mutation rate (2 per chromosome)

  5. - inversion_rate: inversion probability (0.1)

  7. - is_transposition_rate: non-root transposition (0.1)

  9. - is_transposition_length: IS transposition lengths (1,2,3)

  11. - ris_transposition_rate: root transposition rate (0.1)

  13. - ris_transposition_length: RIS transposition lengths (1,2,3)

  15. - gene_transposition_rate: gene transposition rate (0.1)

  17. - crossover_one_point_rate: 1-point crossover rate (0.3)

  19. - crossover_two_point_rate: 2-point crossover rate (0.3)

  21. - crossover_gene_rate: full gene crossover (0.1)
复制代码

__init__()函数：def __init__(self, cls, size, head, genes=1, linker=default_linker)
cls表示种群对应的染色体类型，
size表示种群大小，
head种群头部长度，
genes 基因组数，
linker基因链接的方式
这些就是我们调用population的接口。整个改造过程是有cycle()函数完成的。
类中还定义了一个best属性用以标示当前的最好染色体。

  1. best = property(
  2.         # Gives preference to later individuals tied for best
  3.         lambda self: max(reversed(self.population)),
  4.         doc='The best Chromosome of the current generation'
  5.     )
复制代码


二、染色体
问题的描述主要是在染色体中进行的。所以实质上染色体的构造是我们对模型的构造。而种群是对算法的构造。
之所以说GEP算法在并行上有很大的优势是因为她天生就是并行的。因为对于每个染色体的解码是与其他染色体无关的。
分析chromosome类的外部作用，我们需要知道的是如何生成chromosome类，并且需要定义如何评估（fitness及solved）
chromosome类提供了一个generate函数：
def generate(cls, head, genes=1, linker=default_linker)
head是头部长度，genes是基因组数，linker是指链接方式。
同样还有：
def __init__(self, genes, head, linker=default_linker)同上。
关于fitness和solved，PyGEP提供了两个属性来连接。
fitness = property(lambda self: self._fitness(), doc='Fitness value')
solved  = property(lambda self: self._solved(),  doc='Problem solved')
而所有具体问题的描述都将在fitness和solved中进行。
当然在进行如上操作之前，我们必须要设置符号空间    functions和 terminals ，用以表示表达式。
三、一个例子
这个例子是包中的例子（majority），理解一下就清楚如何使用这个程序包了。
首先定义一个类：

  1. from pygep import *
  2. from pygep.functions.linkers import or_linker
  3. from pygep.functions.logical import LOGIC_ALL

  5. # Demo for the boolean majority function

  7. class Majority(object):
  8.     SAMPLE = []
  9.    
  10.     @staticmethod
  11.     def populate():
  12.         for a in 0, 1:
  13.             for b in 0, 1:
  14.                 for c in 0, 1:
  15.                     Majority.SAMPLE.append(Majority(a, b, c))
  16.         
  17.     def __init__(self, a, b, c):
  18.         self.a = a
  19.         self.b = b
  20.         self.c = c
  21.         self.majority = (a and b) or (a and c) or (b and c)
复制代码



该方程实质上是给了Majority这个bool运算函数的一个运算法则。我们想让gep通过符号运算来得到这法则。
在这个地方我所知道的其实就是不同的参数对应的Majority的值。

  1. class BooleanFunction(Chromosome):
  2.     functions = LOGIC_ALL
  3.     terminals = 'a', 'b', 'c'
  4.    
  5.     def _fitness(self):
  6.         # Fitness function: number of hits
  7.         hits = 0
  8.         for i in Majority.SAMPLE:
  9.             if i.majority == bool(self(i)):
  10.                 hits += 1
  11.         
  12.         return hits
  13.    
  14.     def _solved(self):
  15.         return self.fitness >= len(Majority.SAMPLE)
复制代码



这部分继承chromosome类，初始化了functions和terminals的算符表，并重写了_fitness()和_solved()函数。

  1. if __name__ == '__main__':
  2.     Majority.populate() #生成数据

  4.     # Search for a solution
  5.     p = Population(BooleanFunction, 10, 3, 3, or_linker)#建立种群
  6.     print p

  8.     for _ in xrange(100):
  9.         if p.best.solved:#（检查是否满足要求）
  10.             break
  11.         p.cycle() #不满足就继续
  12.         print
  13.         print p
  14.         
  15.     if p.best.solved:
  16.         print
  17.         print 'SOLVED:', p.best #输出最终计算结果
复制代码

这部分所做的，调用了整个过程。由这个例子可以看出，如果我们要使用该包去计算一个问题时，我们所要做的是
1、建立一个变量列表所对应的数据表（如Majority中给出的是a，b，c，及对应的majority）
2、建立一个chromosome的类的继承，给定算符空间（functions and terminals）在继承中重写_fitness()和_solved()函数，
以满足我们具体问题的需要。
3、生成初始种群P0，并调用_solved()函数检查是否满足要求，调用cycle()函数进行种群改造。
这些就是我们所要做的事，如果你对比如选择算子，变异算子，等的概率设置不满意，同样可以到population类代码中进行修改。
同样对其他不满意也可以自己去修改，因为代码就在你手中。


又一个例子 ：
寻找一个多项式y=x**4+x**3+x**2+x
数据存在data.txt文件里。data.txt文件内容：
9.500366    9103.54918799079
-6.130432    1213.47815867463
3.252685    160.18146840485
7.88797        4432.23529247904
9.090484    7671.79392962917
1.485199    11.8327154232663
-3.950531    193.570365560718
10.003326    11124.3786278048
-.607453    -.326443121119579
5.469299    1093.78835980998
源代码：

  1. from pygep import *
  2. from pygep.functions.mathematical.arithmetic import add_op, subtract_op, multiply_op, divide_op
  3. from pygep.functions.linkers import sum_linker
  4. # example for y=x**4+x**3+x**2+x
  5. HH=add_op, subtract_op, multiply_op, divide_op
  6. class mp(object):
  7.     SAMPLE = []
  8.    
  9.     @staticmethod
  10.     def populate():
  11.                 data=open("data.txt").read().split()               
  12.                 for i in range(0,5):
  13.                         mp.SAMPLE.append(mp(float(data[2*i]),float(data[2*i+1])))                        
  14.         
  15.     def __init__(self, x,y):
  16.         self.x =x      
  17.         self.value = y


  20. class mpf(Chromosome):
  21.     functions = HH
  22.     terminals = tuple('x')
  23.    
  24.     def _fitness(self):
  25.         # Fitness function: number of hits
  26.         hits = 0
  27.         for i in mp.SAMPLE:
  28.             if abs(i.value -float(self(i)))/abs(i.value)<0.01:
  29.                 hits += 1
  30.         
  31.         return hits
  32.    
  33.     def _solved(self):
  34.         return self.fitness >= len(mp.SAMPLE)
  35. if __name__ == '__main__':
  36.     mp.populate()
  37.     # Search for a solution
  38.     p = Population(mpf, 10, 10, 1,sum_linker)
  39.     for _ in xrange(1000):
  40.         if p.best.solved:
  41.             break
  42.         p.cycle()
  43.     if p.best.solved:
  44.         print p
  45.         print 'SOLVED:', p.best
复制代码

result:
+*x++x*x/*xxxxxxxxxxx [2814]: 4
+x*x++x***xxxxxxxxxxx [2951]: 5
+*x++x*x/*xxxxxxxxxxx [2945]: 4
+*+++x*x/*xxxxxxxxxxx [2952]: 0
+*/++x*x/*xxxxxxxxxxx [2953]: 0
+*x++x*x/*xxxxxxxxxxx [2908]: 4
+*x++x*x/*xxxxxxxxxxx [2908]: 4
+*x++x*-+*xxxxxxxxxxx [2954]: 0
+*x++x*x/*xxxxxxxxxxx [2814]: 4
+*x++x*x/*xxxxxxxxxxx [2814]: 4
SOLVED: +x*x++x***xxxxxxxxxxx

so we get the result:x**4+x**3+x**2+x.

  1. 这里需要重新修改一下乘法和除法运算，在pygep.functions.mathematical.arithmetic 这个文件里。做如下修改：
  2. #divide_op   = symbol('/')(lambda i, j: float(i) / j)
  3. #modulus_op  = symbol('%')(lambda i, j: i % j)
  4. def divfun(i,j):
  5.         if abs(j)<10**(-12):
  6.                 return float(i)*10**12
  7.         else:
  8.                 return float(i)/j
  9. def modfun(i,j):
  10.         if abs(j)<10**(-12):
  11.                 return float(i)%10**(-12)
  12.         else:
  13.                 return float(i)/j
  14. divide_op   = symbol('/')(divfun)
  15. modulus_op  = symbol('%')(modfun)
复制代码


代码分析与应用

继续我们的讨论哈!
#--encoding=UTF-8--
#分析PyGEP
一. Chromosome类分析
1.MetaChromosome
这个metachromosome类只重写了type的new方法,并在new方法中为meta实例
添加了symbols字段,用以表示符号集
arity字段,用以识别符号集中的最大操作数目.比如'+'法函数又两个参数'a+b',则arity=2
_fitness字段用以得到适应性函数

2.Chromosome

    __metaclass__ = MetaChromosome #得到类型名
    __next_id = 1  
    gene_type = KarvaGene #表示gene类型是KarvaGene的
#下面重点分析一下代码
    @classmethod
    def generate(cls, head, genes=1, linker=default_linker):


        tail = head * (cls.arity - 1) + 1 #得到尾部长度

         while True:
            new_genes = [None] * genes #新建一个列表
            for i in xrange(genes):
            #调用KarvaGene,新建染色体序列
                new_genes = cls.gene_type(
                    [random.choice(cls.symbols)   for _ in xrange(head)] + \
                    [random.choice(cls.terminals) for _ in xrange(tail)], head
                )

            yield cls(new_genes, head, linker) #返回该染色体,也就是Chromosome的一个实例
            #注意这里使用yield生成的是一个代码序列.
#以下__init__代码:
     def __init__(self, genes, head, linker=default_linker):
        ''
        if head < 0:
            raise ValueError('Head length must be at least 0')
        if not genes:
            raise ValueError('Must have at least 1 gene')

        self.genes  = genes
        self.head   = head
        self.linker = linker

        # Unique number of the organism
        self.__id = type(self).__next_id #每一个实例对应一个id
        type(self).__next_id += 1     #这里设置了一个迭代,这个迭代是在上面的while True:..中实现.
#__call__函数
    def __call__(self, obj):
        return self.linker(*[g(obj) for g in self.genes])
#该函数只有一句,在return中执行,却是得到的我们染色体序列的值.具体需要分析karvaGene.这里不做解释.
#其他是一些操作.暂时不做解释.这些不影响我们的使用.
#现在有一个问题,如果我希望将一些染色体本身作为染色体序列中的一个gene来进行处理,那么该如何去做呢?
#这个问题的背景可以设想为,比如我知道当前的符号表达式是一些三角多项式,我要固定三角函数的形式,比如
#sin(x),那么我希望在生成的序列中要么是一些变量比如x,要么是sin(x)的形式.或者更复杂一些,我希望是
#sin(f(a)*x)的形式,其中f(a)是一些常数的组合,这些常数是需要我们在基因生成过程中来选择的.
#那么就需要考虑刚才的问题,这让我想到了<盗梦空间>中的多层梦境.
#和rainboy约定sin(x**2+x)*x来演示sin(f(x))*h(x),其中f,h需要在GEP中模拟.
#发现分别用种群来模拟,很难加入到种群遗传序列里.
#最终不得不用了一个最为简单的办法.重写__call__()函数.
#这里是两个例子.可见该方法还是有效的,目前看是最简单的,显然也是有很大限制的,这个一方面是由于我用的这个
#勉强使用的办法的缺陷,另一方面,由于需要在整体上进行进化,这也是不得已的办法.希望rainboy和各位朋友能想
#出更好的办法.
例子:sin(f(x))*h(x)
1.f(x)=x**2+x,h(x)=x
代码:(有很多与以上例子重复的地方,但是为了清晰起见还是把他放在这里.

  1. HH=add_op, subtract_op, multiply_op, divide_op
  2. class mp(object):
  3.     SAMPLE = []
  4.    
  5.     @staticmethod
  6.     def populate():
  7.          
  8.      
  9.         for i in range(0,5):
  10.             x=(i+1)/math.e*math.pi
  11.             mp.SAMPLE.append(mp(x,math.sin(x**2+x)*(x**3+x**2)))         
  12.         
  13.     def __init__(self, x,y):
  14.         self.x =x      
  15.         self.value = y


  18. class mpf(Chromosome):
  19.     functions = HH
  20.     terminals = tuple('x')
  21.     def _fitness(self):
  22.         # Fitness function: number of hits
  23.         hits = 0
  24.         for i in mp.SAMPLE:
  25.             if abs(i.value -float(self(i)))/abs(i.value)<0.01:
  26.                 hits += 1
  27.         
  28.         return hits
  29.    
  30.     def _solved(self):
  31.         return self.fitness >= len(mp.SAMPLE)
  32.    
  33.     def __call__(self, obj):
  34.         '''
  35.         Evaluates a given GEP chromosome against some instance.  The
  36.         terminals in the chromosome are assumed to be attributes on
  37.         the object instance provided (unless they are numeric constants).

  39.         @param obj: an object instance with terminal attributes set
  40.         @return:    result of evaluating the chromosome
  41.         '''
  42.         
  43.         return   math.sin(float(self.genes[0](obj)))*float(self.genes[1](obj))

  45. if __name__ == '__main__':
  46.     mp.populate()

  48.     # Search for a solution
  49.     p = Population(mpf, 10, 10, 2,sum_linker)
  50.    

  52.     for _ in xrange(1000):
  53.         if p.best.solved:
  54.             break
  55.         p.cycle()
  56.            
  57.     if p.best.solved:
  58.         print p
  59.         print 'SOLVED:', p.best
复制代码





运行结果:
*++xxxxx+xxxxx [876]: 1
*++xxxxx/xxxxx [902]: 1
*++xxxxxx*xxxx [853]: 1
x+xxxxxx//xxxx [904]: 0
*++xxxxx/xxxxx [878]: 1
*+*xxxxxxxxxxx [901]: 0
-++xxxxxx-xxxx [903]: 0
*++xxxxx++xxxx [897]: 1
*++xxxxx+xxxxx [876]: 1
+x*xxxxx+xxxxx [898]: 5
SOLVED: +x*xxxxx+xxxxx
这个代码就是表示x**2+x 和x
另外也算了一个h(x)=x**3+x**2的例子.
只需要在以上代码里将mp类中的初始数据改一下就可以了.
运行结果:
[Generation: 442  |  Best: #4664 (5)  |  Mean: 0.8]
012345678901234567890123456789012345678901
------------------------------------------
+-*--/+*//xxxxxxxxxxx*-*/+/*-+/xxxxxxxxxxx [4624]: 1
+-*--/+*//xxxxxxxxxxx*-*/+/*-+/xxxxxxxxxxx [4670]: 1
+-*+-/+*//xxxxxxxxxxx*-*/+-*-+/xxxxxxxxxxx [4661]: 0
+-**-/+*//xxxxxxxxxxx*-*/+//-+/xxxxxxxxxxx [4669]: 0
+-*--/+*//xxxxxxxxxxx*-*/+/+-+/xxxxxxxxxxx [4663]: 0
+-*--/--+*xxxxxxxxxxx*-*/+/*-+/xxxxxxxxxxx [4664]: 5
+-*-//+*/-xxxxxxxxxxx--*-*/+**-xxxxxxxxxxx [4666]: 0
+-*--/+*//xxxxxxxxxxx*--/+/*-+/xxxxxxxxxxx [4667]: 0
+-*--/+*//xxxxxxxxxxx*-*/+/*-+/xxxxxxxxxxx [4624]: 1
*-*/+/*-+/xxxxxxxxxxx+-*--/+*//xxxxxxxxxxx [4668]: 0
SOLVED: +-*--/--+*xxxxxxxxxxx*-*/+/*-+/xxxxxxxxxxx 这个符号序列解码出来就是x**2+x和x**3+x**2.
希望有更好的办法,开始的时候一直希望能通过设置算符,来满足要求,结果却没有实现哈.d

评论

我的评论：

发表评论

请 登录 后发表评论。还没有在Zeuux哲思注册吗？现在 注册 ！
王俊亚

回复 王俊亚  2011年04月20日 星期三 14:42

你好，还想问一下就是，在pyGEP中怎样加入数字，例如对一个方程y=3.5x**2+8.9*x的辨识的话怎样在基因中加入随即数字？谢谢

0条回复

王俊亚

回复 王俊亚  2011年03月17日 星期四 10:56

请问楼主一个问题，就是在包中的第二个例子majority，我的问题是，怎么理解这个例子，说是演化出布尔运算的规则，但是我该怎样去理解这个这个例子产生的结果，谢谢

1条回复

  • 张麟

    回复 张麟  2011年03月22日 星期二 21:16

    majority是一个运算符，在布尔运算中是常用的运算。
    他的例子实际上是给出了，majority的全部的变量与应变量的值。
    不知道你不知道怎么理解的是什么

    0条回复

暂时没有评论
反馈意见 帮助中心 服务条款 版权声明 关于哲思

Zeuux © 2024

京ICP备05028076号