少し前まで、確率シミュレーションが好きではなくて、数理モデルを前にして力学系の手法だけで料理を試みていた僕ですが、友人からの熱烈な進言をうけて、C++でシミュレーションを回す練習を始めてみました。
集団ベースモデルによるHamilton and May (1977) シミュレーション
Hamilton and May modelとは、パッチ状の構造をもったメタ個体群における移動分散の進化を論じたモデルです。
たくさんのパッチを考えます。各パッチは1個体だけ収容できるとします。各パッチ内の個体は子供をたくさん産み、死にます。
そして子どもたちは確率 で移動して他のパッチへ移動しますが、残りの
は生まれたパッチに留まるとします。
ただし移動した子どもたちは確率
でしか生き残れないものとします。逆に言うと、移動した個体のうち、割合
は死んでしまうというものです。
さて、進化的に安定な移住率 はいくらでしょう？
移動分散率 をもつ突然変異個体が、移動分散率
をもつ集団で稼げる適応度は
\begin{align}
W = \frac{1-d_\mathrm{m}}{1-d_\mathrm{m}+(1-c)d} + \frac{(1-c) d_\mathrm{m}}{1-cd}
\end{align}
ですので
が任意の
で成立するための
の条件を求めると、
\begin{align}
d^* = \frac{1}{1+c}
\end{align}
と決まります。
個体ベースモデルはいささか時間がかかりますので、次のようなアルゴリズムにしたがって、集団が単型であるという仮定のもとで進化の方向を予測*1してみましょう。
- メルセンヌ・ツイスタによって実一様乱数
をとってくる。この値を、集団における移動分散率の初期値とする。
- メルセンヌ・ツイスタによって、[0,1] の実一様乱数
をとってくる。
という小さい値にたいして、
なら、突然変異は起こらないが、そうでない場合、在来型が一個体だけ変異型になる *2。
- 在来型
は突然変異の結果、
に変異する。ただし
は、メルセンヌ・ツイスタにしたがって、実一様乱数として選び（
）、もちろん
]となるように制限する。
- 変異型の適応度
を計算し、これが1より大きいなら（
）変異型がのっとる。そうでないなら在来型がとどまる。*3
- これを1000世代くりかえし、収束先を調べる。
- なお、パッチの数は無限大とし、移動のコスト
とした（したがって、ESSの予測値は 5/8=0.625になる）。
その結果がこちら。10回ほど走らせた結果をのせます。横軸は世代、縦軸は移動分散率です（つまり移動分散率の変化の時間経過をトラックしている）。初期値がバラバラなのは、ランごとに乱数を振っているためです。
ほぼすべて、0.625に収束してくれていますね。めでたしめでたし。なお、10回のランにかかった時間は、数秒でした。
後記： なお、このアルゴリズムには、収束安定性が保証されないという問題がある。すなわちESSに本当に到達できるかどうかは考えずに、シミュレーションをまわしている。そもそも、収束安定性が成り立たないような場合、「集団は常に単型」であるという仮定が、（実証的にはもちろんのこと）理論的にも極めてアヤシイものになってしまう。
