复现《复杂》二分类元胞自动机

  • 元胞具有一定的计算能力,可以当作最小的非冯诺依曼架构并行计算机。在《复杂》一书当中介绍了采用遗传算法筛选出能进行二分类算法的一维元胞规则——即判断输入的元胞两种状态谁多谁少,哪一种多就全部演变成为这一种状态。
  • 本文尝试复现此操作,由于元胞DNA序列种类为2的128次方种,因此我找到的规则和书中所述不同,但是效果一致 。

先上贡分类器dna序列

[1, 0, 1, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 1, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 0, 0, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 1, 0, 1, 0, 0, 1, 0, 0, 1, 1, 1, 0, 1, 0, 1, 1, 0, 0, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0, 1, 0, 0, 1, 1, 0, 1, 1, 0, 1, 1, 0, 0, 0, 0, 0, 1, 0, 0]

分类效果:

小于一半 大于一半
小于一半数量
大于一半数量
小于一半数量(47)
大于一半数量(51)

配置

jyputer notebook
依赖包:

matrix可视化使用matplotlib 的matshow()函数

from numpy import *;#导入numpy的库函数
import numpy as np; #这个方式使用numpy的函数时,需要以np.开头。
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
%matplotlib inline
import time
import math

partone 元胞自动机dna编码

简单元胞自动机编码:

若一个元胞能观察到自己和左右两个元胞,则能观测2**3种类型,每种类型对应一种生死状态(0死/1生)即形成了这个元胞的生存策略,例如:
⚪ ⚪ ⚪——1
⭕ ⚪ ⚪——0
⚪ ⭕ ⚪——1
⚪ ⚪ ⭕——1
⭕ ⭕ ⚪——0
⭕ ⚪ ⭕——1
⚪ ⭕ ⭕——1
⭕ ⭕ ⭕——0
即可用序列10110110作为此元胞的DNA序列

typelist=[]
for i in range(4):
    typelist.extend(unique_list(permutations(innitlist(3,i), 3)))
print(typelist)
#[(0, 0, 0),
 #(1, 0, 0),
 #(0, 1, 0),
 #(0, 0, 1),
 #(1, 1, 0),
 #(1, 0, 1),
 #(0, 1, 1),
 #(1, 1, 1)]
DNA=[1,0,1,1,0,1,1,0]

迭代环境的搭建,主要是要接入DNA来控制迭代,这里我采用利用排列组合的list:
[(0, 0, 0),
(1, 0, 0),
(0, 1, 0),
(0, 0, 1),
(1, 1, 0),
(1, 0, 1),
(0, 1, 1),
(1, 1, 1)]
对应产生DNA字典

#4领域
class nodcell:
    def __init__(self, i,  list_one,dickey):
        self.__i = i
        self.lod=list_one[i]
###############构成首尾循环################        
        if i+1>=len(list_one):
            righti=0
        else:
            righti=list_one[i+1]          
        if i-1<0:
            lefti=list_one[len(list_one)-1]
        else:
            lefti=list_one[i-1]
###########################################                 
        self.neiboer_mat=(lefti,list_one[i],righti)
#dickey即DNA策略表
        self.dickey=dickey
        #print(self.neiboer_mat)
# 产生下一代
    def relod (self):
        key=self.neiboer_mat
#lod:  live   or  die  
        self.lod=dickey[key]
        return self.lod
                
#迭代一代
def on_generationa(initlist,dickey):
    nodlist=[]
    for i in range(len(initlist)):
        nod=nodcell(i, initlist,dickey)
        nodlist.append(nod)
    ons=[]
    for i in (nodlist):
        ons.append(i.relod())
    return ons

#迭代n代   
def all_gen(initlist,n,dickey):
    finalone=[]
    finalone.append(initlist)
    for i in range(n):
        initlist=on_generationa(initlist,dickey)
        finalone.append(initlist)
    return finalone
#准备好8种类型的list
def prp_originlist():

    originlist=[]
    for i in range(4):
        originlist.extend(unique_list(permutations(innitlist(3,i), 3)))
    len(originlist)
    return originlist
#以类型为key  存放 0/1  (DNA与类型关联)

def prop_key(originlist,n):
    dickey={}
    for i in range(len(originlist)):
        dickey[(originlist[i])]=DNA[i]
    return dickey

    DNA=[1,0,1,1,0,1,1,0]
    #初始化一个200长度的元胞组cc
    cc=list(random.randint(0,2,200))
