Housing

NAMESPACE

@@ -24,6 +24,7 @@ export(set_clinical_treatment)
 export(set_demography)
 export(set_drugs)
 export(set_equilibrium)
+export(set_housing)
 export(set_mass_pev)
 export(set_mda)
 export(set_parameter_draw)

R/correlation.R

@@ -5,7 +5,8 @@ INTS <- c(
   'smc',
   'tbv',
   'bednets',
-  'spraying'
+  'spraying',
+  'housing'
 )
 
 #' Class: Correlation parameters

R/parameters.R

@@ -149,7 +149,9 @@
 #' * k0 - proportion of females bloodfed with no net; default = 0.699
 #' * spraying - boolean for if indoor spraying is enabled; default = FALSE
 #' * phi_indoors - proportion of bites taken indoors; default = 0.90
-#'
+#' * housing - boolean for if housing improvements are enabled; default = FALSE
+#' * phi_housing - fractional increase to outdoor biting due to improved housing; default = 1 simulating no increase
+#' 
 #' treatment parameters:
 #' please set treatment parameters with the convenience functions in
 #' `drug_parameters.R`
@@ -330,6 +332,8 @@ get_parameters <- function(overrides = list()) {
     # indoor spraying
     spraying = FALSE,
     phi_indoors = .90,
+    # housing improvements
+    housing = FALSE,
     # treatment
     drug_efficacy          = numeric(0),
     drug_rel_c             = numeric(0),

R/processes.R

@@ -235,6 +235,18 @@ create_processes <- function(
     )
   }
 
+  if (parameters$housing) {
+    processes <- c(
+      processes,
+      housing_improvement(
+        variables,
+        parameters,
+        correlations
+      ),
+      house_usage_renderer(variables$house_time, renderer)
+    )
+  }
+
   # ======================
   # Progress bar process
   # ======================

R/variables.R

@@ -26,6 +26,7 @@
 #' haven't been any
 #' * net_time - The timestep when a net was last put up (-1 if never)
 #' * spray_time - The timestep when the house was last sprayed (-1 if never)
+#' * house_time - The timestep for adaptation to the house to reduce entry/kill vectors (-1 if never)
 #' * infectivity - The onward infectiousness to mosquitos
 #' * drug - The last prescribed drug
 #' * drug_time - The timestep of the last drug
