Merge pull request #263 from mrc-ide/antimalarial_resistance/etf

Antimalarial resistance/etf

NAMESPACE

@@ -18,6 +18,7 @@ export(rtss_profile)
 export(run_metapop_simulation)
 export(run_simulation)
 export(run_simulation_with_repetitions)
+export(set_antimalarial_resistance)
 export(set_bednets)
 export(set_carrying_capacity)
 export(set_clinical_treatment)

R/antimalarial_resistance.R

@@ -0,0 +1,139 @@
+#' @title Parameterise antimalarial resistance
+#' @description
+#' Parameterise antimalarial resistance
+#' 
+#' @param parameters the model parameters
+#' @param drug the index of the drug which resistance is being set, as set by the set_drugs() function, in the parameter list
+#' @param timesteps vector of time steps for each update to resistance proportion and resistance outcome probability
+#' @param artemisinin_resistance vector of updates to the proportions of infections that are artemisinin resistant at time t
+#' @param partner_drug_resistance vector of updates to the proportions of infections that are partner-drug resistant at time t
+#' @param slow_parasite_clearance_prob vector of updates to the proportion of artemisinin-resistant infections that result in early treatment failure
+#' @param early_treatment_failure_prob vector of updates to the proportion of artemisinin-resistant infections that result in slow parasite clearance
+#' @param late_clinical_failure_prob vector of updates to the proportion of partner-drug-resistant infections that result in late clinical failure
+#' @param late_parasitological_prob vector of updates to the proportion of partner-drug-resistant infections that result in late parasitological failure
+#' @param reinfection_prob vector of updates to the proportion of partner-drug-resistant infections that result in reinfection during prophylaxis
+#' @param slow_parasite_clearance_time single value representing the mean time individual's experiencing slow parasite clearance reside in the treated state
+#' @export
+set_antimalarial_resistance <- function(parameters,
+                                        drug,
+                                        timesteps,
+                                        artemisinin_resistance,
+                                        partner_drug_resistance,
+                                        slow_parasite_clearance_prob,
+                                        early_treatment_failure_prob,
+                                        late_clinical_failure_prob,
+                                        late_parasitological_prob,
+                                        reinfection_prob,
+                                        slow_parasite_clearance_time) {
+
+  if(any(c(length(artemisinin_resistance),
+           length(partner_drug_resistance), 
+           length(slow_parasite_clearance_prob), 
+           length(early_treatment_failure_prob), 
+           length(late_clinical_failure_prob), 
+           length(late_parasitological_prob), 
+           length(reinfection_prob)) != length(timesteps))) {
+    stop("Length of one or more resistance parameter vectors does not match time steps specified for update")
+  }
+
+  if(any(artemisinin_resistance < 0 | artemisinin_resistance > 1 |
+         partner_drug_resistance < 0 | partner_drug_resistance > 1)) {
+    stop("Artemisinin and partner-drug resistance proportions must fall between 0 and 1")
+  }
+
+  if(any(slow_parasite_clearance_prob < 0 | slow_parasite_clearance_prob > 1 |
+         early_treatment_failure_prob < 0 | early_treatment_failure_prob > 1 |
+         late_clinical_failure_prob < 0 | late_clinical_failure_prob > 1 |
+         late_parasitological_prob < 0 | late_parasitological_prob > 1 |
+         reinfection_prob < 0 | reinfection_prob > 1)) {
+    stop("Resistance outcome probabilities must fall between 0 and 1")
+  }
+
+    if(length(slow_parasite_clearance_time) != 1) {
+    stop("Error: length of slow_parasite_clearance_time not equal to 1")
+  }
+
+  if(slow_parasite_clearance_time <= 0) {
+    stop("Error: slow_parasite_clearance_time is non-positive")
+  }
+
+  parameters$antimalarial_resistance <- TRUE
+
+  n_drugs <- length(parameters$drug_efficacy)
+
+  if (drug < 1 | drug > n_drugs) {
+    stop('Drug index is invalid, please set drugs using set_drugs')
+  }
+
+  drug_index <- which(parameters$antimalarial_resistance_drug == drug)
+
+  if (length(drug_index) == 0) {
+    drug_index <- length(parameters$antimalarial_resistance_drug) + 1
+  }
+
+  parameters$antimalarial_resistance_drug[[drug_index]] <- drug
+  parameters$antimalarial_resistance_timesteps[[drug_index]] <- timesteps
+  parameters$prop_artemisinin_resistant[[drug_index]] <- artemisinin_resistance
+  parameters$prop_partner_drug_resistant[[drug_index]] <- partner_drug_resistance
+  parameters$slow_parasite_clearance_prob[[drug_index]] <- slow_parasite_clearance_prob
+  parameters$early_treatment_failure_prob[[drug_index]] <- early_treatment_failure_prob
+  parameters$late_clinical_failure_prob[[drug_index]] <- late_clinical_failure_prob
+  parameters$late_parasitological_failure_prob[[drug_index]] <- late_parasitological_prob
+  parameters$reinfection_during_prophylaxis[[drug_index]] <- reinfection_prob
+  parameters$dt_slow_parasite_clearance[[drug_index]] <- slow_parasite_clearance_time
+
+  return(parameters)
+
+}
+
+#' @title Retrieve resistance parameters
+#' @description
+#' Retrieve the resistance parameters associated with the drug each individual receiving clinical 
+#' treatment has been administered in the current time step.
+#' 
+#' @param parameters the model parameters
+#' @param drugs vector of integers representing the drugs administered to each individual receiving treatment
+#' @param timestep the current time step
+get_antimalarial_resistance_parameters <- function(parameters, drugs, timestep) {
+
+  if(!parameters$antimalarial_resistance) {
+    stop("Error: Antimalarial resistance has not been parameterised; antimalarial_resistance = FALSE")
+  }
+
+  blank_vector <- numeric(length = length(drugs))
+  artemisinin_resistance_proportion <- blank_vector
+  partner_drug_resistance_proportion <- blank_vector
+  slow_parasite_clearance_probability <- blank_vector
+  early_treatment_failure_probability <- blank_vector
+  late_clinical_failure_probability <- blank_vector
+  late_parasitological_failure_probability <- blank_vector
+  reinfection_during_prophylaxis_probability <- blank_vector
+  dt_slow_parasite_clearance <- rep(parameters$dt, length = length(drugs))
+
+  for(i in seq_along(parameters$antimalarial_resistance_drug)) {
+    drug <- parameters$antimalarial_resistance_drug[[i]]
+    treated_with_drug <- which(drugs == drug)
+    resistance_timestep <- match_timestep(ts = parameters$antimalarial_resistance_timesteps[[i]], t = timestep)
+    artemisinin_resistance_proportion[treated_with_drug] <- parameters$prop_artemisinin_resistant[[i]][resistance_timestep]
+    partner_drug_resistance_proportion[treated_with_drug] <- parameters$prop_partner_drug_resistant[[i]][resistance_timestep]
+    slow_parasite_clearance_probability[treated_with_drug] <- parameters$slow_parasite_clearance_prob[[i]][resistance_timestep]
+    early_treatment_failure_probability[treated_with_drug] <- parameters$early_treatment_failure_prob[[i]][resistance_timestep]
+    late_clinical_failure_probability[treated_with_drug] <- parameters$late_clinical_failure_prob[[i]][resistance_timestep]
+    late_parasitological_failure_probability[treated_with_drug] <- parameters$late_parasitological_failure_prob[[i]][resistance_timestep]
+    reinfection_during_prophylaxis_probability[treated_with_drug] <- parameters$reinfection_during_prophylaxis[[i]][resistance_timestep]
+    dt_slow_parasite_clearance[treated_with_drug] <- parameters$dt_slow_parasite_clearance[[i]]
+  }
+
+  resistance_parameters <- list()
+  resistance_parameters$artemisinin_resistance_proportion <- artemisinin_resistance_proportion
+  resistance_parameters$partner_drug_resistance_proportion <- partner_drug_resistance_proportion
+  resistance_parameters$slow_parasite_clearance_probability <- slow_parasite_clearance_probability
+  resistance_parameters$early_treatment_failure_probability <- early_treatment_failure_probability
+  resistance_parameters$late_clinical_failure_probability <- late_clinical_failure_probability
+  resistance_parameters$late_parasitological_failure_probability <- late_parasitological_failure_probability
+  resistance_parameters$reinfection_during_prophylaxis_probability <- reinfection_during_prophylaxis_probability
+  resistance_parameters$dt_slow_parasite_clearance <- dt_slow_parasite_clearance
+
+  return(resistance_parameters)
+
+}

