Prey-Predator系における能力の進化 - 7. LV方程式に個体の能力に関する多様性と突然変異を導入

これまで

Prey-Predator系（被食者-捕食者系）における個体能力の進化のダイナミクスを見たいので，エージェントベースドモデルを作って進化シミュレーションしてみようとしていましたが，失敗しました（前回: d:id:swarm_of_trials:20081221）．
今回は，早々にあきらめていたロトカ・ヴォルテラ方程式で被食者・捕食者の能力の共進化をもう一度考え直してみようと思います．

ロトカ・ヴォルテラ方程式の拡張

個体の能力を考慮したLV方程式の記述

進化の対象として考える個体の能力は第4回進化する能力を考察する（d:id:swarm_of_trials:20080821）と同様に次のようにしたいと思います．

被食者: 捕食者から逃げる能力（回避能力）a，繁殖能力 b
捕食者: 被食者を捕まえる能力（捕食能力）h，繁殖能力 d

能力に関する値（a, b, d, h）を考えてロトカ・ヴォルテラ方程式を書くと，次のようになりました．
$\dot{x} = b x - f(a, h) x y$ (1)
$\dot{y} = -c y + f(a, h) d x y$ (2)
ここで $f(a, h)$ は回避能力 a の被食者が捕食能力 h の捕食者に捕食される確率です．
各項の意味は以下のようになります．

式1の第1項被食者による自然増（生まれた数 - 自然死した数）
式1の第2項被食者が捕食者に食べられる数
式2の第1項捕食者の自然死
式2の第2項被食者を食べることのできた捕食者が生む子供の数

ここで，ひとつ仮定を導入します．
それは各能力 a と b ，h と d の間にはトレードオフが存在し，一方が決まると他方が一意に決まるというものです．
すなわち， $a = g(b)$ （ $d g(b) / d b < 0$ ）， $h = g'(d)$ （ $d g'(d) / d d < 0$ ）（その逆も可）ということです．
なぜ，トレードオフを仮定するかというと，各能力の値は大きければ大きいほど良いため，このような仮定がなければ進化によってどんどん大きくなっていくことが自明だからです．
実際，繁殖行為中は無防備になりがち・出産などによって体力が減少するなどの制約も考えられるので，あながち不自然な仮定でもない気がします．

個体能力の多様性の導入

個体の能力が進化するためにはその能力に多様性が無ければなりません．そうでないと選択が起こらないので．
上で考慮した個体の能力は被食者・捕食者共に 2 つで，それらにトレードオフを仮定したことで，ある能力を持った個体を一つの番号 i でラベル付けすることができます．
すなわち，番号 i の被食者の能力を $a_i, b_i$ と表現できます．
とすると，ロトカ・ヴォルテラ方程式は次のようにできます．
$\dot{x_i} = b_i x_i - \sum_{j=1}^{M} f(a_i, h_j) x_i y_j$ (3)
$\dot{y_i} = -c y_i + \sum_{j=1}^{N} f(a_j, h_i) d_i x_j y_i$ (4)
ここで，N は被食者，M 捕食者の番号の最大値です．
式3の第2項は各能力の捕食者（j が 1 から N ）について出会って（確率は $x_i y_j$ に比例）捕食される確率 $f(a_i, h_j)$ を計算して和を取っています．
N, M をものすごく大きくすれば個体の能力は連続的に変化する（a, b, d, h が実数）として考えてもいいかなと思います．

突然変異の導入

次に突然変異の導入をします．
ここで，突然変異を「能力が少し変化する」と考え， i の個体の子が i + 1 になるいうことにします．
突然変異率を被食者は $\mu_x$ ，捕食者は $\mu_y$ とすると，方程式は以下のようになります．
$\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} + \frac{x_{i+1}}{2} - x_i)$ (5)
$\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} + \frac{y_{i+1}}{2} - y_i)$ (6)

ちなみに式5，6を整理すると以下のようになるので，
$\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}$ (5')
$\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}$ (6')
$1 << N, M$ かつ i と i + 1 の能力の差が非常に小さい場合（例えば $|a_i - a_{i + 1}| << 1$ ）の場合には，次のように式5，6の右辺第3項は偏微分形式にすることができます（ただし偏微分の分子部分の書き方がかなり怪しいです．．）．
$\dot{x_i} = b_i x_i - \sum_{j=1}^{M} f(a_i, h_j) x_i y_j + \frac{\mu_x}{2} \frac{\partial^2 x_{i}}{\partial i^2}$ (5'')
$\dot{y_i} = -c y_i + \sum_{j=1}^{N} f(a_j, h_i) d_i x_j y_i + \frac{\mu_y}{2} \frac{\partial^2 y_{i}}{\partial i^2}$ (6'')
これは反応拡散方程式と呼ばれる方程式の一種です．