degree dependent initial conditions

CovidAlertVaccinationModel/abm_sensitivity_analysis.jl

@@ -8,14 +8,6 @@ const samples = 10
 const len = 10 #number of points to evaluate
 gr()
 
-#add beta
-#run model without vaccination
-#use derivative of log of no infections to calibrate p
-#seed inital infections according to degree
-#total number of notifications
-#age structured vaccination
-#
-
 const univarate_test_list = (
     # (:I_0_fraction, range(0.0, 0.05; length = len)), 
     # (:base_transmission_probability, range(0.0002, 0.002; length = len)),

CovidAlertVaccinationModel/src/ABM/mixing_graphs.jl

@@ -10,7 +10,7 @@ struct GraphEdge
 end
 
 """
-Define a hash on GraphEdge such that ``hash(a,b) = hash(b,a)`` (hash is commutative).
+Define a hash on GraphEdge such that `hash(a,b) = hash(b,a)` (hash is commutative).
 
 This is helpful because then we only need to store (a,b) in the graph edges weights dictionary, rather than both (a,b) and (b,a).
 """
@@ -238,9 +238,7 @@ function fast_chung_lu(pop_i,mixing_dist)
     return GraphEdge.(stubs_i,stubs_j)
 end
 
-
-
-neighbors(g::WeightedGraph,i) = neighbors(g.g,i)
+neighbors(g::WeightedGraph,i) = LightGraphs.neighbors(g.g,i)
 get_weight(g::WeightedGraph,e) = g.weights_dict[e]
 function Base.show(io::IO, g::WeightedGraph) 
     print(io, "WG $(ne(g.g))")

CovidAlertVaccinationModel/src/ABM/model_setup.jl

@@ -3,7 +3,7 @@ function get_parameters()#(0.0000,0.00048,0.0005,0.16,-1.30,-1.24,-0.8,0.35,0.35
         sim_length = 500,
         num_households = 5000,
         I_0_fraction = 0.003,
-        β_y = 0.0011,
+        β_y = 0.001095,
         β_m = 0.00061,
         β_o = 0.04,
         α_y = 0.4,
@@ -15,7 +15,7 @@ function get_parameters()#(0.0000,0.00048,0.0005,0.16,-1.30,-1.24,-0.8,0.35,0.35
         π_base_o = -0.95,
         η = 0.0,
         κ = 0.0,
-        ω = 0.0056,
+        ω = 0.0055,
         ω_en = 0.00,
         Γ = 1/7,
         ξ = 5.0,

CovidAlertVaccinationModel/src/ABM/solve.jl

@@ -142,12 +142,25 @@ function weighted_degree(node,network::TimeDepMixingGraph)
     return weighted_degree
 end
 
+function sample_initial_nodes(nodes,graphs,I_0_fraction)
+    weighted_degrees = zeros(nodes)
+    for v in 1:nodes
+        for g in graphs
+            for w in neighbors(g,v)
+                weighted_degrees[v] += get_weight(g,GraphEdge(v,w))
+            end
+        end
+    end
+    wv = Weights(weighted_degrees ./sum(weighted_degrees))
+    num = round(Int,nodes*I_0_fraction)
+    init_indices = sample(Random.default_rng(Threads.threadid()), 1:nodes,wv, num; replace = false)
+    return init_indices
+end    
 
-function solve!(modelsol,recordings...)
-    init_indices = rand(Random.default_rng(Threads.threadid()), 1:modelsol.nodes, round(Int,modelsol.nodes*modelsol.params.I_0_fraction))
 
-    for t in 1:modelsol.sim_length
-        
+function solve!(modelsol,recordings...)
+    init_indices = sample_initial_nodes(modelsol.nodes, modelsol.inf_network.graph_list[begin], modelsol.params.I_0_fraction)
+    for t in 1:modelsol.sim_length  
         #this also resamples the soc network weights since they point to the same objects, but those are never used
         if t>1
             remake!(t,modelsol.inf_network,modelsol.index_vectors,modelsol.demographics)

CovidAlertVaccinationModel/src/CovidAlertVaccinationModel.jl

@@ -53,7 +53,8 @@ include("ABM/model_setup.jl")
 include("ABM/parameter_optimization.jl")
 include("ABM/output.jl")
 include("ABM/solve.jl")
-include("ABM/abm.jl")
+include("ABM/abm.jl") 
+include("ABM/parameter_planes.jl")
 
 include("IntervalsModel/intervals_model.jl")
 include("IntervalsModel/interval_overlap_sampling.jl")

CovidAlertVaccinationModel/test/ABM/output_test.jl

@@ -25,15 +25,15 @@ using Random
     println("obs postinf vac: $total_postinf_vaccination,target: $target_postinf_vac")
     @test all(abs.(final_size .- target_final_size) .< [75,100,25])
     @test all(abs.(total_preinf_vaccination .- target_preinf_vac) .< [50,100,25])
-    @test all(abs.(total_postinf_vaccination .- target_postinf_vac) .< [120,100,25])
+    @test all(abs.(total_postinf_vaccination .- target_postinf_vac) .< [120,120,40])
 end
 
 @testset "perf" begin
     b = @belapsed CovidAlertVaccinationModel.bench()
-    @test  4.2 < b < 4.85
-    if Threads.nthreads() == 20
+    @test  4.0 < b < 4.5
+    if Threads.nthreads() == 25
         b = @belapsed CovidAlertVaccinationModel.threaded_bench()
-        @test  10.65 < b < 15.00
+        @test  12.0 < b < 16.00
     else
         display("incorrect number of threads, skipping threaded perf test")
     end