2024-11-18 update : adding functionality for prior or fixed sigma ww …

…and sigma hosp.

src/helper_functions.jl

@@ -269,3 +269,4 @@ function calculate_quantiles_without_chain(df, var_prefix, quantiles)
 
     return medians, lower_bounds, upper_bounds
 end
+

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,

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,18 +86,28 @@ 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,
     w_init_sd::Float64=0.04, w_init_mean::Float64=logit(0.35),
     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,