#生成DNA字典
    dickey=prop_key(prp_originlist(),"")
#按照DNA字典迭代100次
    m=np.matrix(all_gen(cc,100,dickey))
#展示
    plt.matshow(m)
    
DNA=[1,0,1,1,0,1,1,0]

所有视觉为1的DNA数量为 2**9=512种

Part Two 升级元胞视野

元胞视野升级为左三右三,类型组合有2**7=128种
所有视觉为3的DNA数量为 -

  • 2**128=340282366920938463463374607431768211456

数不清了,因此一个一个找不可能找到,这里借助遗传算法来解决
遗传算法步骤:

  • 1 进行编码解码(元胞编码如上已变为为二进制)
  • 2 确定适应度函数——(稍后探讨)
  • 3根据轮盘赌选择算子,选取适应度较大的个体。
  • 4确定交叉概率Pc,对上一步得到的种群进行单点交叉。每次交叉点的位置随机。
  • 5确定变异概率Pm,每次变异点的位置随机。

完整程序

from numpy import *;#导入numpy的库函数
import numpy as np; #这个方式使用numpy的函数时,需要以np.开头。
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
%matplotlib inline
import time
import math
#4领域
class nodcell:
    def __init__(self, i,  list_one,thekey):
        
        
        
        self.__i = i
        self.dickey=thekey
        self.lod=list_one[i]
        if i+1>=len(list_one):
            self.neiboer_mat=(list_one[i-3],list_one[i-2],list_one[i-1],list_one[i],list_one[0],list_one[1],list_one[2])
        elif i+2>=len(list_one):
            self.neiboer_mat=(list_one[i-3],list_one[i-2],list_one[i-1],list_one[i],list_one[i+1],list_one[0],list_one[1])
        elif i+3>=len(list_one):
            self.neiboer_mat=(list_one[i-3],list_one[i-2],list_one[i-1],list_one[i],list_one[i+1],list_one[i+2],list_one[0])
        
            
        elif i-3<0:
            self.neiboer_mat=(list_one[len(list_one)-1],list_one[i-2],list_one[i-1],list_one[i],list_one[i+1],list_one[i+2],list_one[i+3])
        elif i+3<len(list_one):
            self.neiboer_mat=(list_one[len(list_one)-2],list_one[len(list_one)-1],list_one[i-1],list_one[i],list_one[i+1],list_one[i+2],list_one[i+3])
        elif i+3<len(list_one):
            self.neiboer_mat=(list_one[len(list_one)-3],list_one[len(list_one)-2],list_one[len(list_one)-1],list_one[i],list_one[i+1],list_one[i+2],list_one[i+3])
        else:
            self.neiboer_mat=(list_one[i-3],list_one[i-2],list_one[i-1],list_one[i],list_one[i+1],list_one[i+2],list_one[i+3])
            

            
          
        
        
        
        #print(self.neiboer_mat)
    def relod (self):
        
        key=self.neiboer_mat
        
        self.lod=self.dickey[key]
        return self.lod
        



#########################################对DNA字典的封装###########################################################################################
from scipy.special import comb, perm
from itertools import combinations, permutations
def innitlist(a,n):
    add=[]
    for i in range(n):
        add.append(1)
    for i in range(a-n):
        add.append(0)
    return add

def unique_list(list1):
    list2 = []
    
    for i in list1:
        if i not in list2:
            list2.append(i)
    return list2
def prp_originlist(n):
#  准备DNA字典
    originlist=[]
    for i in range(n+1):
        originlist.extend(unique_list(permutations(innitlist(n,i), n)))
    len(originlist)
    return originlist
def prop_key(originlist,DNAli):
    dickey={}
    for i in range(len(originlist)):
        dickey[(originlist[i])]=DNAli[i]
    return dickey
##################################################################################################################################
def randomkeylist(neibor=7,size=128,leng=20):
    keylist=[]
    originlist=prp_originlist(neibor)
    for i in range(leng):
        
        DNAli=list(random.randint(0,2,size))
        key=prop_key(originlist,DNAli)
        keylist.append(key)
    return keylist
keylist=randomkeylist()

