- 元胞具有一定的计算能力,可以当作最小的非冯诺依曼架构并行计算机。在《复杂》一书当中介绍了采用遗传算法筛选出能进行二分类算法的一维元胞规则——即判断输入的元胞两种状态谁多谁少,哪一种多就全部演变成为这一种状态。
- 本文尝试复现此操作,由于元胞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]
分类效果:
小于一半 | 大于一半 |
---|---|
配置
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)
所有视觉为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次。
待改进
- 运算速度 太慢,日后有时间会进行向量化操作。
- 元胞的种类太丰富,上面实现的判断很简单,希望以后能设计出能够选择出逻辑运算的适应函数