logo

Agent-based Model Simulation of Mortality 📂Dynamics

Agent-based Model Simulation of Mortality

Simulation

malthusian_growth_simulation2.gif

This post aims to simulate the degrowth of a population from a macroscopic perspective by implementing an action in which generated agents die. Everything related to space or movement in this simulation is merely for visualization and has no actual purpose.

Variables

  • $t$: Represents the current turn.
  • $N(t)$: Represents the number of agents at turn $t$.

Parameters

  • $N_{0} \in \mathbb{N}$: Represents the number of agents at the start of the simulation.
  • $d \in [0,1]$: Represents the probability of an agent dying, referred to as the death rate.
  • $t_{\text{end}} \in \mathbb{N}$: Determines the turn at which the simulation ends.

Actions

(For all agents, every turn)

  • Death: Excluded from the simulation with a probability of d.
  • Movement: Moves according to a vector drawn from a bivariate normal distribution $N \left( \mathbb{0} , \Sigma \right)$.

Code Review

Step 1. Initialization and Agent Creation

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[]

Similar to the Agent-Based Simulation Tutorial, this step involves creating agents. Here, each variable reflects the parameters set above.

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

Random.seed!(0) is set to reproduce results, and time_evolution is a stack for recording the number of entities, which does not impact the simulation itself but is an important technique for repeating and analyzing simulations.

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

Since $N$ represents the number of entities, initially, it must be provided as in $N \gets N_{0}$.

Step 2. Agent Death

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

For the current number of agents, random numbers between 0 and 1 are drawn. If this is less than d, it means that death has occurred with a probability of d. In the screenshot above, 19 agents died, leaving 481 agents. This shows that about 4% of the original 500 agents died, indicating that our parameter settings have been almost accurately reflected.

Step 3. Recording and Agent Movement

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

Since 19 agents were reduced in the Step 2., $N$ was also updated, and the number of entities at that time was recorded. Agents’ coordinates are represented as matrices, so adding a $N \times 2$ matrix, which is stacked by drawing $N$ times from a bivariate normal distribution, implements agent movement. After completing this, it returns to Step 2. Repeat this process until reaching $t = t_{\text{end}}$, and then stop.

By just repeating the above processes, one can observe the phenomenon similar to the gif at the top of the post.

Full Code

Below is the Julia code used in this post.

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