diff --git a/src/uciwweihr_model.jl b/src/uciwweihr_model.jl
index 92b29fb..79c248f 100644
--- a/src/uciwweihr_model.jl
+++ b/src/uciwweihr_model.jl
@@ -43,7 +43,13 @@ The defaults for this fuction will follow those of the default simulation in gen
         epsilon_non_centered ~ Normal()
         # Parameters for wastewater
         rho_gene_non_centered ~ Normal() # gene detection rate
+
+        tau_non_centered ~ Normal() # for showing identifyability issue
+
         # Parameters for hospital
+
+        sigma_hosp_non_centered ~ Normal() # for showing identifyability issue
+
         # Non-constant Rt
         Rt_params_non_centered ~ MvNormal(zeros(l_param_change_times + 2), I) # +2 for sigma and init
         sigma_Rt_non_centered = Rt_params_non_centered[1]
@@ -67,7 +73,13 @@ 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 wastewater
         rho_gene = exp(rho_gene_non_centered * params.rho_gene_sd + params.log_rho_gene_mean)
+
+        tau = exp(tau_non_centered * params.tau_sd + params.log_tau_mean) # for showing identifyability issue
+
         # 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
+
         # Non-constant Rt
         Rt_init = exp(Rt_init_non_centered * params.Rt_init_sd + params.Rt_init_mean)
         sigma_Rt = exp(sigma_Rt_non_centered * params.sigma_Rt_sd + params.sigma_Rt_mean)
@@ -120,10 +132,12 @@ The defaults for this fuction will follow those of the default simulation in gen
         # E - 1 // I - 2 // H - 3 // R - 4
         # Likelihood calculations------------
         for i in 1:l_obs_ww
-            data_wastewater[i] ~ Normal(log_W_means[i], params.sigma_wastewater)
+            #data_wastewater[i] ~ Normal(log_W_means[i], params.sigma_wastewater)
+            data_wastewater[i] ~ Normal(log_W_means[i], tau)
         end
         for i in 1:l_obs_hosp
-            data_hosp[i] ~ NegativeBinomial2(H_means[i], params.sigma_hosp)
+            #data_hosp[i] ~ NegativeBinomial2(H_means[i], params.sigma_hosp)
+            data_hosp[i] ~ NegativeBinomial2(H_means[i], sigma_hosp)
         end
 
 
@@ -147,6 +161,10 @@ The defaults for this fuction will follow those of the default simulation in gen
             rt_vals = rt_vals,
             sigma_Rt = sigma_Rt,
             rho_gene = rho_gene,
+
+            tau = tau, # for showing identifyability issue
+            sigma_hosp = sigma_hosp, # for showing identifyability issue
+            
             H = H_means,
             log_genes_mean = log_W_means,
             rt_init = rt_init,
diff --git a/src/uciwweihr_model_params.jl b/src/uciwweihr_model_params.jl
index b9db889..646d3f7 100644
--- a/src/uciwweihr_model_params.jl
+++ b/src/uciwweihr_model_params.jl
@@ -44,8 +44,14 @@ struct uciwweihr_model_params
     log_epsilon_mean::Float64
     rho_gene_sd::Float64
     log_rho_gene_mean::Float64
-    sigma_wastewater::Float64
-    sigma_hosp::Float64
+    #sigma_wastewater::Float64
+    #sigma_hosp::Float64
+
+    tau_sd::Float64
+    log_tau_mean::Float64
+    sigma_hosp_sd::Float64
+    sigma_hosp_mean::Float64
+
     Rt_init_sd::Float64
     Rt_init_mean::Float64
     sigma_Rt_sd::Float64
@@ -75,8 +81,12 @@ 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::Float64=0.1,
+    #sigma_hosp::Float64=500.0,
+    
+    tau_sd::Float64=0.02, log_tau_mean::Float64=log(0.1),
+    sigma_hosp_sd::Float64=50.0, sigma_hosp_mean::Float64=500.0,
+
     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),
@@ -92,7 +102,11 @@ 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,
+
+        tau_sd, log_tau_mean,
+        sigma_hosp_sd, sigma_hosp_mean,
+
         Rt_init_sd, Rt_init_mean,
         sigma_Rt_sd, sigma_Rt_mean,
         w_init_sd, w_init_mean,
diff --git a/src/uciwweihr_visualizer.jl b/src/uciwweihr_visualizer.jl
index e468292..e5f115f 100644
--- a/src/uciwweihr_visualizer.jl
+++ b/src/uciwweihr_visualizer.jl
@@ -54,6 +54,7 @@ function uciwweihr_visualizer(
         ["gamma", "nu", "epsilon"],
         ["rt_init", "w_init"],
         ["sigma_w", "sigma_Rt"],
+        ["tau", "sigma_hosp"], # For identifyability issue checking
         ["rho_gene"]
     ],
     time_varying_params = ["rt_vals", "w_t"],
@@ -169,6 +170,7 @@ function uciwweihr_visualizer(
         ["gamma", "nu", "epsilon"],
         ["rt_init", "w_init"],
         ["sigma_w", "sigma_Rt"],
+        ["tau", "sigma_hosp"], # For identifyability issue checking
         ["rho_gene"]
     ],
     time_varying_params = ["rt_vals", "w_t"],