Prey-Predator系における能力の進化 - 8. 個体の多様性・突然変異を導入したLV方程式の数値計算と考察

また，えらく間が開いてしまいました．

これまで

Prey-Predator系（被食者-捕食者系）における個体能力の進化のダイナミクスを見るために，前回（d:id:swarm_of_trials:20090111）はロトカ・ヴォルテラ方程式に集団内の能力の多様性と突然変異を導入しました．
それが以下の方程式です．
$\dot{x_i} = b_i x_i - \sum_{j=1}^{M} f(a_i, h_j) x_i y_j + \mu_x \frac{x_{i-1} - 2 x_i + x_{i+1}}{2}$ (1)
$\dot{y_i} = -c y_i + \sum_{j=1}^{N} f(a_j, h_i) d_i x_j y_i + \mu_y \frac{y_{i-1} - 2 y_i + y_{i+1}}{2}$ (2)
x は被食者の数で，a は回避能力，b は繁殖能力を表します．
y は捕食者の数で，h は捕食能力，d は繁殖能力を表します．
詳しいことは前回を見てください．

今回は適当な仮定を置いて，数値計算によってこの系の性質を調べ，最後にその結果について少し考察します．

仮の関数の具体化

前回，そこで作った方程式系には適切な具体的な関数がわからないものがあり，そこは，仮に未知の関数として $f(x)$ とか置いて話を進めました．
数値計算をするためには具体的な方程式がないとできませんので，まず，適当な仮定を置いて具体的な方程式を作ります．

実際の生態系の知見から仮定するなどいろいろなやり方があると思いますが，ここでの方針は，想定できるものの内，できるだけ簡単なものを仮定するということにします．

具体的な関数の形の作成

まず， $f(a, h)$ について考えます．これは回避能力 a の被食者が捕食能力 h の捕食者に捕食される確率（以下捕食率）を表します．
この関数がちゃんと捕食率を表すためには，少なくとも以下の三点を満たしている必要があります．

確率なので $f(a, h) \in [0,1]$
回避能力 a が高いことは捕食率 f を減少させるはずなので， $\partial f(a, h) / \partial a \leq 0$
捕食能力 h が高いことは捕食率 f を増加させるはずなので， $\partial f(a, h) / \partial h \geq 0$

この三点を満たす関数として $f(a, h) = 1 / (1 + \exp(-\gamma (h - a)))$ と仮定します（ $\gamma \in \mathcal{R^+}$ ）．

次に， $a = g(b)$ ， $h = g'(d)$ について考えます．
これらは被食者の能力 a と b，捕食者の能力 h と d にトレードオフ関係があることを表しています．
トレードオフなので，これらの関数の条件はそれぞれ以下のようになります．

b の増加は a の減少なので， $d g(b) / d b < 0$ ．
d の増加は h の減少なので， $d g'(d) / d d < 0$ ．

そこで次のようにしました．

$a = \frac{A}{B} (B - b)$
$d = \frac{D}{H} (H - h)$

大文字は対応する小文字で表される能力の最大値を表します．
ここから，上の用にインデックス i を用いて各能力を表すと次のようになります．

$a_i = A \frac{N - i}{N}$
$b_i = B \frac{i}{N}$
$d_i = D \frac{M - i}{M}$
$h_i = H \frac{i}{M}$

これで未知の関数は無くなりました．

数値計算

パラメータはたくさんありますが適当に設定して数値計算しました．
ここでは，突然変異率 $\mu_x, \mu_y$ のみ変えて方程式の挙動を調べます．
以下で示す図は式 1, 2 を数値計算したもので，初期状態からある程度時間が経過した後の軌道を書いています．
横軸は被食者の総数 $\sum_{i=1}^N x_i$ ，縦軸は捕食者の総数 $\sum_{i=1}^N y_i$ にして示しています．
この経過時間は各図が見やすいように設定ごとに調整しています．

まず，個体の多様性はあるが，突然変異は全く起こらない場合（ $\mu_x = 0, \mu_y = 0$ ）場合についてやってみました．

被食者数の初期値を変えて 4 つの軌道を書いています．
sum(x) = p は $\sum_{i=1}^N x_i = p$ であることを表します．
また，初期値において $x_i = x_j$ $\forall i, j (i \not= j)$ としました．

同図より初期値によって軌道が変わっていることがわかります．
このことから $\mu_x = 0, \mu_y = 0$ の式1，2は中立安定な周期解を持つことがわかります．
この性質はロトカ・ヴォルテラ方程式系と同様です．
すなわち，被食者・捕食者系の解軌道の性質は個体の多様性を導入しても変わらないということが予想されます．

次に，突然変異率が非常に小さい場合（ $\mu_x = 0.01, \mu_y = 0.01$ ）の場合についてやってみました．

