logo

エージェントベースモデルシミュレーションにおける繁殖 📂動力学

エージェントベースモデルシミュレーションにおける繁殖

シミュレーション

malthusian_growth_simulation.gif

このポストでは、生成されたエージェントに自身を複製するアクションを与え、マクロレベルでの集団成長を模倣しようとする。このシミュレーションで空間や移動に関連するすべては、視覚化のためだけのもので、実際の目的とは何の関係もない。

変数

  • $t$: 現在のターンを意味する。
  • $N(t)$: $t$ターンでのエージェントの数を示す。

パラメータ

  • $N_{0} \in \mathbb{N}$: シミュレーションが始まるときのエージェントの数を示す。
  • $b \in [0,1]$: 繁殖率で、エージェントが繁殖する確率を示す。
  • $t_{\text{end}} \in \mathbb{N}$: シミュレーションが終わるターンを定める。

アクション

（すべてのエージェントが、毎ターンごとに）

  • 繁殖: b の確率で自分の位置に新しいエージェントを作る。
  • 移動: 2変量正規分布 $N \left( \mathbb{0} , \Sigma \right)$から引いたベクトル通りに移動する。

コードレビュー

ステップ1. 初期化とエージェント生成

julia> N0 = 5 # 초기 인구수
5

julia> b = 0.05 # 번식률
0.05

julia> max\_iteration = 180 # 시뮬레이션 기간
180

julia> gaussian2 = MvNormal([0.0; 0.0], 0.02I) # 2차원 정규분포
IsoNormal(
dim: 2
μ: [0.0, 0.0]
Σ: [0.02 0.0; 0.0 0.02]
)

julia> Random.seed!(2);

julia> Time\_evolution = [] # 인구수를 기록하기 위한 스택
Any[]

エージェントベースのシミュレーションチュートリアルと同様にエージェントを生成する。ここでの各変数は上のパラメータ設定を反映している。

  • N0 $\leftarrow N_{0} = 5$
  • b $\leftarrow b = 0.05$
  • max_iteration $\leftarrow t_{\text{end}} = 180$

Random.seed!(2)は再現性のためのシード番号設定であり、time_evolutionは個体数を記録するためのスタックだ。シミュレーション自体には影響を与えないが、シミュレーションを繰り返し分析するために重要なテクニックである。

julia> coordinate = rand(gaussian2, N0)'
5×2 Adjoint{Float64,Array{Float64,2}}:
  0.104598   -0.105289
 -0.086056   -0.243734
 -0.0955465   0.0787217
 -0.121846    0.0732683
  0.176526    0.165338

julia> N = N0
5

$N$は個体数を示しているので、最初に$N \gets N_{0}$のような初期個体数を受け取る必要がある。

20201112_193047.png

エージェントの位置を散布図にして見ると、上のように密集していることがわかる。

ステップ2. エージェントの繁殖

julia> replicated = (rand(N) .< b)
5-element BitArray{1}:
 0
 0
 1
 0
 0

julia> new\_coordinate = coordinate[replicated,:]
1×2 Array{Float64,2}:
 -0.0955465  0.0787217

julia> coordinate = cat(coordinate, new\_coordinate, dims = 1)
6×2 Array{Float64,2}:
  0.104598   -0.105289
 -0.086056   -0.243734
 -0.0955465   0.0787217
 -0.121846    0.0732683
  0.176526    0.165338
 -0.0955465   0.0787217

現在のエージェントの数だけ、0と1の間で乱数を引く。これがbより小さいということは、bの確率で繁殖が起こったという意味だ。上のスクリーンショットでは、三番目のエージェントが$5\%$の確率で繁殖したことが、**ステップ1.**の座標の三番目の行がそのまま新しい行として追加されたことでわかる。

ステップ3. 記録とエージェント移動

julia> N = size(coordinate, 1)
6

julia> push!(time\_evolution, N)
1-element Array{Any,1}:
 6

julia> coordinate = coordinate + rand(gaussian2, N)'
6×2 Array{Float64,2}:
  0.0678177  -0.0642068
 -0.137464   -0.0672691
 -0.0571189   0.06278
 -0.0635334   0.226498
  0.18355     0.180084
  0.0920502   0.0682156

**ステップ2.**でエージェントがもう一つ増えたので、$N$もアップデートされ、この時の個体数を記録した。各エージェントの座標は行列で表されているので、$N$回引いた$N \times 2$行列を足すことで、エージェントが移動したことを実装できる。実際に散布図を描くと、**ステップ1.**とは異なり僅かに移動したことがわかる。

20201112_192927.png

これを終えて、**ステップ2.**に戻る。$t = t_{\text{end}}$に達するまで続けて繰り返す。

このプロセスを繰り返すだけで、ポスト上部のgifと同様の現象を観察できる。

全コード

以下は、このポストで使用されたJuliaのコードだ。

cd(@__DIR__) # 파일 저장 경로

@time using Plots
@time using Random
@time using Distributions
@time using LinearAlgebra

N0 = 5 # 초기 인구수
b = 0.05 # 번식률
max_iteration = 180 # 시뮬레이션 기간
gaussian2 = MvNormal([0.0; 0.0], 0.02I) # 2차원 정규분포

Random.seed!(2);
time_evolution = [] # 인구수를 기록하기 위한 스택
let
  coordinate = rand(gaussian2, N0)'
  N = N0

  anim = @animate for t = (0:max_iteration)/100
    plot(coordinate[:,1], coordinate[:,2], Seriestype = :scatter,
      markercolor = RGB(1.,94/255,0.), markeralpha = 0.4, markerstrokewidth	= 0.1,
      title = "t = $t", aspect_ratio = 1, size = [400,400],
      xaxis=true,yaxis=true,axis=nothing, legend = false)
      xlims!(-10.,10.)
      ylims!(-10.,10.)

    replicated = (rand(N) .< b)
    new_coordinate = coordinate[replicated,:]
    coordinate = cat(coordinate, new_coordinate, dims = 1)

    N = size(coordinate, 1)
    push!(time_evolution, N)
    coordinate = coordinate + rand(gaussian2, N)'
  end
  gif(anim, "malthusian_growth_simulation.gif", fps = 18)
end