@@ -196,6 +197,7 @@ create_variables <- function(parameters) {
   # Init vector controls
   net_time <- individual::IntegerVariable$new(rep(-1, size))
   spray_time <- individual::IntegerVariable$new(rep(-1, size))
+  house_time <- individual::IntegerVariable$new(rep(-1, size))
 
   variables <- list(
     state = state,
@@ -219,7 +221,8 @@ create_variables <- function(parameters) {
     pev_profile = pev_profile,
     tbv_vaccinated = tbv_vaccinated,
     net_time = net_time,
-    spray_time = spray_time
+    spray_time = spray_time,
+    house_time = house_time
   )
 
   # Add variables for individual mosquitoes

R/vector_control.R

@@ -11,7 +11,7 @@ prob_bitten <- function(
   parameters
   ) {
   n <- parameters$human_population
-  if (!(parameters$bednets || parameters$spraying)) {
+  if (!(parameters$bednets || parameters$spraying || parameters$housing)) {
     return(
       list(
         prob_bitten_survives = rep(1, n),
@@ -65,37 +65,84 @@ prob_bitten <- function(
       js_prime,
       parameters$k0
     )
+    spray_on = 1
     rs_comp <- 1 - rs
     ss <- rep(1, n)
     ss[protected_index] <- prob_survives_spraying(
       ks_prime,
       parameters$k0
     )
   } else {
-    phi_indoors <- 0
+    spray_on = 0
     rs <- 0
     rs_comp <- 1
     ss <- 1
   }
+
+  if (parameters$housing) {
+    phi_housing <- parameters$phi_housing[[species]]  ## if a change in behaviour is caused by the housing changes, we can increase outdoor biting proportion
+    phi_indoors <- parameters$phi_indoors[[species]]
+    house_time <- variables$house_time$get_values()
+    since_housing <- timestep - house_time
+    matches <- match(house_time, parameters$housing_timesteps)
+    rh <- prob_repelled_house(matches, since_housing, species, parameters) ## if housing prevents entry to house, we increase proportion needing to repeat foraging
+    sh <- prob_survives_house(rh, matches, since_housing, species, parameters)
+    unused <- house_time == -1
+    sh[unused] <- 1
+    rh[unused] <- 0
+  } else {
+    phi_housing <- 1
+    rh <- 0
+    sh <- 1
+  }
 
+  if ((!parameters$housing & !parameters$spraying)) {
+    phi_indoors <- 0 ## we want phi_indoors to be applied if housing is on
+  }
+
   list(
     prob_bitten_survives = (
-      1 - phi_indoors +
-      phi_bednets * rs_comp * sn * ss +
-      (phi_indoors - phi_bednets) * rs_comp * ss
+      1 - phi_indoors * phi_housing + ##
+        (1 - rh) * (phi_bednets * phi_housing * rs_comp * sn * ss * sh^2) +  ## * sh if some mortality from housing
+        (1 - rh) * ((phi_indoors * phi_housing - phi_bednets * phi_housing) * rs_comp * ss * sh^2) ## * sh if some mortality from housing
     ),
     prob_bitten = (
-      1 - phi_indoors +
-      phi_bednets * rs_comp * sn +
-      (phi_indoors - phi_bednets) * rs_comp
+      1 - phi_indoors * phi_housing +
+        (1 - rh) * (phi_bednets * phi_housing * rs_comp * sn * sh) + ## * sh if some mortality from housing
+        (1 - rh) * ((phi_indoors * phi_housing - phi_bednets * phi_housing) * rs_comp * sh) ## * sh if some mortality from housing
     ),
     prob_repelled = (
-      phi_bednets * rs_comp * rn +
-      phi_indoors * rs
+      phi_bednets * phi_housing * rs_comp * sh * rn + 
+        phi_indoors * phi_housing * rs * sh * spray_on + 
+        phi_indoors * phi_housing * rh
     )
   )
 }
 
+#' @title Simulate housing improvements
+#' @description simulates improved housing so that harder for vectors to get indoors
+#' from `set_housing` and correlation parameters from
+#' `get_correlation_parameters`
+#'
+#' @param variables list of variables in the model
+#' @param parameters the model parameters
+#' @param correlations correlation parameters
+#' @noRd
+housing_improvement <- function(variables, parameters, correlations) {
+  function(timestep) {
+    matches <- timestep == parameters$housing_timesteps
+    if (any(matches)) {
+      target <- which(sample_intervention(
+        seq(parameters$human_population),
+        'housing',
+        parameters$housing_coverages[matches],
+        correlations
+      ))
+      variables$house_time$queue_update(timestep, target)
+    }
+  }
+}
+
 #' @title Indoor spraying
 #' @description models indoor residual spraying according to the strategy
 #' from `set_spraying` and correlation parameters from
@@ -106,7 +153,7 @@ prob_bitten <- function(
 #' @param correlations correlation parameters
 #' @noRd
 indoor_spraying <- function(spray_time, parameters, correlations) {
-  function(timestep) {
+      function(timestep) {
     matches <- timestep == parameters$spraying_timesteps
     if (any(matches)) {
       target <- which(sample_intervention(
@@ -187,6 +234,22 @@ spraying_decay <- function(t, theta, gamma) {
   1 / (1 + exp(-(theta + gamma * t)))
 }
 
+house_decay <- function(t, gamma) {
+  exp(-t / gamma)
+}
+
+prob_repelled_house <- function(matches, dt, species, parameters) {
+  rhm <- parameters$house_rhm[matches, species]
+  gammah <- parameters$house_gammah[matches]
+  (parameters$house_rh[matches, species] - rhm) * house_decay(dt, gammah) + rhm
+}
+
+prob_survives_house <- function(rh, matches, dt, species, parameters) {
+  dh0 <- parameters$house_dh0[matches, species]
+  dh <- dh0 * house_decay(dt, parameters$house_gammah[matches])
+  1 - rh - dh
+}
+
 net_usage_renderer <- function(net_time, renderer) {
   function(t) {
     renderer$render(
@@ -196,3 +259,13 @@ net_usage_renderer <- function(net_time, renderer) {
     )
   }
 }
+
+house_usage_renderer <- function(house_time, renderer) {
+  function(t) {
+    renderer$render(
+      'n_use_house_adaptation',
+      house_time$get_index_of(-1)$not(TRUE)$size(),
+      t
+    )
+  }
+}

R/vector_control_parameters.R

@@ -179,3 +179,61 @@ get_init_carrying_capacity <- function(parameters){
   names(init_cc) <- parameters$species
   return(init_cc)
 }
+
+#' @title Parameterise a housing improvement strategy
+#'
+#' @description The model will simulate improved housing at `timesteps` to a random
+#' sample of the entire human population. The sample size will be a proportion
+#' of the human population taken from the corresponding `coverages`.
+#' The sample _can_ contain humans who have already benefited from housing.
+#'
+#' If a human in the sample lives in a house that has been improved to reduce
+#' mosquito entry, the efficacy of the
+#' housing improvement will be high - as determined by the parameter rn_house
+#'
+#' The structure for the housing model will be documented in a publication
+#' Sherrard-Smith et al in prep
+#'
+#' @param parameters a list of parameters to modify
+#' @param timesteps the timesteps at which to distribute housing adaptations
+#' @param coverages the proportion of the human population who reside in protected housing
+#' @param dh0 a matrix of death probabilities for each species over time.
+#' With nrows=length(timesteps), ncols=length(species)
+#' @param rh a matrix of repelling probabilities for each species over time
+#' With nrows=length(timesteps), ncols=length(species)
+#' @param rhm a matrix of minimum repelling probabilities for each species over time
+#' With nrows=length(timesteps), ncols=length(species)
+#' @param gammah a vector of house adaptation half-lives for each distribution timestep
+#' @export
+set_housing <- function(
+    parameters,
+    timesteps,
+    coverages,
+    phi_housing,
+    dh0,
+    rh,
+    rhm,
+    gammah
+) {
+  stopifnot(all(coverages >= 0) && all(coverages <= 1))
+  lengths <- vnapply(list(coverages, gammah), length)
+  if (!all(lengths == length(timesteps))) {
+    stop('timesteps and time-varying parameters must align')
+  }
+  for (x in list(dh0, rh, rhm)) {
+    if (ncol(x) != length(parameters$species)) {
+      stop('death and repelling probabilities rows need to align with species')
+    }
+    if (nrow(x) != length(timesteps)) {
+      stop('death and repelling probabilities columns need to align with timesteps')
+    }
+  }
+  parameters$housing <- TRUE
+  parameters$housing_timesteps <- timesteps
+  parameters$housing_coverages <- coverages
+  parameters$house_dh0 <- dh0
+  parameters$house_rh <- rh
+  parameters$house_rhm <- rhm
+  parameters$house_gammah <- gammah
+  parameters
+}

R/vector_parameters.R

@@ -6,6 +6,7 @@
 #' Q0: 0.92
 #' phi_bednets: 0.85
 #' phi_indoors: 0.90
+#' phi_housing: 1
 #' mum: 0.132
 #'
 #' parameters from:
@@ -18,6 +19,7 @@ gamb_params <- list(
   Q0 = .92,
   phi_bednets = .85,
   phi_indoors = .90,
+  phi_housing = 1,
   mum = .132
 )
 
@@ -29,6 +31,7 @@ gamb_params <- list(
 #' Q0: 0.71
 #' phi_bednets: 0.8
 #' phi_indoors: 0.86
+#' phi_housing: 1
 #' mum: 0.132
 #'
 #' parameters from:
@@ -41,6 +44,7 @@ arab_params <- list(
   Q0 = .71,
   phi_bednets = .8,
   phi_indoors = .86,
+  phi_housing = 1,
   mum = .132
 )
 
@@ -52,6 +56,7 @@ arab_params <- list(
 #' Q0: 0.94
 #' phi_bednets: 0.78
 #' phi_indoors: 0.87
+#' phi_housing: 1
 #' mum: 0.112
 #'
 #' parameters from:
@@ -64,6 +69,7 @@ fun_params <- list(
   Q0 = .94,
   phi_bednets = .78,
   phi_indoors = .87,
+  phi_housing = 1,
   mum = .112
 )
 
@@ -75,6 +81,7 @@ fun_params <- list(
 #' Q0: 0.21
 #' phi_bednets: 0.57
 #' phi_indoors: 0.37
+#' phi_housing: 1
 #' mum: 0.112
 #'
 #' parameters reference:
@@ -91,6 +98,7 @@ steph_params <- list(
   Q0 = 0.21,
   phi_bednets = 0.52186,
   phi_indoors = 0.4776,
+  phi_housing = 1,
   mum = .112
 )
 
@@ -114,6 +122,7 @@ set_species <- function(parameters, species, proportions) {
     'Q0',
     'phi_bednets',
     'phi_indoors',
+    'phi_housing',
     'mum'
   )
   for (key in keys) {