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.
It is also available from the Models
module as Models.wright_fisher
.
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
numagents = 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:
model = ABM(Haploid)
AgentBasedModel with 0 agents of type Haploid no space scheduler: fastest
Create n
random individuals:
for i in 1:numagents
add_agent!(model, rand(model.rng))
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!(model, nagents(model))
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!(model::ABM) = sample!(model, nagents(model))
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!(model, dummystep, modelstep_neutral!, 20; adata = [(:trait, mean)])
data
step | mean_trait | |
---|---|---|
Int64 | Float64 | |
1 | 0 | 0.540674 |
2 | 1 | 0.555865 |
3 | 2 | 0.534324 |
4 | 3 | 0.570363 |
5 | 4 | 0.596332 |
6 | 5 | 0.582346 |
7 | 6 | 0.5772 |
8 | 7 | 0.584075 |
9 | 8 | 0.590567 |
10 | 9 | 0.58364 |
11 | 10 | 0.578553 |
12 | 11 | 0.62567 |
13 | 12 | 0.606685 |
14 | 13 | 0.592396 |
15 | 14 | 0.577259 |
16 | 15 | 0.551849 |
17 | 16 | 0.540901 |
18 | 17 | 0.547915 |
19 | 18 | 0.561568 |
20 | 19 | 0.561245 |
21 | 20 | 0.57563 |
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.
model = ABM(Haploid)
for i in 1:numagents
add_agent!(model, rand(model.rng))
end
modelstep_selection!(model::ABM) = sample!(model, nagents(model), :trait)
data, _ = run!(model, dummystep, modelstep_selection!, 20; adata = [(:trait, mean)])
data
step | mean_trait | |
---|---|---|
Int64 | Float64 | |
1 | 0 | 0.525768 |
2 | 1 | 0.684621 |
3 | 2 | 0.736674 |
4 | 3 | 0.794236 |
5 | 4 | 0.837069 |
6 | 5 | 0.845882 |
7 | 6 | 0.863521 |
8 | 7 | 0.889436 |
9 | 8 | 0.89822 |
10 | 9 | 0.912791 |
11 | 10 | 0.927714 |
12 | 11 | 0.921683 |
13 | 12 | 0.927152 |
14 | 13 | 0.92976 |
15 | 14 | 0.923554 |
16 | 15 | 0.915291 |
17 | 16 | 0.922749 |
18 | 17 | 0.934277 |
19 | 18 | 0.936009 |
20 | 19 | 0.93479 |
21 | 20 | 0.932559 |
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".