def on_generationa(initlist,dickey):
    nodlist=[]
    for i in range(len(initlist)):
        
        nod=nodcell(i, initlist,dickey)
        nodlist.append(nod)
    ons=[]
    for i in (nodlist):
        ons.append(i.relod())
    return ons 




def all_gen(initlist,n,dickey):
    
    
    
    finalone=[]
    finalone.append(initlist)
    for i in range(n):
        initlist=on_generationa(initlist,dickey)
        finalone.append(initlist)
    return finalone





def plot_generate(init,generate,DNAli):
    dickey=prop_key(prp_originlist(7),DNAli)
    m=np.matrix(all_gen(init,generate,dickey))
    plt.matshow(m)

    plt.title(1)
def m_generate(init,generate,DNAli):
    dickey=prop_key(prp_originlist(7),DNAli)
    m=np.matrix(all_gen(init,generate,dickey))
    return m



适应度函数:元胞自动机目标为二分类器,不同于函数

  • 第一次尝试:仿照softmax函数归一,-log(x)函数及-log(1-x)作为两种情况进行判断
    结果:经过上百代迭代不论初始,最终元胞全为0
    这里没有让适应度函数同时考虑两种情况,一种情况成立 就被判断为高适应度
  • 第二次尝试:没有随机初始元胞,以用同一个元胞尝试DNA,得到的结果换一个元胞序列就变了

  • 第三次尝试:同时随机正反两例,设计正反适应度,最后结合到一起作为最终适应度。
    如下贴出:

############################################
#遗传算法的实现

import numpy as np
from scipy.optimize import fsolve, basinhopping

import timeit
#获取初始种群:
def getIntialPopulation(length,huge):
    Populationlist=[]
    for i in range(huge):
        init=list(random.randint(0,2,length))
        Populationlist.append(init)
    return Populationlist


dnak=prp_originlist(7)
#计算单个元胞的适应度
def get_one_fit_level(DNAli,lenss,dnak,generate=200):
    
    for i in range(1000):
        init=list(random.randint(0,2,lenss))
        if sum(init)>int(lenss/2):
            initp=init
            break

    for i in range(1000):
        init=list(random.randint(0,2,lenss))
        if sum(init)<int(lenss/2):
            initn=init
            break

    
    dickey=prop_key(dnak,DNAli)
    m1=np.matrix(all_gen(initp,generate,dickey))
    m2=np.matrix(all_gen(initn,generate,dickey))
    lel=aa
    fit=0
    xp=sum(m1[-1])
    xn=sum(m2[-1])
    if xp<=(lel*3)/4:
        fitp=0.001*(xp)/lel
    else:
        fitp=(xp)/lel
    
    
    if xn>=lel/4:
        fitn=0.001*(lel-xn)/lel
    else:
        fitn=(lel-xn)/lel
    fit=fitn*fitp
    return fit
        
 # 计算所有元胞的适应度       
def all_fit_mean(Populationlist,initenvir,generate=200):
    listfit=[]
    for i in Populationlist:
        listfit.append(math.exp(get_one_fit_level(i,lenss,dnak,generate)))
        
    fitsum=sum(listfit)
    li=[]
    for i in listfit:
        li.append(i/fitsum)
    return fitsum,li
    
        

    

#####轮盘赌注
#####轮盘赌注

def geration_bat(fitlist,Populationlist):
    licol=[]
    count=0
# 累计概率的计算
    for i in range(len(fitlist)):
        count+=fitlist[i]
        licol.append(count)
    choselist=[]    
#轮盘选取序号
    licol=[0]+licol
    for i in range(len(licol)-1):
        aim=random.random(1)[0]
        for i in range(len(licol)-1):
            if (aim>=licol[i])&(aim<licol[i+1]):
                choselist.append(i)
                #print(i)
            
                
    #print(len(choselist))
    newgeneration    =[]   
    for i in choselist :
        newgeneration.append(Populationlist[i])
    #print(len(newgeneration))
    return newgeneration
#交叉

def cross_dna(lista,listb,pc):
    np.random.seed(0)
   
    p = np.array([1-pc, pc])
    index = np.random.choice([0, 1], p = p.ravel())
    if index==1:
        
        pos=random.randint(0,128)
        mida=lista[pos]
        midb=listb[pos]
        #print(listb[pos])
        #print(lista[pos])
        newa=[]
        newb=[]
        for k,l,m in zip(lista,listb,range(len(listb))):
            if m==pos:
                newa.append(l) 
                newb.append(k) 
            else:
                newa.append(k) 
                newb.append(l) 

     
    return (newa,newb)
