Fixed some bugs

seroanalytics · May 17, 2024 · 05c9b6a · 05c9b6a
1 parent ea9a3db
commit 05c9b6a
Show file tree

Hide file tree

Showing 5 changed files with 15 additions and 13 deletions.
diff --git a/R/helpers.R b/R/helpers.R
@@ -77,7 +77,7 @@ get_n_alive_group <- function(antibody_data, times, demographics=NULL, melt_data
     sample_mask <- times[create_sample_mask(antibody_data, times)]
     masks <- data.frame(cbind(age_mask, sample_mask))
     DOBs <- cbind(DOBs, masks)
-    n_alive <- plyr::ddply(DOBs, ~population_group, function(y) sapply(seq(1, length(times)), function(x)
+    n_alive <- plyr::ddply(DOBs, ~population_group, function(y) sapply(times, function(x)
         nrow(y[y$age_mask <= x & y$sample_mask >= x, ])))
   }
     n_alive <- as.matrix(n_alive[, 2:ncol(n_alive)])

diff --git a/R/plot_antibody_model.R b/R/plot_antibody_model.R
@@ -207,7 +207,6 @@ plot_model_fits <- function(chain, infection_histories,
     measurement_ranges <- par_tab %>% dplyr::filter(names %in% c("min_measurement","max_measurement")) %>% dplyr::select(names,values,biomarker_group) %>%
       pivot_wider(names_from=names,values_from=values)
   }
-
   max_x <- max(inf_hist_densities$variable) + 5
   time_range <- range(inf_hist_densities$variable)
   ## If provided, add true infection histories
@@ -217,7 +216,7 @@ plot_model_fits <- function(chain, infection_histories,
     known_infection_history <- reshape2::melt(known_infection_history)
     colnames(known_infection_history) <- c("individual","variable","inf")
     known_infection_history <- known_infection_history[known_infection_history$inf == 1,]
-    known_infection_history$variable <- possible_exposure_times[known_infection_history$variable]
+    known_infection_history$variable <- possible_exposure_times[as.numeric(as.factor(known_infection_history$variable))]
     known_infection_history$individual <- individuals[known_infection_history$individual]
     expand_samples <- expand_grid(individual=individuals,sample_time=unique(to_use$sample_time))
     known_infection_history <- known_infection_history %>% left_join(expand_samples,by="individual",relationship="many-to-many") %>% filter(variable <= sample_time)

diff --git a/R/posteriors.R b/R/posteriors.R
@@ -43,7 +43,7 @@ create_posterior_func <- function(par_tab,
                                   antibody_data,
                                   antigenic_map=NULL,
                                   possible_exposure_times=NULL,
-                                  prior_version = 1,
+                                  prior_version = 2,
                                   solve_likelihood = TRUE,
                                   age_mask = NULL,
                                   measurement_bias = NULL,

diff --git a/R/simulate_data.R b/R/simulate_data.R
@@ -52,7 +52,7 @@ simulate_data <- function(par_tab,
                           age_group_bounds = NULL,
                           attack_rates,
                           repeats = 1,
-                          measurement_indices = NULL,
+                          measurement_bias = NULL,
                           data_type = NULL,
                           demographics=NULL,
                           verbose=FALSE) {
@@ -100,6 +100,7 @@ simulate_data <- function(par_tab,
         timevarying_demographics <- demographics %>% expand_grid(time=possible_exposure_times)
         timevarying_demographics$age <- timevarying_demographics$time - timevarying_demographics$birth
         timevarying_demographics$age_group <- as.numeric(cut(timevarying_demographics$age, breaks=c(0, age_group_bounds)))
+        timevarying_demographics$age_group <- as.numeric(as.factor(timevarying_demographics$age_group))
       } else {
         timevarying_demographics <- NULL
       }
@@ -151,10 +152,8 @@ simulate_data <- function(par_tab,
     #########################################################
     message("Simulating data\n")    
 
-    measurement_bias <- NULL
-    if (!is.null(measurement_indices)) {
+    if (!is.null(measurement_bias)) {
         message(cat("Measurement bias\n"))
-        #measurement_bias <- pars[measurement_indices_par_tab]
     }
 
     ## Simulate infection histories
@@ -166,7 +165,11 @@ simulate_data <- function(par_tab,
     }
     ## Merge with final sample time and any relevant population group keys
     use_demo <- use_demo %>% left_join(final_sample_times,by="individual")
-    use_demo <- suppressMessages(use_demo %>% left_join(population_groups))
+    if(is.null(population_groups)){
+      use_demo$population_group <- 1
+    } else {
+      use_demo <- suppressMessages(use_demo %>% left_join(population_groups))
+    }
 
     ## Get simulated infection histories and attack rates
     tmp <- simulate_infection_histories(
@@ -191,11 +194,11 @@ simulate_data <- function(par_tab,
     ## Correct arrangement
     antibody_data <- antibody_data %>% 
       arrange(individual, biomarker_group, sample_time, biomarker_id, repeat_number)
-
     ## Simulate data!
     f <- create_posterior_func(par_tab,antibody_data,antigenic_map,function_type=3,
                                possible_exposure_times = possible_exposure_times,
                                demographics=timevarying_demographics,
+                               measurement_bias=measurement_bias,
                                start_level="none")
     antibody_data$measurement <- f(par_tab$values, infection_history)
     ## Add noise, but need to be specific to the data type

diff --git a/README.Rmd b/README.Rmd
@@ -36,12 +36,12 @@ devtools::install_github("seroanalytics/serosolver",ref="published")
 * Moving more options and inputs behind the scenes to streamline the user interface
 * Generalization to consider multiple biomarker types per sample (e.g., antibody titre and avidity)
 * Support for continuous as well as discrete observations (e.g., can now fit to ELISA data as well as HAI titres)
+* Model infection histories and antibody kinetics as a function of demographic variables
+* Allow some infection states to be fixed during fitting
+* Ways to fix or estimate starting/baseline titres
 * _IN PROGRESS_ Some small improvements to the MCMC sampler and parameter transformations
 * _IN PROGRESS_ Improved guidance and support for using priors
-* _IN PROGRESS_ Model infection histories and antibody kinetics as a function of demographic variables
-* _IN PROGRESS_ Allow some infection states to be fixed during fitting
 * _IN PROGRESS_ Inclusion of explicit immunity model
-* _IN PROGRESS_ Ways to fix or estimate starting/baseline titres
 * _IN PROGRESS_ Added tests
 </details>