In pahanc/malariasimulation_Ace_params: An individual based model for malaria

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

suppressPackageStartupMessages(library(ggplot2))
library(malariasimulation)
library(malariaEquilibrium)
library(reshape2)

Parameterisation

We are going to set the default parameters to run the simulation from an equilibrium.

year <- 365
month <- 30
sim_length <- 3 * year
human_population <- 1000
starting_EIR <- 50

simparams <- get_parameters(
  list(
    human_population = human_population,
    model_seasonality = TRUE, # Let's try a bi-modal model
    g0 = 0.28605,
    g = c(0.20636, -0.0740318, -0.0009293),
    h = c(0.173743, -0.0730962, -0.116019),
    severe_incidence_rendering_min_ages = 2*year,
    severe_incidence_rendering_max_ages = 10*year
  )
)

simparams <- set_equilibrium(simparams, starting_EIR)

# Plotting functions
plot_prevalence <- function(output) {
  ggplot(output) + geom_line(
    aes(x = timestep, y = (n_inc_severe_730_3650 / n_730_3650))) + 
    labs(x = "timestep", y = "Severe incidence")
}

add_intervention_lines <- function(plot, events) {
  plot + geom_vline(
    data = events,
    mapping = aes(xintercept=timestep),
    color="blue"
  ) + geom_text(
    data = events,
    mapping = aes(x = timestep, y = 0, label = name),
    size = 4,
    angle = 90,
    vjust = -0.4,
    hjust = 0
  )
}

Then we can run the simulation for a variety of vector control strategies:

Bed nets

We can distribute bed nets once a year for two years, changing the characteristics of the bed nets on each distribution...

bednetparams <- simparams

bednet_events = data.frame(
  timestep = c(1, 2) * year,
  name=c("Bednets 1", "Bednets 2")
)

bednetparams <- set_bednets(
  bednetparams,
  timesteps = bednet_events$timestep,
  coverages = c(.8, .8),
  retention = 5 * year,
  dn0 = matrix(c(.533, .45), nrow=2, ncol=1),
  rn = matrix(c(.56, .5), nrow=2, ncol=1),
  rnm = matrix(c(.24, .24), nrow=2, ncol=1),
  gamman = rep(2.64 * 365, 2)
)

output <- run_simulation(sim_length, bednetparams)

add_intervention_lines(plot_prevalence(output), bednet_events)

Indoor spraying

We can do the same for IRS...

sprayingparams <- simparams

peak <- peak_season_offset(sprayingparams)
spraying_events = data.frame(
  timestep = c(1, 2) * year + peak - 3 * month,
  name=c("Spraying 1", "Spraying 2")
)

sprayingparams <- set_spraying(
  sprayingparams,
  timesteps = spraying_events$timestep,
  coverages = rep(.8, 2),
  ls_theta = matrix(c(2.025, 2.025), nrow=2, ncol=1),
  ls_gamma = matrix(c(-0.009, -0.009), nrow=2, ncol=1),
  ks_theta = matrix(c(-2.222, -2.222), nrow=2, ncol=1),
  ks_gamma = matrix(c(0.008, 0.008), nrow=2, ncol=1),
  ms_theta = matrix(c(-1.232, -1.232), nrow=2, ncol=1),
  ms_gamma = matrix(c(-0.009, -0.009), nrow=2, ncol=1)
)

output <- run_simulation(sim_length, sprayingparams)

add_intervention_lines(plot_prevalence(output), spraying_events)