logo

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

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

シミュレーション

malthusian_growth_simulation2.gif

このポストでは、生成されたエージェントが死亡するアクションを与えることで、集団の逆成長をマクロの観点から模倣しようとしている。このシミュレーションの中で、空間や移動に関連するすべては、単に視覚化のためのものであり、実際の目的とは一切関係がない。

変数

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

パラメーター

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

アクション

（全エージェントが、毎ターン）

  • 死亡: d の確率でシミュレーションから除外される。
  • 移動: 二変量正規分布 $N \left( \mathbb{0} , \Sigma \right)$ から引き出されたベクトルに従って移動する。

コードレビュー

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

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

julia> d = 0.04 # 사망률
0.04

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

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

julia> Random.seed!(0);

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

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

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

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

julia> coordinate = rand(gaussian2, N0)'
500×2 Adjoint{Float64,Array{Float64,2}}:
  0.135821    0.165683
 -0.0706015  -0.0269708
  0.117323    0.0594672
  0.0129895  -0.0218035
 -0.102842    0.314866
 -0.137781   -0.152561
  0.0794965   0.162326
 -0.0692709  -0.0375145
  ⋮
 -0.22697     0.178074
  0.036792   -0.35028
 -0.253518    0.418463
 -0.226243   -0.111246
  0.120505   -0.0111552
  0.0882067   0.0526977
  0.398744   -0.114504
 -0.215111    0.259829

julia> N = N0
500

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

ステップ2. エージェントの死亡

julia> coordinate = coordinate[(rand(N) .> d),:]
481×2 Array{Float64,2}:
  0.135821    0.165683
 -0.0706015  -0.0269708
  0.117323    0.0594672
  0.0129895  -0.0218035
 -0.102842    0.314866
  0.0794965   0.162326
 -0.0692709  -0.0375145
  ⋮
 -0.253518    0.418463
 -0.226243   -0.111246
  0.120505   -0.0111552
  0.0882067   0.0526977
  0.398744   -0.114504
 -0.215111    0.259829

現在のエージェントの数だけ、0と1の間で乱数を引く。これが d より小さいということは、d の確率で 死亡 が起きたという意味である。上のスクリーンショットでは、19のエージェントが死亡して481のエージェントが残っていることが分かる。これは元のエージェント500個の約4%が死亡したことを示しており、私たちのパラメータ設定がほぼ正確に反映されたことを確認できる。

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

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

julia> coordinate = coordinate + rand(gaussian2, N)'
481×2 Array{Float64,2}:
  0.399141     0.492689
 -0.151015    -0.00161691
 -0.00843121   0.10119
  0.111592     0.00865223
  0.0591775    0.150437
  0.0060825    0.0383255
 -0.0652535   -0.0688566
 -0.702124    -0.485473
  ⋮
  0.437486    -0.354416
 -0.149827     0.156366
 -0.0737965   -0.168831
  0.344358    -0.136392
 -0.231089    -0.0695535
  0.446688     0.143493
  0.133667     0.156873

ステップ2. でエージェントが19個減ったので、$N$ も更新して、その時の個体数を記録した。各エージェントの座標は行列で表されており、二変量正規分布から $N$ 回引いて積み上げた $N \times 2$ 行列を足すことで、エージェントの移動が実装される。これを終えて ステップ2. に戻る。$t = t_{\text{end}}$ に達するまでこのプロセスを繰り返し、停止する。

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

全コード

以下は、このポストに使用されるジュリアのコードである。

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

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

N0 = 500 # 초기 인구수
d = 0.04 # 사망률
max_iteration = 180 # 시뮬레이션 기간
gaussian2 = MvNormal([0.0; 0.0], 0.04I) # 2차원 정규분포

Random.seed!(0);
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.)

    coordinate = coordinate[(rand(N) .> d),:]

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