Skip to content

Commit

Permalink
Fixed group subset error for timevarying demographics in proposal func
Browse files Browse the repository at this point in the history
  • Loading branch information
jameshay218 committed May 13, 2024
1 parent a6fd2da commit 4bd3acd
Show file tree
Hide file tree
Showing 4 changed files with 44 additions and 30 deletions.
1 change: 1 addition & 0 deletions R/helpers.R
Original file line number Diff line number Diff line change
Expand Up @@ -882,6 +882,7 @@ transform_parameters <- function(pars, scale_table, theta_indices,scale_par_indi

add_scale_pars <- function(par_tab, antibody_data, timevarying_demographics=NULL){
if(!is.null(timevarying_demographics)){
timevarying_demographics <- timevarying_demographics %>% arrange(individual, time)
timevarying_demographics <- as.data.frame(timevarying_demographics)
demographics <- create_demographic_table(timevarying_demographics,par_tab)
stratification_pars <- setup_stratification_table(par_tab, demographics)
Expand Down
6 changes: 3 additions & 3 deletions R/posteriors.R
Original file line number Diff line number Diff line change
Expand Up @@ -131,8 +131,10 @@ create_posterior_func <- function(par_tab,
antibody_data[, c("individual", "sample_time", "biomarker_group","biomarker_id")],
antibody_data_unique[, c("individual", "sample_time", "biomarker_group","biomarker_id")]
)

## Setup data vectors and extract
if(!is.null(demographics)){
demographics <- demographics %>% arrange(individual, time)
demographics <- as.data.frame(demographics)
timevarying_demographics <- TRUE
demographic_groups <- create_demographic_table(demographics,par_tab)
Expand All @@ -158,9 +160,7 @@ create_posterior_func <- function(par_tab,
population_group_id_vec <- setup_dat$population_group_id_vec
indiv_group_indices <- setup_dat$indiv_group_indices
n_demographic_groups <- setup_dat$n_demographic_groups
demographics <- setup_dat$demographics


demographics_groups <- setup_dat$demographics

## List of melted antigenic maps, one entry for each observation type
antigenic_map_melted <- setup_dat$antigenic_map_melted
Expand Down
11 changes: 10 additions & 1 deletion src/antibody_models.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -161,11 +161,20 @@ NumericVector antibody_model(const NumericMatrix theta,
// If trying to solve how demographic groups change over time, then need to extract the demographic group indices for this individuals
if(timevarying_groups){
// Treat the indiv_theta_groups vector as a flattened matrix -- this individual's demography vector is from [i-1,2] to [i-1, number_possible_exposures + 1] (the +1 is for birth group)
groups = indiv_theta_groups[Range((number_possible_exposures+1)*(i-1)+1, (number_possible_exposures+1)*i - 1)];
groups = indiv_theta_groups[Range((number_possible_exposures+1)*(i-1) + 1, (number_possible_exposures+1)*i - 1)];
//Rcpp::Rcout << "Indiv: " << i << std::endl;
//Rcpp::Rcout << "Groups: " << groups << std::endl;

// Then subset to groups with infections
groups = groups[use_indices];

//Rcpp::Rcout << "Groups subset: " << groups << std::endl;
// Birth group is [i-1, 1]
birth_group = indiv_theta_groups[(number_possible_exposures+1)*(i-1)];
//Rcpp::Rcout << "Infection times: " << infection_times << std::endl;

//Rcpp::Rcout << "Birth group: " << birth_group << std::endl;
//Rcpp::Rcout << "Theta: " << theta << std::endl;
} else {
// Otherwise, the individual's group is unchanged
group = indiv_theta_groups[i-1];
Expand Down
56 changes: 30 additions & 26 deletions src/proposal.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -416,10 +416,14 @@ List inf_hist_prop_prior_v2_and_v4(
// If trying to solve how demographic groups change over time, then need to extract the demographic group indices for this individuals
if(timevarying_groups){
// Treat the indiv_group_indices vector as a flattened matrix -- this individual's demography vector is from [i-1,2] to [i-1, number_possible_exposures + 1] (the +1 is for birth group)
groups = indiv_group_indices[Range((number_possible_exposures+1)*(indiv)+1, (number_possible_exposures+1)*(indiv+1) - 1)];
groups = indiv_group_indices[Range((number_possible_exposures+1)*(indiv) + 1, (number_possible_exposures+1)*(indiv+1) - 1)];

// Birth group is [i-1, 1]
birth_group = indiv_group_indices[(number_possible_exposures+1)*(indiv)];

//Rcpp::Rcout << "Indiv (prop): " << indiv << std::endl;
//Rcpp::Rcout << "Groups (prop): " << groups << std::endl;
//Rcpp::Rcout << "Birth group (prop): " << birth_group << std::endl;
} else {
// Otherwise, the individual's group is unchanged
group = indiv_group_indices[indiv];
Expand Down Expand Up @@ -614,30 +618,30 @@ List inf_hist_prop_prior_v2_and_v4(
//if(TRUE){
//Rcpp::Rcout << "Infection history after change: " << new_infection_history<< std::endl;

// Calculate likelihood!
indices = new_infection_history > 0;
use_indices =infection_history_mat_indices[indices];
infection_times = possible_exposure_times[use_indices];
infection_times_indices_tmp = possible_exposure_times_indices[use_indices];

// Subset group indices to those times with infections
if(timevarying_groups){
groups_subset = groups[use_indices];
}
//Rcpp::Rcout << infection_times_indices_tmp << std::endl;

// Start end end location of the type_data matrix
type_start = type_data_start(indiv);
type_end = type_data_start(indiv+1)-1;

//Rcpp::Rcout << "Type start: " << type_start << std::endl;
//Rcpp::Rcout << "Type end: " << type_end << std::endl << std::endl;

// ====================================================== //
// =============== CHOOSE MODEL TO SOLVE =============== //
// ====================================================== //
// For each observation type solved for this individual
new_prob = 0;
// Calculate likelihood!
indices = new_infection_history > 0;
use_indices =infection_history_mat_indices[indices];
infection_times = possible_exposure_times[use_indices];
infection_times_indices_tmp = possible_exposure_times_indices[use_indices];
// Subset group indices to those times with infections
if(timevarying_groups){
groups_subset = groups[use_indices];
}
//Rcpp::Rcout << infection_times_indices_tmp << std::endl;
// Start end end location of the type_data matrix
type_start = type_data_start(indiv);
type_end = type_data_start(indiv+1)-1;
//Rcpp::Rcout << "Type start: " << type_start << std::endl;
//Rcpp::Rcout << "Type end: " << type_end << std::endl << std::endl;
// ====================================================== //
// =============== CHOOSE MODEL TO SOLVE =============== //
// ====================================================== //
// For each observation type solved for this individual
new_prob = 0;

for(int index = type_start; index <= type_end; ++index){
//Rcpp::Rcout << "index: " << index << std::endl;
Expand Down Expand Up @@ -686,7 +690,7 @@ List inf_hist_prop_prior_v2_and_v4(
wane_long_parameters,
antigenic_seniority_parameters,
infection_times,
groups,
groups_subset,
birth_group,
infection_times_indices_tmp,
biomarker_id_indices,
Expand Down

0 comments on commit 4bd3acd

Please sign in to comment.