02 離散ロジスティック成長モデル

指数増殖モデルでは個体数が制限なく増え続けるが，現実の個体群には様々な制約がある．

密度依存性 (density dependence)を取り入れた離散ロジスティック成長モデル (discrete logistic growth model)を実装し，平衡点 (equilibrium point)と局所安定性 (local stability)を調べる．

離散ロジスティック成長モデル¶

離散ロジスティック成長モデルは以下で定義される：

X_{t+1} = X_t + r \left(1 - \frac{X_t}{K}\right) X_t

(1)

ここで，

$r$ ：内的自然増加率 (intrinsic rate of increase)　個体数が十分小さいときの1世代あたりの増殖率
$K$ ：環境収容力 (carrying capacity)　ある環境で維持されうる個体数

$1 - X_t / K$ が密度依存性を表現し， $X_t \to K$ で増加率 $\to 0$ となる．

実装¶

# 02-01. 離散ロジスティック成長モデル
import matplotlib.pyplot as plt

Solution to Exercise 1

r = 0.5
K = 100
x = 1
t = 0

t_list = [t]
x_list = [x]
for i in range(50):
    t = t + 1
    x = x + r * x * (1 - x / K)
    t_list.append(t)
    x_list.append(x)

plt.figure(dpi=200)
plt.plot(t_list, x_list, "-", t_list, x_list, "r.")
plt.title("Logistic growth", fontsize="xx-large")
plt.xlabel("Time (t)", fontsize="x-large")
plt.ylabel("Pop. size (x)", fontsize="x-large")

平衡点と局所安定性¶

平衡点 $\bar X$ とは $X_{t+1} = X_t$ を満たす個体数で， $f(X) = X + r X(1 - X/K)$ とおくと $\bar X = f(\bar X)$ の解として求まる．

離散ロジスティック成長モデルの平衡点は：

$\bar X = 0$ （自明な平衡点：個体群が絶滅した状態）
$\bar X = K$ （非自明な平衡点：環境収容力に達した状態）

平衡点からの微小なずれ（ $n_t$ ）が生じた際にどうなるか（局所安定性）を $f'(\bar X)$ の絶対値で判定しよう（詳しくはスライド参照）．

$|f'(\bar X)| < 1$ ならば $n_t \to 0$ となり，平衡点 $\bar N$ は局所安定．
$|f'(\bar X)| > 1$ ならば $|n_t| \to \infty$ となり，平衡点 $\bar N$ は不安定．

離散ロジスティック成長モデルの平衡点でどうなるか，場合分けしてみよう．

パラメータと動態の多様性¶

$r$ の値に応じて離散ロジスティックモデルは安定・周期振動・カオスなど多様な動態を示す．

まずは $r$ を変えて時間発展がどう変化するかを観察する．

Solution to Exercise 2

K = 100
x0 = 1
t_end = 100

for r in [0.5, 1.5, 2.0, 2.5, 2.9]:
    x = x0
    t_list = [0]
    x_list = [x]
    for i in range(t_end):
        x = x + r * x * (1 - x / K)
        t_list.append(i + 1)
        x_list.append(x)

    plt.figure(dpi=120)
    plt.plot(t_list, x_list, "-")
    plt.title(f"Discrete logistic, r = {r}")
    plt.xlabel("Time (t)")
    plt.ylabel("Pop. size (x)")
    plt.show()

$r$ が小さいうちは $K$ への収束だが， $r$ が大きくなるにつれて振動が現れ，さらに不規則な動態（カオス）が現れる．

分岐図¶

時間発展を $r$ ごとに別々に描くことで動態の違いは分かる．

さらに，パラメータ全体での挙動を一望するには分岐図 (bifurcation diagram)が役に立つ．

素朴な作り方として，各 $r$ で十分長く回した最終数世代の $X$ をプロットする方法を試す．

# 02-02. 分岐図

Solution to Exercise 3

K = 100
r_list = []
x_list = []
for i in range(150, 300):
    r = i / 100
    x = 10
    for t in range(1900):
        x = x + r * x * (1 - x / K)

    for t in range(100):
        x = x + r * x * (1 - x / K)
        r_list.append(r)
        x_list.append(x)

plt.figure(dpi=200)
plt.plot(r_list, x_list, ".", markersize=1)
plt.title("Bifurcation diagram")
plt.xlabel("r")
plt.ylabel("x")

$r$ が小さい範囲では平衡点 $\bar X = K$ が見える． $r$ が大きくなるにつれ周期解やカオスに入る．