new sims

IntervalsModel/src/IntervalsModel.jl

@@ -25,7 +25,7 @@ const rng = Xoroshiro128Plus()
 #consts that let give us nicer names for the indices
 const YOUNG, MIDDLE,OLD = 1,2,3
 const durmax = 144
-const sample_repeat = 10
+const sample_repeat = 100
 """
 Number of start times to sample for a given set of distribution parameters. 
 """
@@ -44,9 +44,9 @@ include("plotting_functions.jl")
 Runs parameter estimation for the three scenarios, hh, ws, and rest. 
 """
 function main()
-    # do_hh(400)
-    do_ws(400)
-    do_rest(400)
+    do_hh(500)
+    do_ws(500)
+    do_rest(500)
     plot_all()
 end
 
@@ -88,7 +88,7 @@ Number of particles to obtain representing posterior distribution of the distrib
 
 Threshold adaptive rate, see `KissABC.smc` for more details.
 """
-function bayesian_estimation(fname, err_func, priors_list, dists, particles; alpha = 0.99)
+function bayesian_estimation(fname, err_func, priors_list, dists, particles; alpha = 0.98)
     
     data_pairs = map(zip(dists,priors_list)) do (dist,priors)
         init = rand(priors)
@@ -112,11 +112,12 @@ function do_hh(particles)
     ]
 
     poisson_priors = filter_priors("home")
+    display(poisson_priors)
     # Set parameter priors for fitting
     priors_list = [
         Factored(poisson_priors...),
     ]
-    bayesian_estimation("hh",err_hh,priors_list,dists,particles; alpha = 0.99)
+    bayesian_estimation("hh",err_hh,priors_list,dists,particles; alpha = 0.98)
 end
 
 
@@ -136,7 +137,7 @@ function do_ws(particles)
     priors_list = [
         Factored(poisson_priors...)
     ]
-    bayesian_estimation("ws",err_ws,priors_list,dists,particles; alpha = 0.99)
+    bayesian_estimation("ws",err_ws,priors_list,dists,particles; alpha = 0.98)
 end
 
 """

IntervalsModel/src/hh_durations_model.jl

@@ -38,7 +38,7 @@ function err_hh(p,dist)
     # display(p_matrix)
     AGERESP =  dat.AGERESP #age of the respondent
     errsum = 0
-    @inbounds for i = 1:length(duration_subarray) #loop through entire file
+    @inbounds for i = 1:length(AGERESP) #loop through entire file
         age_sample = AGERESP[i]
         @inbounds for age_j in YOUNG:OLD #for a given age_sample loop over possible contact ages
 

IntervalsModel/src/plotting_functions.jl

 using LaTeXStrings
 using Plots
-const color_palette = palette(:seaborn_bright) #color theme for the plots
+const color_palette = palette(:seaborn_pastel) #color theme for the plots
+const color_palette_bright = palette(:seaborn_bright) #bright color theme for the plots
 
 default(dpi = 300)
 default(framestyle = :box)
@@ -34,7 +35,12 @@ function plot_dists(fname,dist_constructors,data)
             hasnans = any(any.(map(l -> isnan.(l),dist_pts)))
             err_down = hasnans ? 0 : quantile.(dist_pts,0.05)
             err_up = hasnans ? 0 : quantile.(dist_pts,0.95)
-            plot!(p_matrix[i,j] ,x_range,mean_dat; ribbon = ( mean_dat .- err_down,err_up .- mean_dat),legend = false,label = string(dist_constructor),seriescolor = color_palette[k])
+            plot!(p_matrix[i,j] ,x_range,mean_dat;
+                ribbon = ( mean_dat .- err_down,err_up .- mean_dat),
+                legend = false,
+                label = string(dist_constructor),
+                seriescolor = color_palette_bright[k]
+            )
         end
         annotate!(p_matrix[i,j],compute_x_pos(p_matrix[i,j]),compute_y_pos(p_matrix[i,j]), Plots.text("$(ymo[i])→$(ymo[j])", :left, 10))
     end
@@ -54,7 +60,7 @@ function plot_posteriors(fname,data)
         hist = fit(Histogram,data.P[i].particles; nbins = 40)
         kde_est = kde(data.P[i].particles)
         kernel_data = [pdf(kde_est,x) for x in minimum(data.P[i].particles):maximum(data.P[i].particles)]
-        plot!(p_list[i],hist;legend = false, xlabel = L"\lambda_%$i")
+        plot!(p_list[i],hist;legend = false, xlabel = L"\lambda_%$i", color = color_palette[1]  )
         # display(kernel_data)
         # vline!(p_list[i],[argmax(kernel_data)]; seriescolor = color_palette[2])