Scaling for Malaria Metrics

October 14, 2024

Load the required packages:

library(ramp.xds)
library(ramp.work)
library(deSolve)

An important question in malaria is the relationship between various metrics, especially the relationship between the annual entomological inoculation rate (aEIR) and the parasite rate (PR). To facilitate analysis of malaria data, we have developed some functions to compute scaling relationships.

xde_scaling

The function xde_scaling() defines the relationship between the EIR and the PR. It analyzes stable orbits and outputs the average annual EIR and average annual PR for an even mesh of log(aEIR) values running from $10^{-1}$ up to $10^{3}$ The code is in mob_library/Work

To illustrate, we pick a function describing a seasonal pattern using ramp.xds::make_F_sin

tt <- seq(0, 730, by=5) 
p1 <- makepar_F_sin(floor=0.1)
Fsin <- make_function(p1)
plot(tt, Fsin(tt), type="l")

Next, we set up a cohort model:

xds_setup_cohort(Xname = "SIS", F_season=Fsin) -> sis

xds_solve_cohort(sis) -> sis

The function xde_scaling_eir runs the model over a mesh of N=25 values:

xde_scaling_eir(sis, 25) -> sis

The results are attached as sis$outputs$eirpr

plot_eirpr(sis)

require(viridis)

## Loading required package: viridis

## Loading required package: viridisLite

clrs = turbo(25)

plot_eirpr(sis)

with(sis$output$eirpr,{
  points(aeir, pr, col = clrs)
  lines(scaling[[5]]$aeir, scaling[[5]]$pr, col = clrs[5])
  lines(scaling[[10]]$aeir, scaling[[10]]$pr, col = clrs[10])
  lines(scaling[[15]]$aeir, scaling[[15]]$pr, col = clrs[15])
  lines(scaling[[20]]$aeir, scaling[[20]]$pr, col = clrs[20])
})

xde_pr2eir()

Since xde_scaling defines the relationship between the EIR and the PR, we can now run xde_pr2eir() to get the predicted value of the eir, for any given value of the pr. The code is in mob_library/Work

We can run this for 50 randomly chosen values of the PfPR.

preir_i = xde_pr2eir(c(0.001, runif(25, 0, 1), 0.999), sis)

The function flags any values that are outside of the acceptable range. This may not seem important for the SIS model, but the range of other models can be bounded, so we don’t want to return nonsense values.

preir_i$errors

##         pr1         pr2         pr3 
## 0.001000000 0.007399441 0.999000000

We can plot the others:

plot_eirpr(sis)
with(sis$outputs$eirpr, points(aeir, pr, pch = 15))
with(preir_i, points(365*eir, pr, pch = 19, col = "red"))