前提・実現したいこと
初めて質問させていただきます。
N状態K近傍のセルオートマトン(CA)を実装して
ラングトンが提唱したλパラメータと相互情報量のグラフを出力し、カオスの縁を確認したいと思っています。
しかし、いざコードを実行してみると、状態4あたりでλの値が後半で急に上昇してしまいます。
おそらく一般化したCAのプログラムが違うと思うのですが、どこをどう直せがいいかさっぱり分かりません。
もし、間違っているところがあれば指摘して欲しいです。よろしくお願いいたします。
該当のソースコード
Python
1import matplotlib.pyplot as plt 2import numpy as np 3import random as rnd 4 5#一次元データの系列をリストで受け取り、そのエントロピーを返す 6def calcEntropy(data): 7 dic = {} 8 for d in data: 9 if d in dic: 10 dic[d] = dic[d] + 1 11 else: 12 dic[d] = 1 13 14 probdist = np.array(list(dic.values()), dtype = float) / float(len(data)) 15 return(np.sum([- p * np.log2(p) for p in probdist])) 16 17#2種類のリストの結合エントロピーを返す 18def calcJointEntropy(x,y): 19 data = [] 20 21 for m in range(len(x)): 22 data1 = x[m] 23 data2 = y[m] 24 data.append((data1,data2)) 25 26 dic = {} 27 28 for d in data: 29 if d in dic: 30 dic[d] = dic[d] + 1 31 else: 32 dic[d] = 1 33 34 probdist = np.array(list(dic.values()), dtype = float) / float(len(data)) 35 return(np.sum([- p * np.log2(p) for p in probdist])) 36 37#2種類のリストの相互情報量を返す 38def calcMI(x,y): 39 H_X = calcEntropy(x) 40 H_Y = calcEntropy(y) 41 H_XY = calcJointEntropy(x,y) 42 return(H_X + H_Y - H_XY) 43 44#CAのデータから各セルの相互情報量を算出し、リストとして返す 45def calcCAMIList(data): 46 MIList = [] 47 48 for l in range(len(data)): 49 List = [] 50 dd = [row[l] for row in data] 51 52 for t in range(0,len(dd)): 53 x = [] 54 y = [] 55 for i in range(6): 56 x.append(dd[(t+i)%len(dd)]) 57 y.append(dd[(t+1+i)%len(dd)]) 58 a = calcMI(x,y) 59 List.append(a) 60 61 MIList.append(List) 62 63 return(MIList) 64 65#全てのセルの相互情報量を平均を返す 66def calcCAMI(data): 67 d = calcCAMIList(data) 68 List = [] 69 darray = np.array(d) 70 for l in range(len(darray)): 71 L = np.sum(darray[l]) / darray.shape[0] 72 List.append(L) 73 74 List_n = np.array(List) / len(List) 75 return(np.sum(List_n)) 76 77#K状態N近傍のCAを計算する 78def NK_CA(l, t, cell_i, K, N, n): 79 cell = cell_i 80 data = [cell] 81 lambda1 = (K**N - n)/K**N 82 83#K状態でのRULEの振る舞い 84 np.random.seed(10000) 85 RULE = [] 86 for i in range(K**N): 87 if np.random.rand()< lambda1: 88 r = np.random.randint(1,K) 89 RULE.append(r) 90 91 else:RULE.append(0) 92 93#N近傍でのneighboringstateの振る舞い。N/2の切り落としでどこから見ていくのかを割り出し、一つずつ足して行きました。 94 for t in range(t): 95 cell_next = [0 for i in range(l)] 96 for j in range(l): 97 neighboringstate = 0 98 for m in range(N): 99 neighboringstate = cell[int(j+(m-(N-1)/2)+l)%l]*(K**(N-1-m)) 100 101 cell_next[j] = RULE[neighboringstate] 102 103 cell = cell_next 104 data.append(cell) 105 106 return(data,lambda1) 107 108#横軸λ、縦軸相互情報量の平均としたグラフに表す 109L = 101 110T = 100 111cell_init = [0 for i in range(L)] 112cell_init[L//2] = 1 113 114lambda2 = [] 115mi = [] 116 117for k in range(2,5): 118 for n in range(3,5): 119 for n_n in range(1, k**n+1): 120 dataXY, lam = KN_CA(L, T, cell_init, k, n, n_n) 121 calc_mi,mi_sum = calcCAMI(dataXY) 122 lambda2.append(lam) 123 mi.append(mi_sum) 124 125 126fig = plt.figure(figsize=(5,5)) 127ax = fig.add_subplot(1,1,1) 128ax.plot(lambda2, mi,'.') 129plt.xlabel('lambda') 130plt.ylabel('averageMI') 131plt.show()
試したこと
NK_CAのneighboringstateでjから離れるほどかける値を小さくなるように組みましたが、余計値がばらついてしまいました。
###環境
MacOS Catalina バージョン10.15.4
Jupyter Notebook
あなたの回答
tips
プレビュー