Skip to content

Commit

Permalink
MCMC amendments
Browse files Browse the repository at this point in the history
MCMC code simplified (always estimate FOI and R0 coefficients, no other options to estimate FOI/R0 directly or estimate FOI coefficients only) and other simplifying edits carried out. Documentation amended with line breaks.
  • Loading branch information
KeithJF82 committed Jun 18, 2024
1 parent aa5a3c5 commit 7c39c6d
Show file tree
Hide file tree
Showing 18 changed files with 714 additions and 1,107 deletions.
2 changes: 1 addition & 1 deletion DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
Package: YEP
Type: Package
Title: Yellow Fever Epidemic Prediction
Version: 0.1.0
Version: 1.0.0
Authors@R: person("Keith", "Fraser", email = "[email protected]", role = c("aut", "cre"))
Maintainers@R: person("Keith", "Fraser", email = "[email protected]", role = c("aut", "cre"))
Description: A package for running a dynamic, deterministic SEIRV model of yellow fever, creating input data and processing output data.
Expand Down
1 change: 0 additions & 1 deletion NAMESPACE
Original file line number Diff line number Diff line change
@@ -1,7 +1,6 @@
# Generated by roxygen2: do not edit by hand

export(Generate_Dataset)
export(Generate_Dataset2)
export(MCMC)
export(Model_Run)
export(Model_Run_Many_Reps)
Expand Down
561 changes: 271 additions & 290 deletions R/main.R

Large diffs are not rendered by default.

648 changes: 305 additions & 343 deletions R/mcmc.R

Large diffs are not rendered by default.