R/human_infection.R

@@ -263,23 +263,27 @@ update_severe_disease <- function(
 #' @param renderer simulation renderer
 #' @noRd
 calculate_treated <- function(
-  variables,
-  clinical_infections,
-  parameters,
-  timestep,
-  renderer
-  ) {
+    variables,
+    clinical_infections,
+    parameters,
+    timestep,
+    renderer
+) {
+
+  if(clinical_infections$size() == 0) {
+    return(individual::Bitset$new(parameters$human_population))
+  }
+
   treatment_coverages <- get_treatment_coverages(parameters, timestep)
   ft <- sum(treatment_coverages)
-
+  
   if (ft == 0) {
     return(individual::Bitset$new(parameters$human_population))
   }
-
+  
   renderer$render('ft', ft, timestep)
   seek_treatment <- sample_bitset(clinical_infections, ft)
   n_treat <- seek_treatment$size()
-
   renderer$render('n_treated', n_treat, timestep)
 
   drugs <- as.numeric(parameters$clinical_treatment_drugs[
@@ -290,29 +294,51 @@ calculate_treated <- function(
       replace = TRUE
     )
   ])
-
-  successful <- bernoulli_multi_p(parameters$drug_efficacy[drugs])
-  treated_index <- bitset_at(seek_treatment, successful)
-
-  # Update those who have been treated
-  if (treated_index$size() > 0) {
-    variables$state$queue_update('Tr', treated_index)
+
+  if(parameters$antimalarial_resistance) {
+    resistance_parameters <- get_antimalarial_resistance_parameters(
+      parameters = parameters,
+      drugs = drugs, 
+      timestep = timestep
+    )
+    unsuccessful_treatment_probability <- resistance_parameters$artemisinin_resistance_proportion * resistance_parameters$early_treatment_failure_probability
+    susceptible_to_treatment_index <- bernoulli_multi_p(p = 1 - unsuccessful_treatment_probability)
+    susceptible_to_treatment <- bitset_at(seek_treatment, susceptible_to_treatment_index)
+    drugs <- drugs[susceptible_to_treatment_index]
+  } else {
+    susceptible_to_treatment <- seek_treatment
+
+  }
+
+  n_early_treatment_failure <- n_treat - susceptible_to_treatment$size()
+  successfully_treated_index <- bernoulli_multi_p(parameters$drug_efficacy[drugs])
+  successfully_treated <- bitset_at(susceptible_to_treatment, successfully_treated_index)
+  successfully_treated_drugs <- drugs[successfully_treated_index]
+  n_treat_eff_fail <- susceptible_to_treatment$size() - length(successfully_treated_index)
+  renderer$render('n_early_treatment_failure', n_early_treatment_failure, timestep)
+  renderer$render('n_treat_eff_fail', n_treat_eff_fail, timestep)
+  renderer$render('n_treat_success', successfully_treated$size(), timestep)
+
+  # Update variables of those who have been successfully treated:
+  if (successfully_treated$size() > 0) {
+    variables$state$queue_update("Tr", successfully_treated)
     variables$infectivity$queue_update(
-      parameters$cd * parameters$drug_rel_c[drugs[successful]],
-      treated_index
+      parameters$cd * parameters$drug_rel_c[successfully_treated_drugs],
+      successfully_treated
     )
     variables$drug$queue_update(
-      drugs[successful],
-      treated_index
+      successfully_treated_drugs,
+      successfully_treated
     )
     variables$drug_time$queue_update(
       timestep,
-      treated_index
+      successfully_treated
     )
   }
-  treated_index
+  successfully_treated
 }
 
+
 #' @title Schedule infections
 #' @description
 #' Schedule infections in humans after the incubation period

R/model.R

@@ -74,6 +74,9 @@
 #' susceptible
 #'  * net_usage: the number people protected by a bed net
 #'  * mosquito_deaths: number of adult female mosquitoes who die this timestep
+#'  * n_early_treatment_failure: number of clinically treated individuals who experienced early treatment failure in this timestep
+#'  * n_treat_eff_fail: number of clinically treated individuals who's treatment failed due to drug efficacy
+#'  * n_treat_success: number of successfully treated individuals in this timestep
 #'
 #' @param timesteps the number of timesteps to run the simulation for (in days)
 #' @param parameters a named list of parameters to use

R/parameters.R

@@ -192,6 +192,22 @@
 #' * tbv_tra_mu - transmission reduction parameter; default = 12.63
 #' * tbv_gamma1 - transmission reduction parameter; default = 2.5
 #' * tbv_gamma2 - transmission reduction parameter; default = 0.06
+#' 
+#' Antimalarial resistance parameters:
+#' please set antimalarial resistance parameters with the convenience functions in
+#' `antimalarial_resistance.R:set_antimalarial_resistance`
+#' 
+#' * antimalarial_resistance - boolean for if antimalarial resistance is enabled; default = FALSE
+#' * antimalarial_resistance_drug - vector of drugs for which resistance can be parameterised; default = NULL
+#' * antimalarial_resistance_timesteps - vector of time steps on which resistance updates occur; default = NULL
+#' * prop_artemisinin_resistant - vector of proportions of infections resistant to the artemisinin component of a given drug; default = NULL
+#' * prop_partner_drug_resistant - vector of proportions of infections resistant to the parter drug component of a given drug; default = NULL
+#' * slow_parasite_clearance_prob - vector of probabilities of slow parasite clearance for a given drug; default = NULL
+#' * early_treatment_failure_prob - vector of probabilities of early treatment failure for a given drug; default = NULL
+#' * late_clinical_failure_prob - vector of probabilities of late clinical failure for a given drug; default = NULL
+#' * late_parasitological_failure_prob - vector of probabilities of late parasitological failure for a given drug; default = NULL
+#' * reinfection_during_prophylaxis - vector of probabilities of reinfection during prophylaxis for a given drug; default = NULL
+#' * dt_slow_parasite_clearance - the delay for humans experiencing slow parasite clearance to move from state Tr to S; default = NULL
 #'
 #' rendering:
 #' All values are in timesteps and all ranges are inclusive
@@ -377,6 +393,17 @@ get_parameters <- function(overrides = list()) {
     tbv_timesteps = NULL,
     tbv_coverages = NULL,
     tbv_ages = NULL,
+    # antimalarial resistance
+    antimalarial_resistance = FALSE,
+    antimalarial_resistance_drug = NULL,
+    antimalarial_resistance_timesteps = NULL,
+    prop_artemisinin_resistant = NULL,
+    prop_partner_drug_resistant = NULL,
+    slow_parasite_clearance_prob = NULL,
+    early_treatment_failure_prob = NULL,
+    late_clinical_failure_prob = NULL,
+    late_parasitological_failure_prob = NULL,
+    reinfection_during_prophylaxis = NULL,
     # flexible carrying capacity
     carrying_capacity = FALSE,
     carrying_capacity_timesteps = NULL,