Merge pull request #138 from dd-harp/dev

pre-trip

NAMESPACE

@@ -429,7 +429,7 @@ export(xde_setup_cohort)
 export(xde_setup_human)
 export(xde_setup_mosy)
 export(xde_solve)
-export(xde_steady)
+export(xde_stable_orbit)
 importFrom(MASS,ginv)
 importFrom(deSolve,lagderiv)
 importFrom(deSolve,lagvalue)

R/adult-RM.R

@@ -183,7 +183,7 @@ make_MYZpar_RM = function(pars, MYZopts=list(), calK,
     pars = MBionomics.RM(0, 0, pars)
 
     return(pars)
-})}
+  })}
 
 #' @title Make inits for RM adult mosquito model
 #' @param pars a [list]
@@ -205,7 +205,7 @@ make_MYZinits_RM = function(pars, MYZopts = list(),
 
     pars$MYZinits = inits
     return(pars)
-})}
+  })}
 
 #' @title Add indices for adult mosquitoes to parameter list
 #' @description Implements [make_indices_MYZ] for the RM model.

R/adult-Ztrace.R

@@ -85,7 +85,7 @@ make_MYZpar_Ztrace = function(pars, MYZopts,
     MYZpar$f0 <- checkIt(f, pars$nPatches)
     MYZpar$q0 <- checkIt(q, pars$nPatches)
     MYZpar$Zm <- checkIt(Zm, pars$nPatches)
-    if(is.null(Zf)) Zf = function(t, y, pars){return(1)}
+    if(is.null(Zf)) Zf = function(t, pars){return(0*t + 1)}
     MYZpar$Zf <- Zf
 
     pars$MYZpar = MYZpar

R/aquatic-trace.R

@@ -15,7 +15,7 @@ LBionomics.trace <- function(t, y, pars) {
 #' @return a [numeric] vector of length `nHabitats`
 #' @export
 F_alpha.trace <- function(t, y, pars) {
-  with(pars$Lpar, Lambda*Lt(t, pars))
+  pars$Lpar$Lt(t, pars)
 }
 
 #' @title Derivatives for aquatic stage mosquitoes
@@ -62,7 +62,7 @@ make_Lpar_trace = function(pars, Lopts=list(),
     Lpar = list()
     class(Lpar) <- "trace"
     Lpar$Lambda = checkIt(Lambda, pars$nHabitats)
-    if(is.null(Lt)) Lt = function(t, pars){1}
+    if(is.null(Lt)) Lt = function(t, pars){pars$Lpar$Lambda}
     Lpar$Lt = Lt
     pars$Lpar <- Lpar
     return(pars)
@@ -94,7 +94,7 @@ update_inits_L.trace<- function(pars, y0) {
 #' @param Lt is a [function] of the form Lt(t,pars) that computes temporal fluctuations
 #' @return none
 #' @export
-make_parameters_L_trace <- function(pars, Lambda, Lt=function(t,pars){1}) {
+make_parameters_L_trace <- function(pars, Lambda, Lt=function(t,pars){pars$Lpar$Lambda}) {
   stopifnot(is.numeric(Lambda))
   Lpar <- list()
   class(Lpar) <- 'trace'

R/diffeqn.R

@@ -55,6 +55,7 @@ xDE_diffeqn <- function(t, y, pars) {
   if(pars$fqZ_out==TRUE) {fqZt = F_fqZ(t, y, pars)}else{fqZt = numeric(0)}
   if(pars$NI_out==TRUE) {NIt = F_X(t, y, pars)/F_H(t, y, pars)}else{NIt = numeric(0)}
   if(pars$kappa_out==TRUE) {kappat = kappa}else{kappat = numeric(0)}
+
   return(list(c(dL, dMYZ, dX, EIRt, fqZt, NIt, kappat)))
 }
 
@@ -77,7 +78,6 @@ xDE_diffeqn_human <- function(t, y, pars) {
   pars <- EIP(t, pars)
 
   # set and modify the baseline mosquito bionomic parameters
-  pars <- EIP(t, pars)
   pars <- MBionomics(t, y, pars)
   pars <- VectorControlEffectSizes(t, y, pars)
 
@@ -95,8 +95,9 @@ xDE_diffeqn_human <- function(t, y, pars) {
 
   if(pars$eir_out==TRUE) {EIRt = EIR}else{EIRt = numeric(0)}
   if(pars$fqZ_out==TRUE) {fqZt = F_fqZ(t, y, pars)}else{fqZt = numeric(0)}
-  if(pars$NI_out==TRUE) {NIt = F_X(t, y, pars)/F_H(t, y, pars)}else{NIt = numeric(0)}
-  if(pars$kappa_out==TRUE) {kappat = kappa}else{kappat = numeric(0)}
+  if(pars$NI_out==TRUE) {NIt = as.numeric(F_X(t, y, pars)/F_H(t, y, pars))}else{NIt = numeric(0)}
+  if(pars$kappa_out==TRUE) {kappat <- F_kappa(t, y, pars, beta)}else{kappat = numeric(0)}
+
   return(list(c(dX, EIRt, fqZt, NIt, kappat)))
 }
 
@@ -153,7 +154,7 @@ xDE_diffeqn_mosy <- function(t, y, pars) {
 xDE_diffeqn_cohort <- function(a, y, pars, F_eir) {
 
   # EIR: entomological inoculation rate trace
-  EIR <- F_eir(a, pars)
+  EIR <- F_eir(a, pars)*pars$Hpar$wts_f
 
   # FoI: force of infection
   FoI <- Exposure(a, y, pars, EIR)
@@ -163,6 +164,7 @@ xDE_diffeqn_cohort <- function(a, y, pars, F_eir) {
 
   if(pars$eir_out==TRUE) {EIRt = EIR}else{EIRt = numeric(0)}
   if(pars$NI_out==TRUE) {NIt = F_X(a, y, pars)/F_H(a, y, pars)}else{NIt = numeric(0)}
+
   return(list(c(dX, EIRt, NIt)))
 }
 

R/human-SIP.R

@@ -78,7 +78,7 @@ dXdt.SIPdXdH <- function(t, y, pars, FoI) {
 setup_X.SIP = function(pars, Xname, Xopts=list()){
 
   pars$Xname = "SIP"
-  pars = make_Xpar_SIP(pars, Xopts)
+    pars = make_Xpar_SIP(pars, Xopts)
   pars = make_Xinits_SIP(pars, Xopts)
 
   return(pars)

R/solving.R

@@ -15,7 +15,7 @@ xde_solve = function(pars, Tmax=365, dt=1){
 #' @param Ymax the number of years to burn-in
 #' @return a [list]
 #' @export
-xde_steady = function(pars, Ymax=10){
+xde_stable_orbit = function(pars, Ymax=10){
   pars <- xde_solve(pars, Tmax = Ymax*365, dt=1)
   deout = tail(pars$orbits$deout, 366)
   deout[,1] = c(0:365)

R/utils.R

@@ -126,31 +126,35 @@ parse_deout <- function(deout, pars){
     if(pars$nStrata>1){
       varslist$eir=apply(cum_eir, 2, diff)/dtime
     } else {
-      varslist$eir= diff(cum_eir)/dtime
+      eir = diff(cum_eir)/dtime
+      varslist$eir= c(eir, eir[1])
     }
   }
   if (pars$fqZ_out == TRUE) {
     cum_fqZ = deout[,1+pars$fqZ_ix]
     if(pars$nPatches>1){
       varslist$fqZ=apply(cum_fqZ, 2, diff)/dtime
     } else {
-      varslist$fqZ= diff(cum_fqZ)/dtime
+      fqZ= diff(cum_fqZ)/dtime
+      varslist$fqZ = c(fqZ, fqZ[1])
     }
   }
   if (pars$NI_out == TRUE) {
-    cum_ni = deout[,1+pars$ni_ix]
+    cum_ni = deout[,1+pars$NI_ix]
     if(pars$nStrata>1){
       varslist$ni=apply(cum_ni, 2, diff)/dtime
     } else {
-      varslist$ni= diff(cum_ni)/dtime
+      ni= diff(cum_ni)/dtime
+      varslist$ni = c(ni, ni[1])
     }
   }
   if (pars$kappa_out == TRUE) {
     cum_kappa = deout[,1+pars$kappa_ix]
     if(pars$nPatches>1){
       varslist$kappa=apply(cum_kappa, 2, diff)/dtime
     } else {
-      varslist$kappa= diff(cum_kappa)/dtime
+      kappa = diff(cum_kappa)/dtime
+      varslist$kappa = c(kappa, kappa[1])
     }
   }
   return(varslist)

_pkgdown.yml

@@ -694,7 +694,7 @@ reference:
   - xde_solve.mosy_dde
   - xde_solve.human
   - xde_solve.cohort
-  - xde_steady
+  - xde_stable_orbit
 - title: Utilities
   contents:
   - make_parameters_xde