##############################变异################################
#配对
def cors_new_gen(pc,newgeneration):
    lens=len(newgeneration)
    for i in range(1000):
    
        c=random.randint(0,2,lens)
        if (sum(c)==int(lens/2)):
            break
    
    manl=[]
    fman=[]
    for i ,j in zip(newgeneration,c):
        if j==0:
            manl.append(i)
        else:
            fman.append(i)
    new_cors_generation=[]
    for i ,j in zip(manl,fman):
        a,b=cross_dna(i,j,pc)
        new_cors_generation.append(a)
        new_cors_generation.append(b)
    return new_cors_generation

cors_gen=cors_new_gen(1,newgeneration) 
        
len(cors_gen)
#突变
def R_change(cors_gen,pm=0.0001):
    new_list=[]
    
    for i in cors_gen:
        p = np.array([1-pm, pm])
        index = np.random.choice([0, 1], p = p.ravel())
        
        if index==1:
            pos=random.randint(0,128)
            newone=[]
            for j in range(len(i)):
                if j==pos:
                    if i[j]==0:
                        newone.append(1)
                    else:
                        newone.append(0)
                else:
                    newone.append(i[j])
            new_list.append(newone)
        else:
            new_list.append(i)
    return new_list

    

开始种群繁衍500代

#初始化种群
Populationlist=getIntialPopulation(128,200)
costlist=[]

for ge in range(500):
    aa=40
    init=list(random.randint(0,2,aa))
    #获取初始适应度
    fitsum,fitlist=all_fit_mean(Populationlist,init,40)
    #轮盘赌住选择夫妻
    newgeneration=geration_bat(fitlist,Populationlist)
    #进行交叉及编译
    cors_gen=cors_new_gen(0.7,newgeneration) 
    newlist=R_change(cors_gen,pm=0.05)  
    Populationlist=newlist
    costlist.append(fitsum)
    print(ge)
    print(fitsum)
##防止陷入局部最优点
    if costlist[-1]==costlist[-2]:
        break
    

计算机为i5四代,这里用的元胞的长度为40,计算速度为9秒/代,若用长度为200的元胞,计算速度为50秒/代,实际迭代效果适合的长度就行,但是元胞的自动繁殖时间需要足够。我取的是200次。

待改进

  • 运算速度 太慢,日后有时间会进行向量化操作。
  • 元胞的种类太丰富,上面实现的判断很简单,希望以后能设计出能够选择出逻辑运算的适应函数
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 212,332评论 6 493
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 90,508评论 3 385
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 157,812评论 0 348
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 56,607评论 1 284
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 65,728评论 6 386
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 49,919评论 1 290
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 39,071评论 3 410
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 37,802评论 0 268
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 44,256评论 1 303
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,576评论 2 327
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,712评论 1 341
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 34,389评论 4 332
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 40,032评论 3 316
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,798评论 0 21
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,026评论 1 266
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 46,473评论 2 360
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 43,606评论 2 350

推荐阅读更多精彩内容

  • 书名贪婪的大脑:为何人类会无止境地寻求意义作者(英)丹尼尔·博尔(Daniel Bor)译者林旭文豆瓣http:/...
    xuwensheng阅读 15,206评论 8 54
  • 氨基酸是蛋白质的基本结构单位,参与合成蛋白质的氨基酸有20种,可作为原料在核糖体工厂通过肽键连接形成多肽链,都有密...
    官敏慧阅读 5,642评论 0 7
  • 今天第一次翻到你的微博,那是完全没有我的一个领域。你的微博,是我们分手那天你申请的。或许,你也想着重新开始,在一个...
    罗罗哒罗阅读 204评论 0 0
  • 一直听身边人提起简书平台,今天晚上终于来啦 以后就要记录我的生活和故事了,今天有点兴奋呢。 晚安。
    clare922阅读 104评论 0 0
  • 简书于2013年4月22日晚上10点半上线了公测版本,正式接受用户的检验。 简书应该是在9月份开始内测的,掐指算来...
    简书阅读 1,958评论 5 24