SCJ-PredatorPrey/predator_prey.jl

350 lines
12 KiB
Julia

# # Predator-prey dynamics
# ```@raw html
# <video width="auto" controls autoplay loop>
# <source src="../sheepwolf.mp4" type="video/mp4">
# </video>
# ```
# The predator-prey model emulates the population dynamics of predator and prey animals who
# live in a common ecosystem and compete over limited resources. This model is an
# agent-based analog to the classic
# [Lotka-Volterra](https://en.wikipedia.org/wiki/Lotka%E2%80%93Volterra_equations)
# differential equation model.
# This example illustrates how to develop models with
# heterogeneous agents (sometimes referred to as a *mixed agent based model*),
# incorporation of a spatial property in the dynamics (represented by a standard
# array, not an agent, as is done in most other ABM frameworks),
# and usage of [`GridSpace`](@ref), which allows multiple agents per grid coordinate.
# ## Model specification
# The environment is a two dimensional grid containing sheep, wolves and grass. In the
# model, wolves eat sheep and sheep eat grass. Their populations will oscillate over time
# if the correct balance of resources is achieved. Without this balance however, a
# population may become extinct. For example, if wolf population becomes too large,
# they will deplete the sheep and subsequently die of starvation.
# We will begin by loading the required packages and defining two subtypes of
# `AbstractAgent`: `Sheep`, `Wolf`. Grass will be a spatial property in the model. All three agent types have `id` and `pos`
# properties, which is a requirement for all subtypes of `AbstractAgent` when they exist
# upon a `GridSpace`. Sheep and wolves have identical properties, but different behaviors
# as explained below. The property `energy` represents an animals current energy level.
# If the level drops below zero, the agent will die. Sheep and wolves reproduce asexually
# in this model, with a probability given by `reproduction_prob`. The property `Δenergy`
# controls how much energy is acquired after consuming a food source.
# Grass is a replenishing resource that occupies every position in the grid space. Grass can be
# consumed only if it is `fully_grown`. Once the grass has been consumed, it replenishes
# after a delay specified by the property `regrowth_time`. The property `countdown` tracks
# the delay between being consumed and the regrowth time.
# ## Making the model
# First we define the agent types
# (here you can see that it isn't really that much
# of an advantage to have two different agent types. Like in the [Rabbit, Fox, Wolf](@ref)
# example, we could have only one type and one additional filed to separate them.
# Nevertheless, for the sake of example, we will use two different types.)
using Agents, Random
@agent struct Sheep(GridAgent{2})
energy::Float64
reproduction_prob::Float64
Δenergy::Float64
perception::Int32
nearby_agents
nearby_grass
#speed::Float64
#endurance::Float64
end
function perceive!(sheep::Sheep,model)
sheep.nearby_agents = nearby_agents(sheep, model, model.sheep_perception)#sheep.perception)
sheep.nearby_grass = nearby_fully_grown(sheep, model)
end
function move!(sheep::Sheep,model)
wolves = filter(x -> isa(x, Wolf), collect(sheep.nearby_agents))
if !isempty(wolves)
closest_wolf = findmin(wolf -> sqrt(sum((sheep.pos .- wolf.pos) .^ 2)), wolves)[2]
move_away!(sheep, wolves[closest_wolf].pos, model)
elseif !isempty(sheep.nearby_grass)
pos = random_empty_fully_grown(sheep.nearby_grass, model)
move_towards!(sheep, pos, model)
else
randomwalk!(sheep, model)
end
sheep.energy -= 1
end
function eat!(sheep::Sheep, model)
if model.fully_grown[sheep.pos...]
sheep.energy += model.Δenergy_sheep#sheep.Δenergy
model.fully_grown[sheep.pos...] = false
end
return
end
function reproduce!(sheep::Sheep, model)
print(model.sheep_reproduce)
if rand(abmrng(model)) model.sheep_reproduce#sheep.reproduction_prob
sheep.energy /= 2
replicate!(sheep, model)
end
end
function Agents.agent2string(agent::Sheep)
"""
Sheep
ID = $(agent.id)
energy = $(agent.energy)
"""
end
function move_away!(agent, pos, model)
direction = agent.pos .- pos
direction = clamp.(direction,-1,1)
walk!(agent,direction,model)
end
function move_towards!(agent, pos, model; ifempty=true)
direction = pos .- agent.pos
direction = clamp.(direction,-1,1)
walk!(agent,direction,model; ifempty=ifempty)
end
function nearby_fully_grown(sheep::Sheep, model)
nearby_pos = nearby_positions(sheep.pos, model, sheep.perception)
fully_grown_positions = filter(x -> model.fully_grown[x...], collect(nearby_pos))
return fully_grown_positions
end
function random_empty_fully_grown(positions, model)
n_attempts = 2*length(positions)
while n_attempts != 0
pos_choice = rand(positions)
isempty(pos_choice, model) && return pos_choice
n_attempts -= 1
end
return positions[1]
end
@agent struct Wolf(GridAgent{2})
energy::Float64
reproduction_prob::Float64
Δenergy::Float64
perception::Int32
nearby_agents
#speed::Float64
#endurance::Float64
end
function perceive!(wolf::Wolf,model)
wolf.nearby_agents = nearby_agents(wolf, model, model.wolf_perception)#wolf.perception)
end
function move!(wolf::Wolf,model)
sheeps = filter(x -> isa(x, Sheep), collect(wolf.nearby_agents))
if !isempty(sheeps)
closest_sheep = findmin(sheep -> sqrt(sum((wolf.pos .- sheep.pos) .^ 2)), sheeps)[2]
move_towards!(wolf, sheeps[closest_sheep].pos, model; ifempty=false)
else
randomwalk!(wolf, model; ifempty=false)
end
wolf.energy -= 1
end
function eat!(wolf::Wolf, model)
dinner = first_sheep_in_position(wolf.pos, model)
if !isnothing(dinner)
remove_agent!(dinner, model)
wolf.energy += model.Δenergy_wolf#wolf.Δenergy
end
end
function reproduce!(wolf::Wolf, model)
if rand(abmrng(model)) model.wolf_reproduce#wolf.reproduction_prob
wolf.energy /= 2
replicate!(wolf, model)
end
end
function Agents.agent2string(agent::Wolf)
"""
Wolf
ID = $(agent.id)
energy = $(agent.energy)
"""
end
function first_sheep_in_position(pos, model)
ids = ids_in_position(pos, model)
j = findfirst(id -> model[id] isa Sheep, ids)
isnothing(j) ? nothing : model[ids[j]]::Sheep
end
# The function `initialize_model` returns a new model containing sheep, wolves, and grass
# using a set of pre-defined values (which can be overwritten). The environment is a two
# dimensional grid space, which enables animals to walk in all
# directions.
function initialize_model(;
events = [],
n_sheep = 100,
n_wolves = 50,
dims = (20, 20),
regrowth_time = 30,
Δenergy_sheep = 4,
Δenergy_wolf = 20,
sheep_reproduce = 0.04,
wolf_reproduce = 0.05,
sheep_perception = 0,
wolf_perception = 0,
seed = 23182,
)
rng = MersenneTwister(seed)
space = GridSpace(dims, periodic = true)
## Model properties contain the grass as two arrays: whether it is fully grown
## and the time to regrow. Also have static parameter `regrowth_time`.
## Notice how the properties are a `NamedTuple` to ensure type stability.
## define as dictionary(mutable) instead of tuples(immutable) as per https://github.com/JuliaDynamics/Agents.jl/issues/727
properties = Dict(
:events => events,
:fully_grown => falses(dims),
:countdown => zeros(Int, dims),
:regrowth_time => regrowth_time,
:Δenergy_sheep => Δenergy_sheep,
:Δenergy_wolf => Δenergy_wolf,
:sheep_reproduce => sheep_reproduce,
:wolf_reproduce => wolf_reproduce,
:sheep_perception => sheep_perception,
:wolf_perception => wolf_perception
)
model = StandardABM(Union{Sheep, Wolf}, space;
agent_step! = sheepwolf_step!, model_step! = custom_model_step!,
properties, rng, scheduler = Schedulers.Randomly(), warn = false
)
## Add agents
for _ in 1:n_sheep
energy = rand(abmrng(model), 1:(Δenergy_sheep*2)) - 1
add_agent!(Sheep, model, energy, sheep_reproduce, Δenergy_sheep, sheep_perception, [], [])
end
for _ in 1:n_wolves
energy = rand(abmrng(model), 1:(Δenergy_wolf*2)) - 1
add_agent!(Wolf, model, energy, wolf_reproduce, Δenergy_wolf, wolf_perception, [])
end
## Add grass with random initial growth
for p in positions(model)
fully_grown = rand(abmrng(model), Bool)
countdown = fully_grown ? regrowth_time : rand(abmrng(model), 1:regrowth_time) - 1
model.countdown[p...] = countdown
model.fully_grown[p...] = fully_grown
end
return model
end
# ## Defining the stepping functions
# Sheep and wolves behave similarly:
# both lose 1 energy unit by moving to an adjacent position and both consume
# a food source if available. If their energy level is below zero, they die.
# Otherwise, they live and reproduce with some probability.
# They move to a random adjacent position with the [`randomwalk!`](@ref) function.
# Notice how the function `sheepwolf_step!`, which is our `agent_step!`,
# is dispatched to the appropriate agent type via Julia's Multiple Dispatch system.
function sheepwolf_step!(sheep::Sheep, model)
perceive!(sheep, model)
move!(sheep, model)
if sheep.energy < 0
remove_agent!(sheep, model)
return
end
eat!(sheep, model)
reproduce!(sheep, model)
end
function sheepwolf_step!(wolf::Wolf, model)
perceive!(wolf, model)
move!(wolf, model)
if wolf.energy < 0
remove_agent!(wolf, model)
return
end
eat!(wolf, model)
reproduce!(wolf, model)
end
function custom_model_step!(model)
event_handler!(model)
grass_step!(model)
end
# The behavior of grass function differently. If it is fully grown, it is consumable.
# Otherwise, it cannot be consumed until it regrows after a delay specified by
# `regrowth_time`. The dynamics of the grass is our `model_step!` function.
function grass_step!(model)
@inbounds for p in positions(model) # we don't have to enable bound checking
if !(model.fully_grown[p...])
if model.countdown[p...] 0
model.fully_grown[p...] = true
model.countdown[p...] = model.regrowth_time
else
model.countdown[p...] -= 1
end
end
end
end
# Check current step and start event at step t
function event_handler!(model)
for event in model.events
if event.timer == event.t_start # start event
if event.name == "Drought"
model.regrowth_time = event.value
model.wolf_perception += 1
model.sheep_perception += 1
elseif event.name == "Flood"
model.regrowth_time = event.value
model.Δenergy_wolf = model.Δenergy_wolf - 1
model.Δenergy_sheep = model.Δenergy_sheep - 1
elseif event.name == "PreyReproduceSeasonal"
model.sheep_reproduce = event.value
elseif event.name == "PredatorReproduceSeasonal"
model.wolf_reproduce = event.value
end
end
if event.timer == event.t_end # end event
if event.name == "Drought"
model.regrowth_time = event.pre_value
model.wolf_perception -= 1
model.sheep_perception -= 1
elseif event.name == "Flood"
model.regrowth_time = event.pre_value
model.Δenergy_wolf = model.Δenergy_wolf + 1
model.Δenergy_sheep = model.Δenergy_sheep + 1
elseif event.name == "PreyReproduceSeasonal"
model.sheep_reproduce = event.pre_value
elseif event.name == "PredatorReproduceSeasonal"
model.wolf_reproduce = event.pre_value
end
end
if event.timer == event.t_cycle # reset cycle
event.timer = 1
else
event.timer += 1
end
end
end
mutable struct RecurringEvent
name::String
value::Float64
pre_value::Float64
t_start::Int64
t_end::Int64
t_cycle::Int64
timer::Int64
end