Antibody model fix (#27)

* swap antibody adjustment between infection/vaccination

* perhaps better: let the model decide scaling

* update prior

Author: sbfnk

R/model.R

@@ -12,7 +12,7 @@ i2p_data <- function(prev, ab, vacc, init_ab,
                      prob_detectable, unobserved_time = 14, horizon = 0,
                      inf_ab_delay = c(rep(0, 7 * 4), rep(1 / 7, 7)),
                      vacc_ab_delay = c(rep(0, 7 * 4), rep(1 / 7, 7)),
-                     prop_dont_seroconvert = c(-2, 1), # 10%
+                     prop_seroconvert = c(2, 1), # 90%
                      inf_waning_rate = c(-9, 4),
                      vac_waning_rate = c(-9, 4), # 0.1%
                      vaccine_efficacy = c(3, 1), # 95%
@@ -123,7 +123,7 @@ i2p_data <- function(prev, ab, vacc, init_ab,
     prob_detect_sd = rev(prob_detectable$sd),
     pbt = max(prob_detectable$time) + 1,
     init_inc_mean = logit(init_inc_mean),
-    pbeta = prop_dont_seroconvert,
+    pbeta = prop_seroconvert,
     pgamma_mean = c(inf_waning_rate[1], vac_waning_rate[1]),
     pgamma_sd = c(inf_waning_rate[2], vac_waning_rate[2]),
     pdelta = vaccine_efficacy,

stan/functions/ab.stan

@@ -4,7 +4,7 @@
 //
 // @param vacc Vector of new vaccinations
 // 
-// @param beta Proportion (0 - 1) that don't seroconvert
+// @param beta Proportion (0 - 1) that seroconvert
 //
 // @param gamma Vector of daily rates of antibody waning (0 - 1) for both
 // infection and vaccination. Assumes to exponential waning.
@@ -21,22 +21,23 @@
 //
 // infections <- rep(0.01, 100)
 // vacc <-  rep(0.001, 100)
-// beta <- 0.05
+// beta <- 0.95
 // gamma <- c(0.01, 0.001)
 // delta <- 0.9
 // init_pop_ab <- 0.2
 // t <- 100
 //
 // detectable_antibodies(infections, vacc, beta, gamma, delta, init_pop_ab, t)
 vector detectable_antibodies(vector infections, vector vacc,
-                             real beta, vector gamma, real delta, real init_pop_ab, int t) {
+                             real beta, vector gamma, real delta, real k, real l,
+                             real init_pop_ab, int t) {
   vector[t] pop_ab;
   vector[t] inf_ab;
   vector[t] vac_ab;
 
   // Infection antibodies (assumes all previous from infection)
   inf_ab[1] = init_pop_ab // initial antibodies
-    + (1 - beta) * infections[1] * (1 - init_pop_ab) // new seroconversion
+    + beta * infections[1] // new seroconversion
     - init_pop_ab * gamma[1];
   // Vaccination antibodies
   vac_ab[1] = delta * vacc[1]; // vaccination
@@ -45,11 +46,11 @@ vector detectable_antibodies(vector infections, vector vacc,
   for (i in 2:t) {
     // Infection antibodies
     inf_ab[i] = inf_ab[i - 1] 
-      + (1 - beta) * infections[i] * (1 - pop_ab[i - 1]) // new seroconversion
+      + beta * infections[i] * pow(1 - pop_ab[i - 1], k) // new seroconversion
       - inf_ab[i - 1] * gamma[1]; // waning
     // Vaccination antibodies
     vac_ab[i] = vac_ab[i - 1] 
-      + delta * vacc[i] // new seroconversion
+      + delta * vacc[i]  * pow(1 - pop_ab[i - 1], l) // new seroconversion
       - vac_ab[i - 1] * gamma[2]; // vaccination waning
     // Population antibodies
     pop_ab[i] = inf_ab[i] + vac_ab[i];

stan/inc2prev.stan

@@ -72,6 +72,8 @@ parameters {
   vector<lower = 0, upper = 1>[ab_obs ? 2 : 0] gamma; // antibody waning (inf & vac)
   vector<lower = 0, upper = 1>[ab_obs ? 1 : 0] delta; // vaccine efficacy
   vector<lower = 0, upper = 1>[ab_obs ? 1 : 0] init_dab; // initial proportion with antibodies
+  vector<lower = 0>[n_ab > 0 ? 1 : 0] k; // Potential loss of efficacy from new infections in already seropositive people
+  vector<lower = 0>[n_ab > 0 ? 1 : 0] l; // Potential loss of efficacy from new doses being administered to already seropositive people
 }
 
 transformed parameters {
@@ -110,7 +112,7 @@ transformed parameters {
     infs_with_potential_abs = convolve(infections, inf_ab_delay);
     // calculate detectable antibodies
     dab = detectable_antibodies(infs_with_potential_abs, vacc_with_ab, beta[1],
-                                gamma, delta[1], init_dab[1], t);
+                                gamma, delta[1], k[1], l[1], init_dab[1], t);
     // calculate observed detectable antibodies
     odab = observed_in_window(dab, ab_stime, ab_etime, ut, ab_obs);
     //combined standard error