Wright-Fisher model of evolution
This is one of the simplest models of population genetics that demonstrates the use of sample!
. We implement a simple case of the model where we study haploids (cells with a single set of chromosomes) while for simplicity, focus only on one locus (a specific gene). In this example we will be dealing with a population of constant size.
A neutral model
- Imagine a population of
n
haploid individuals. - At each generation,
n
offsprings replace the parents. - Each offspring chooses a parent at random and inherits its genetic material.
using Agents
n = 100
Let's define an agent. The genetic value of an agent is a number (trait
field).
mutable struct Haploid <: AbstractAgent
id::Int
trait::Float64
end
And make a model without any spatial structure:
m = ABM(Haploid)
AgentBasedModel with 0 agents of type Haploid
no space
scheduler: fastest
Create n
random individuals:
for i in 1:n
add_agent!(m, rand())
end
To create a new generation, we can use the sample!
function. It chooses random individuals with replacement from the current individuals and updates the model. For example:
sample!(m, nagents(m))
The model can be run for many generations and we can collect the average trait value of the population. To do this we will use a model-step function (see step!
) that utilizes sample!
:
modelstep_neutral!(m) = sample!(m, nagents(m))
We can now run the model and collect data. We use dummystep
for the agent-step function (as the agents perform no actions).
using Statistics: mean
data, _ = run!(m, dummystep, modelstep_neutral!, 20; adata = [(:trait, mean)])
data
step | mean_trait | |
---|---|---|
Int64 | Float64 | |
1 | 0 | 0.42748 |
2 | 1 | 0.466301 |
3 | 2 | 0.491889 |
4 | 3 | 0.519205 |
5 | 4 | 0.538416 |
6 | 5 | 0.519753 |
7 | 6 | 0.541751 |
8 | 7 | 0.527809 |
9 | 8 | 0.538966 |
10 | 9 | 0.532131 |
11 | 10 | 0.537809 |
12 | 11 | 0.538059 |
13 | 12 | 0.539744 |
14 | 13 | 0.604315 |
15 | 14 | 0.636157 |
16 | 15 | 0.648527 |
17 | 16 | 0.631555 |
18 | 17 | 0.667634 |
19 | 18 | 0.664524 |
20 | 19 | 0.637829 |
21 | 20 | 0.639602 |
As expected, the average value of the "trait" remains around 0.5.
A model with selection
We can sample individuals according to their trait values, supposing that their fitness is correlated with their trait values.
m = ABM(Haploid)
for i in 1:100
add_agent!(m, rand())
end
modelstep_selection!(m::ABM) = sample!(m, nagents(m), :trait)
data, _ = run!(m, dummystep, modelstep_selection!, 20; adata = [(:trait, mean)])
data
step | mean_trait | |
---|---|---|
Int64 | Float64 | |
1 | 0 | 0.511277 |
2 | 1 | 0.694014 |
3 | 2 | 0.785328 |
4 | 3 | 0.828866 |
5 | 4 | 0.868187 |
6 | 5 | 0.869305 |
7 | 6 | 0.880079 |
8 | 7 | 0.901108 |
9 | 8 | 0.897571 |
10 | 9 | 0.911525 |
11 | 10 | 0.907961 |
12 | 11 | 0.925238 |
13 | 12 | 0.934203 |
14 | 13 | 0.944933 |
15 | 14 | 0.948984 |
16 | 15 | 0.952078 |
17 | 16 | 0.961232 |
18 | 17 | 0.96482 |
19 | 18 | 0.972801 |
20 | 19 | 0.977566 |
21 | 20 | 0.978035 |
Here we see that as time progresses, the trait becomes closer and closer to 1, which is expected - since agents with higher traits have higher probability of being sampled for the next "generation".