From 05b37f683b3a30d55f8ecb1e886e00183eda4362 Mon Sep 17 00:00:00 2001 From: cbernalz Date: Mon, 18 Nov 2024 15:19:10 -0800 Subject: [PATCH] 2024-11-18 update : adding functionality for prior or fixed sigma ww and sigma hosp. --- src/helper_functions.jl | 1 + src/uciwweihr_model.jl | 36 ++++++++++++++++++++++-------- src/uciwweihr_model_params.jl | 41 ++++++++++++++++++++++++----------- 3 files changed, 56 insertions(+), 22 deletions(-) diff --git a/src/helper_functions.jl b/src/helper_functions.jl index 08426bd..fe2f7d6 100644 --- a/src/helper_functions.jl +++ b/src/helper_functions.jl @@ -269,3 +269,4 @@ function calculate_quantiles_without_chain(df, var_prefix, quantiles) return medians, lower_bounds, upper_bounds end + diff --git a/src/uciwweihr_model.jl b/src/uciwweihr_model.jl index ce77d75..67a6de4 100644 --- a/src/uciwweihr_model.jl +++ b/src/uciwweihr_model.jl @@ -44,11 +44,19 @@ The defaults for this fuction will follow those of the default simulation in gen # Parameters for wastewater rho_gene_non_centered ~ Normal() # gene detection rate - sigma_ww_non_centered ~ Normal() # for showing identifyability issue + if isnothing(params.sigma_ww_sd) + sigma_ww = params.sigma_ww + else + sigma_ww_non_centered ~ Normal() + end # Parameters for hospital - sigma_hosp_non_centered ~ Normal() # for showing identifyability issue + if isnothing(params.sigma_hosp_sd) + sigma_hosp = params.sigma_hosp + else + sigma_hosp_non_centered ~ Normal() + end # Non-constant Rt Rt_params_non_centered ~ MvNormal(zeros(l_param_change_times + 2), I) # +2 for sigma and init @@ -74,11 +82,15 @@ The defaults for this fuction will follow those of the default simulation in gen # Parameters for wastewater rho_gene = exp(rho_gene_non_centered * params.rho_gene_sd + params.log_rho_gene_mean) - sigma_ww = exp(sigma_ww_non_centered * params.sigma_ww_sd + params.log_sigma_ww_mean) # for showing identifyability issue + if !isnothing(params.sigma_ww_sd) + sigma_ww = exp(sigma_ww_non_centered * params.sigma_ww_sd + params.log_sigma_ww_mean) + end # Parameters for hospital - sigma_hosp = clamp.(sigma_hosp_non_centered * params.sigma_hosp_sd + params.sigma_hosp_mean, min_neg_bin_sigma, max_neg_bin_sigma) # for showing identifyability issue + if !isnothing(params.sigma_hosp_sd) + sigma_hosp = clamp.(sigma_hosp_non_centered * params.sigma_hosp_sd + params.sigma_hosp_mean, min_neg_bin_sigma, max_neg_bin_sigma) + end # Non-constant Rt Rt_init = exp(Rt_init_non_centered * params.Rt_init_sd + params.Rt_init_mean) @@ -163,8 +175,8 @@ The defaults for this fuction will follow those of the default simulation in gen sigma_Rt = sigma_Rt, rho_gene = rho_gene, - sigma_ww = sigma_ww, # for showing identifyability issue - sigma_hosp = sigma_hosp, # for showing identifyability issue + sigma_ww = sigma_ww, + sigma_hosp = sigma_hosp, H = H_comp_sol, I = I_comp_sol, @@ -209,7 +221,11 @@ The defaults for this fuction will follow those of the default simulation in gen epsilon_non_centered ~ Normal() # Parameters for hospital - sigma_hosp_non_centered ~ Normal() # for showing identifyability issue + if isnothing(params.sigma_hosp_sd) + sigma_hosp = params.sigma_hosp + else + sigma_hosp_non_centered ~ Normal() + end # Non-constant Rt Rt_params_non_centered ~ MvNormal(zeros(l_param_change_times + 2), I) # +2 for sigma and init @@ -234,7 +250,9 @@ The defaults for this fuction will follow those of the default simulation in gen epsilon = exp(epsilon_non_centered * params.epsilon_sd + params.log_epsilon_mean) # Parameters for hospital - sigma_hosp = clamp.(sigma_hosp_non_centered * params.sigma_hosp_sd + params.sigma_hosp_mean, min_neg_bin_sigma, max_neg_bin_sigma) # for showing identifyability issue + if !isnothing(params.sigma_hosp_sd) + sigma_hosp = clamp.(sigma_hosp_non_centered * params.sigma_hosp_sd + params.sigma_hosp_mean, min_neg_bin_sigma, max_neg_bin_sigma) + end # Non-constant Rt Rt_init = exp(Rt_init_non_centered * params.Rt_init_sd + params.Rt_init_mean) @@ -304,7 +322,7 @@ The defaults for this fuction will follow those of the default simulation in gen rt_vals = rt_vals, sigma_Rt = sigma_Rt, - sigma_hosp = sigma_hosp, # for showing identifyability issue + sigma_hosp = sigma_hosp, H = H_comp_sol, I = I_comp_sol, diff --git a/src/uciwweihr_model_params.jl b/src/uciwweihr_model_params.jl index f5a19d8..fcc57e5 100644 --- a/src/uciwweihr_model_params.jl +++ b/src/uciwweihr_model_params.jl @@ -18,8 +18,12 @@ Struct for holding parameters used in the UCIWWEIHR ODE compartmental model. Us - `log_epsilon_mean::Float64=log(1/5)`: Mean for the rate of hospitalization recovery on the log scale. - `rho_gene_sd::Float64=0.02`: Standard deviation for the rho prior. - `log_rho_gene_mean::Float64=log(0.011)`: Mean for the row prior on log scale. -- `sigma_ww::Float64=log(0.1)`: Standard deviation for the normal distribution for wastewater data. Not infered. -- `sigma_hosp::Float64=500.0`: Standard deviation for the negative binomial distribution for hospital data. Not infered. +- `sigma_ww::Union{Float64, Nothing}=log(0.1)`: Standard deviation for the normal distribution for wastewater data. Not infered. +- `sigma_hosp::Union{Float64, Nothing}=500.0`: Standard deviation for the negative binomial distribution for hospital data. Not infered. +- `sigma_ww_sd::Union{Float64, Nothing}=nothing`: Standard deviation for the normal prior of the log standard deviation of the wastewater data. If `nothing`, the sigma_ww is used. +- `log_sigma_ww_mean::Union{Float64, Nothing}=nothing`: Mean for the normal prior of the log standard deviation of the wastewater data. If `nothing`, the sigma_ww is used. +- `sigma_hosp_sd::Union{Float64, Nothing}=nothing`: Standard deviation for the normal prior of the log standard deviation of the hospital data. If `nothing`, the sigma_hosp is used. +- `sigma_hosp_mean::Union{Float64, Nothing}=nothing`: Mean for the normal prior of the log standard deviation of the hospital data. If `nothing`, the sigma_hosp is used. - `Rt_init_sd::Float64=0.3`: Standard deviation for the initial value of the time-varying reproduction number. - `Rt_init_mean::Float64=0.2`: Mean for the initial value of the time-varying reproduction number. - `sigma_Rt_sd::Float64=0.2`: Standard deviation for normal prior of log time-varying reproduction number standard deviation. @@ -44,13 +48,14 @@ struct uciwweihr_model_params log_epsilon_mean::Float64 rho_gene_sd::Float64 log_rho_gene_mean::Float64 - #sigma_ww::Float64 - #sigma_hosp::Float64 - sigma_ww_sd::Float64 - log_sigma_ww_mean::Float64 - sigma_hosp_sd::Float64 - sigma_hosp_mean::Float64 + sigma_ww::Union{Float64, Nothing} + sigma_hosp::Union{Float64, Nothing} + + sigma_ww_sd::Union{Float64, Nothing} + log_sigma_ww_mean::Union{Float64, Nothing} + sigma_hosp_sd::Union{Float64, Nothing} + sigma_hosp_mean::Union{Float64, Nothing} Rt_init_sd::Float64 Rt_init_mean::Float64 @@ -81,11 +86,15 @@ function create_uciwweihr_model_params(; nu_sd::Float64=0.02, log_nu_mean::Float64=log(1/7), epsilon_sd::Float64=0.02, log_epsilon_mean::Float64=log(1/5), rho_gene_sd::Float64=0.02, log_rho_gene_mean::Float64=log(0.011), - #sigma_wastewater::Float64=0.1, - #sigma_hosp::Float64=500.0, + + sigma_wastewater::Union{Float64, Nothing}=0.1, + sigma_hosp::Union{Float64, Nothing}=500.0, - sigma_ww_sd::Float64=0.02, log_sigma_ww_mean::Float64=log(0.1), - sigma_hosp_sd::Float64=50.0, sigma_hosp_mean::Float64=500.0, + #sigma_ww_sd::Float64=0.02, log_sigma_ww_mean::Float64=log(0.1), + #sigma_hosp_sd::Float64=50.0, sigma_hosp_mean::Float64=500.0, + sigma_ww_sd::Union{Float64, Nothing}=nothing, log_sigma_ww_mean::Union{Float64, Nothing}=nothing, + sigma_hosp_sd::Union{Float64, Nothing}=nothing, sigma_hosp_mean::Union{Float64, Nothing}=nothing, + Rt_init_sd::Float64=0.3, Rt_init_mean::Float64=0.2, sigma_Rt_sd::Float64=0.2, sigma_Rt_mean::Float64=-3.0, @@ -93,6 +102,12 @@ function create_uciwweihr_model_params(; sigma_w_sd::Float64=0.2, sigma_w_mean::Float64=-3.5 ) + if sigma_ww_sd == nothing + println("sigma_ww will be fixed") + end + if sigma_hosp_sd == nothing + println("sigma_hosp will be fixed") + end return uciwweihr_model_params( E_init_sd, log_E_init_mean, @@ -102,7 +117,7 @@ function create_uciwweihr_model_params(; nu_sd, log_nu_mean, epsilon_sd, log_epsilon_mean, rho_gene_sd, log_rho_gene_mean, - #sigma_wastewater, sigma_hosp, + sigma_wastewater, sigma_hosp, sigma_ww_sd, log_sigma_ww_mean, sigma_hosp_sd, sigma_hosp_mean,