58 changes: 21 additions & 37 deletions R/outputs_generate.R
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
#' @description Generate annual serological and/or case/death data
#'
#' @details This function is used to generate annual serological and/or case/death data based on templates;
#' it is normally used by the single_like_calc() function.
#' it is normally used by the single_posterior_calc() function.
#'
#' [TBA - Explanation of breakdown of regions to model and how to set lengths of FOI_values and R0_values]
#'
Expand All @@ -21,28 +21,26 @@
#' @param p_death_severe_inf Probability of a severe infection resulting in death
#' @param p_rep_severe Probability of reporting of a severe but non-fatal infection
#' @param p_rep_death Probability of reporting of a fatal infection
#' @param mode_start Flag indicating how to set initial population immunity level in addition to vaccination
#' If mode_start=0, only vaccinated individuals
#' If mode_start=1, shift some non-vaccinated individuals into recovered to give herd immunity (uniform by age, R0 based only)
#' If mode_start=2, use SEIRV input in list from previous run(s)
#' @param mode_start Flag indicating how to set initial population immunity level in addition to vaccination \cr
#' If mode_start=0, only vaccinated individuals \cr
#' If mode_start=1, shift some non-vaccinated individuals into recovered to give herd immunity (uniform by age, R0 based only) \cr
#' If mode_start=2, use SEIRV input in list from previous run(s) \cr
#' If mode_start=3, shift some non-vaccinated individuals into recovered to give herd immunity (stratified by age)
#' @param start_SEIRV SEIRV data from end of a previous run to use as input (list of datasets, one per region)
#' @param dt Time increment in days to use in model (should be either 1.0, 2.5 or 5.0 days)
#' @param n_reps number of stochastic repetitions
#' @param deterministic TRUE/FALSE - set model to run in deterministic mode if TRUE
#' @param mode_parallel Set mode for parallelization, if any:
#' If mode_parallel="none", no parallelization
#' If mode_parallel="clusterMap", all regions run in parallel with different time periods and output types
#' @param cluster Cluster of threads to use if mode_parallel="clusterMap"
#' @param output_frame Flag indicating whether to output a complete data frame of results in template format (if TRUE)
#' @param mode_parallel TRUE/FALSE - set model to run in parallel using cluster if TRUE
#' @param cluster Cluster of threads to use if mode_parallel=TRUE
#' @param output_frame TRUE/FALSE - indicate whether to output a complete data frame of results in template format (if TRUE)
#' or calculated values only (if FALSE)
#' '
#' @export
#'
Generate_Dataset <- function(input_data = list(),FOI_values = c(),R0_values = c(),sero_template = NULL,case_template = NULL,
vaccine_efficacy = 1.0, p_severe_inf = 0.12, p_death_severe_inf = 0.39, p_rep_severe = 1.0,
p_rep_death = 1.0,mode_start = 1,start_SEIRV = NULL, dt = 1.0,n_reps = 1, deterministic = FALSE,
mode_parallel = "none",cluster = NULL,output_frame=FALSE){
mode_parallel = FALSE,cluster = NULL,output_frame=FALSE){

assert_that(input_data_check(input_data),msg=paste("Input data must be in standard format",
" (see https://mrc-ide.github.io/YEP/articles/CGuideAInputs.html)"))
Expand All @@ -58,8 +56,8 @@ Generate_Dataset <- function(input_data = list(),FOI_values = c(),R0_values = c(
assert_that(p_rep_severe>=0.0 && p_rep_severe<=1.0,msg="Severe infection reporting probability must be between 0-1")
assert_that(p_rep_death>=0.0 && p_rep_death<=1.0,msg="Fatal infection reporting probability must be between 0-1")
}
assert_that(mode_parallel %in% c("none","clusterMap"))
if(mode_parallel=="clusterMap"){assert_that(is.null(cluster)==FALSE)}
assert_that(is.logical(mode_parallel))
if(mode_parallel){assert_that(is.null(cluster)==FALSE)}

#Prune input data based on regions
regions=regions_breakdown(c(sero_template$region,case_template$region))
Expand Down Expand Up @@ -102,29 +100,15 @@ Generate_Dataset <- function(input_data = list(),FOI_values = c(),R0_values = c(
model_case_values=model_death_values=rep(0,nrow(case_template))
}

#Set up vector of output types to get from model if needed
#if(mode_parallel %in% c("none","clusterMap")){
output_types=rep(NA,n_regions)
for(n_region in 1:n_regions){
if(is.na(case_line_list[[n_region]][1])==FALSE){
if(is.na(sero_line_list[[n_region]][1])==FALSE){output_types[n_region]="case+sero"} else{output_types[n_region]="case"}
} else {output_types[n_region]="sero"}
}
#}

#Model all regions in parallel if parallel modes in use
# if(mode_parallel=="pars_multi"){
# years_data_all=c(min(year_data_begin):max(year_end))
# if(is.null(sero_template)==FALSE){if(is.null(case_template)==FALSE){output_type="case+sero"} else {output_type="sero"}
# } else {output_type="case"}
# model_output_all=Model_Run_Multi_Input(FOI_spillover = FOI_values,R0 = R0_values,
# vacc_data = input_data$vacc_data, pop_data = input_data$pop_data,
# years_data = years_data_all, start_SEIRV=start_SEIRV,output_type = output_type,
# year0 = input_data$years_labels[1],mode_start = mode_start,
# vaccine_efficacy = vaccine_efficacy, dt = dt, n_particles = n_reps,
# n_threads = n_reps*n_regions,deterministic = deterministic)
# }
if(mode_parallel=="clusterMap"){
#Set up vector of output types to get from model
output_types=rep(NA,n_regions)
for(n_region in 1:n_regions){
if(is.na(case_line_list[[n_region]][1])==FALSE){
if(is.na(sero_line_list[[n_region]][1])==FALSE){output_types[n_region]="case+sero"} else{output_types[n_region]="case"}
} else {output_types[n_region]="sero"}
}

if(mode_parallel){
vacc_data_subsets=pop_data_subsets=years_data_sets=list() #TODO - change input data?
for(n_region in 1:n_regions){
vacc_data_subsets[[n_region]]=input_data$vacc_data[n_region,,]
Expand All @@ -144,7 +128,7 @@ Generate_Dataset <- function(input_data = list(),FOI_values = c(),R0_values = c(
for(n_region in 1:n_regions){

#Run model if not using parallelization
if(mode_parallel=="none"){
if(mode_parallel==FALSE){
#cat("\n\t\tBeginning modelling region ",input_data$region_labels[n_region])
model_output = Model_Run(FOI_spillover = FOI_values[n_region],R0 = R0_values[n_region],
vacc_data = input_data$vacc_data[n_region,,],pop_data = input_data$pop_data[n_region,,],
Expand Down
221 changes: 0 additions & 221 deletions R/outputs_generate2.R

This file was deleted.

Loading

0 comments on commit 7c39c6d

Please sign in to comment.