ロジスティック写像

2021年12月11日

はじめに

ロジスティック写像について Python で実装する。

ロジスティック写像

理論

ロジスティック写像は、次式で表される。

$$x_{t+1} = \lambda x_t (1 - x_t)$$

$x_t$ は時刻 $t$ のときの値で、0 から 1 の間の値を取る。$\lambda$ はパラメータである。

これは、生物の個体数の変動を表したモデルで、増えすぎると減り、減りすぎると増えるといったフィードバックの機構が表されている。

実装

『カオスの素顔』(ニーナ・ホール編、講談社) によると、ロジスティック写像は BASIC で次のようなコードで表される。

FOR I=1 TO 15
X=LAMB*X*(1-X)
PRINT X
NEXT

いまさら BASIC もないので、Python で書いてみる。

def func(x, lamb):
    for i in range(15):
        x = lamb*x*(1. - x)
        print(x)

あとのためにパッケージを読み込んでおく。

from matplotlib import pyplot as plt
import numpy as np

初期値を x = 0.4 とし、lamb を変えてみる。

lamb = 2

func(0.4, 2.)

結果

0.48
0.4992
0.49999872
0.4999999999967232
0.5
0.5
0.5
0.5
0.5
0.5
0.5
0.5
0.5
0.5
0.5

値が 0.5 に落ち着いている。

lamb = 1 + √5

func(0.4, 1. + np.sqrt(5.))

結果

0.7766563145999497
0.5613325051679958
0.7968439545141753
0.5238665497858146
0.8071736905792098
0.5036756121803967
0.8089732746923843
0.5000874331796714
0.8090169696366287
0.5000000494766355
0.8090169943749396
0.5000000000000158
0.8090169943749475
0.49999999999999994
0.8090169943749473

値が 0.5 と 0.809 あたりを繰り返している。つまり、振動している。

lamb = 4

func(0.4, 4.)

結果

0.96
0.15360000000000013
0.5200281600000003
0.9983954912280576
0.006407737294172653
0.02546671278776609
0.09927263731020608
0.3576703231667294
0.918969052370147
0.297859732624244
0.8365572492210314
0.5469168719870904
0.9911952284917879
0.034908990027601214
0.1347614097714162

デタラメな値に見える。

可視化

matplotlib で簡単にグラフを描けるので、これで可視化する。

関数を変える。値を直接出力するのではなくて、配列で値を返すようにする。

def func2(x, lamb):
    xs = [x]
    for i in range(15):
        x = lamb*x*(1. - x)
        xs.append(x)
    return xs

lamb = 2

xs1 = func2(0.4, 2.)
plt.plot(range(len(xs1)), xs1, "o-")
plt.ylabel("x")
plt.show()

lamb = 1 + √5

xs2 = func2(0.4, 1. + np.sqrt(5.))
plt.plot(range(len(xs2)), xs2, "o-")
plt.ylabel("x")
plt.show()

lamb = 4

xs3 = func2(0.4, 4.)
plt.plot(range(len(xs3)), xs3, "o-")
plt.ylabel("x")
plt.show()

lamb = 4 でずっと計算し続けても、ほとんどデタラメな値に見える。

def func3(x, lamb):
    xs = [x]
    for i in range(1000):
        x = lamb*x*(1. - x)
        xs.append(x)
    return xs

xs = func3(0.4, 4.)
plt.plot(range(len(xs)), xs, "-")
plt.ylabel("x")
plt.show()

しかし、もともとは 2 次式なので、$(x_t, x_{t+1})$ をプロットすると放物線になる。

for i in range(len(xs)-1):
    x1 = xs[i]
    x2 = xs[i+1]
    plt.plot(x1, x2, ".")

plt.xlabel("$x_1$")
plt.ylabel("$x_2$")
plt.show()

つまり、明確に規則があるところから、デタラメに見える値が出てき得るということである。これが「カオス」である。

初期値依存性

カオスには初期値依存性があるといわれる。たとえば、初期値を $10^{-15}$ だけずらしてみると、次のようになる。

xs2 = func3(0.4 + 1e-15, 4.)
plt.plot(range(len(xs)), xs, "-")
plt.plot(range(len(xs2)), xs2, "-")
plt.xlim(0, 100)
plt.ylabel("x")
plt.show()

青がもとのもので、オレンジが初期値をずらしたものだが、途中から大幅にずれてしまい、数値誤差が「たかが誤差」とは言えない事態になってしまう。

分岐図

各 lamb の値に対して、どういう値に落ち着くかをプロットしてみる。

lambs = np.linspace(2., 4., 1001)
for lamb in lambs:
    xs = func3(0.4, lamb)
    plt.scatter([lamb]*500, xs[-500:], marker=".", s=1, linewidths=0)
    
plt.xlim(2., 4.)
plt.ylim(0., 1.)
plt.xlabel("lamb")
plt.ylabel("x")
plt.show()

1000 ステップ計算して、最後の 500 ステップの値をプロットしている。lamb が 2 から 3 あたりでは 1 つの値に落ち着いているが、3 を超えたところから 2 つの値に分かれる。lamb = 1 + √5 で振動していたのがこれにあたる。そのあとは 4 つに分かれ、8 つに分かれ…最終的にぐちゃぐちゃになる。この部分が「カオス」である。

lamb = 3.5 あたりの分岐している部分を拡大すると、次のようになる。

lambs = np.linspace(3.5, 3.6, 1001)
for lamb in lambs:
    xs = func3(0.4, lamb)
    plt.scatter([lamb]*500, xs[-500:], marker=".", s=1, linewidths=0)
    
plt.xlim(3.5, 3.6)
plt.ylim(0.4, 0.6)
plt.xlabel("lamb")
plt.ylabel("x")
plt.show()

さらに、lamb = 3.566 付近を拡大すると、次のようになる。

lambs = np.linspace(3.566, 3.571, 1001)
for lamb in lambs:
    xs = func3(0.4, lamb)
    plt.scatter([lamb]*500, xs[-500:], marker=".", s=1, linewidths=0)
    
plt.xlim(3.566, 3.571)
plt.ylim(0.485, 0.515)
plt.xlabel("lamb")
plt.ylabel("x")
plt.show()

このように、分岐部分を拡大していくと、ずっと同じようなパターンが現れ続けるという特徴がある。このように複雑なものが、ただ 1 つの単純な式から出てくることが不思議なところである。

参考

  • ニーナ・ホール編: カオスの素顔, 講談社, 1994