こちらの設定では初期値を変えても，最終的には同様の軌道になりました．
このことから $\mu_x = 0.01, \mu_y = 0.01$ の式1，2は漸近安定な周期解（リミットサイクル）を持つことがわかります．

次に，被食者の突然変異率が小さく，捕食者の突然変異率が大きい場合（ $\mu_x = 0.01, \mu_y = 1.0$ ）場合についてやってみました．

こちらの設定でも初期値に最終的な軌道は依存せず，同じ周期軌道に収束しました．
このことから $\mu_x = 0.01, \mu_y = 1.0$ の式1，2も $\mu_x = 0.01, \mu_y = 0.01$ の場合と同様に漸近安定な周期解（リミットサイクル）を持つことが予想されます．

次に， $\mu_x = 1.0, \mu_y = 0.01$ の場合についてやってみました．

この設定では初期値に最終的な軌道は依存せず，同じ点に収束しました．
このことから $\mu_x = 1.0, \mu_y = 0.01$ の式1，2は漸近安定な不動点（固定点アトラクタ）を持つことがわかります．

次に， $\mu_x = 1.0, \mu_y = 1.0$ の場合についてやってみました．

この設定でも初期値に最終的な軌道は依存せず，同じ点に収束しました．
このことから $\mu_x = 1.0, \mu_y = 1.0$ の式1，2も $\mu_x = 1.0, \mu_y = 0.01$ の場合と同様に漸近安定な不動点（固定点アトラクタ）を持つことがわかります．

結果のまとめ

Prey-Predator系（被食者・捕食者系）における個体能力の進化のダイナミクスを見るために，ロトカ・ヴォルテラ方程式系に集団内の能力の多様性と突然変異を導入（前回（d:id:swarm_of_trials:20090111））し，それの数値計算をしました．
その結果，被食者・捕食者系は個体の多様性・突然変異を考慮に入れることで，ロトカ・ヴォルテラ方程式系とは大きく異なった解軌道の性質を持つことがわかりました．
第一に突然変異を考えず，個体の多様性のみを考慮に入れた場合では，本記事の方程式系は中立安定な周期解を持ち，その点でロトカ・ヴォルテラ方程式系と同様の性質を持っていました．
ただし，少しでも突然変異が存在すると，本記事の方程式系は安定な軌道（アトラクタ）を持ちました．
これは元のロトカ・ヴォルテラ方程式系とは大きく異なる性質です．
さらに，被食者の突然変異率に依存してアトラクタの種類が変化し，（やってみた範囲では）捕食者の突然変異率には依存しませんでした．
具体的には，被食者の突然変異率が小さいと本記事の方程式系はリミットサイクルを持ちます．
そして，被食者の突然変異率が大きい場合，固定点アトラクタを持ちます．
分岐に関してはちゃんと調べていませんが，被食者の突然変異率に関してホップ分岐が存在していることがわかります．

この結果の意味するところ

突然変異を考慮に入れる / 入れない，突然変異率の大小の違いが被食者・捕食者系の振る舞いを大きく変えることは何を意味するのでしょうか？
あまり深く考えず，思ったことを書いてみます．

ロトカ・ヴォルテラ方程式系で考えている個体数動力学に多様性や突然変異を導入することは考えているスケールを大きく変えることであると言えます．
おそらく元のロトカ・ヴォルテラ方程式系は，個体数が多いことを想定して，個体の多様性を無視しています．
そのため個体の多様性を考慮することは考えている系の大きさを変えることと言えます．

また，元のロトカ・ヴォルテラ方程式系は，被食者と捕食者の数の変化は突然変異が起こる頻度よりずっと早いと想定して，突然変異の存在を無視していると思われます．
そのため突然変異を考慮することは考えている系のタイムスケールを変えることと言えます．
突然変異が頻繁に起こるか，または，個体数の変動がすごくゆっくりな場合に相当すると考えています．
一方の突然変異が，他方に比べて非常に大きい場合というのも，世代交代などの早さが両者で異なることを表していると言えます．

すなわち本記事の結果は，被食者・捕食者系では見ているスケールを変えると系の性質は全く変わってしまうということを示唆しています．
また，本記事ではいくつかの仮定を置いて，ロトカ・ヴォルテラ方程式系のスケールを変えてみましたが，他のスケールの変え方もあり，それはまた別の結果を導き出す可能性があります．
フィールド・実験室において被食者・捕食者系にも関わらずロトカ・ヴォルテラ方程式系とは全くことなった挙動を示す場合*1，このような想定しているスケールの違いがあるのかもしれません．

今後は

もうちょっとちゃんと解析できたらいいと思います（せめて分岐に関する相図ぐらい．．．）．

また，上の「この結果の意味するところ」はかなり適当に考察しているので，もうちょっとちゃんと考えてみたいです．

*1:とはいうものの，そういうことがあるのか私は知らないのですが

Reports

Twitter: @mtknnktm.