# Example

## Gaussian random walk

This is a simple implementation to simulate a 2-D Gaussian random walk.

```
using RecordedArrays
using Plots
using Random
Random.seed!(1)
c = DiscreteClock(10000) # define a clock, the particle will walk 10000 epoch
pos = recorded(DynamicEntry, c, [0.0, 0.0]) # create a position vector of the particle
for _ in c
pos .+= randn(2) # walk randomly at each epoch
end
# plot path of particle
phaseportrait(pos; frame=:none, grid=false, legend=false)
```

## Logistic growth

This is a simple implementation of Gillespie algorithm with direct method to simulate a Logistic growth population with growth rate $r=0.5$ and carrying capacity $K=100$.

```
using RecordedArrays
using Plots
using Random
Random.seed!(1)
c = ContinuousClock(100.0) # define a clock, the population will growth for 100 time unit
n = recorded(DynamicEntry, c, 10) # define a scalar to record population size
const r = 0.5
const K = 100
for _ in c
# evaluate a_i
grow = r * n # intrinsic growth
comp = r * n * n / K # resource competition
summed = grow + comp # sum a_i
τ = -log(rand()) / summed # compute time interval
increase!(c, τ) # update current time
# sample a reaction and adjust population size
if rand() * summed < grow
n[1] += 1
else
n[1] -= 1
end
state(n) <= 0 && break # break if population extinct
end
# plot population dynamics
timeseries(n; frame=:box, grid=false, legend=false)
```

## Stochastic Predator–prey Dynamics

This is a simple implementation of Gillespie algorithm with direct method to simulate a Predator–prey Dynamics.

```
using RecordedArrays
using Plots
using Random
Random.seed!(1)
c = ContinuousClock(100.0) # define a clock, the population will growth for 100 time unit
n = recorded(DynamicEntry, c, [100, 100]) # define a vector to record population size (n[1] for prey, n[2] for predator)
const α = 0.5
const β = 0.001
const δ = 0.001
const γ = 0.5
for _ in c
n[2] == 0 && break
# evaluate a_i
grow = α * n[1] # intrinsic growth of prey
predation_prey = β * n[1] * n[2] # predation cause death of prey
predation_pred = δ * n[1] * n[2] # predation cause reproduction of predator
death = γ * n[2] # intrinsic death of prey
summed = grow + predation_prey + predation_pred + death
summed <= 0 && break
τ = -log(rand()) / summed # compute time interval
increase!(c, τ) # update current time
# sample a reaction and adjust population size
r1 = rand() * summed
if (r1 -= grow; r1 < 0)
n[1] += 1
elseif (r1 -= predation_prey; r1 < 0)
n[1] -= 1
elseif (r1 -= predation_pred; r1 < 0)
n[2] += 1
else
n[2] -= 1
end
end
ts = timeseries(n; frame=:box, grid=false, legend=false)
pp = phaseportrait(n; frame=:box, grid=false, legend=false)
plot(ts, pp; titles=["Population dynamics" "Phase portrait"], size=(800, 400))
```