From cd6462b39ad6da1e4c6fd302188b041599929da2 Mon Sep 17 00:00:00 2001
From: pjentsch <pjentsch@uwaterloo.ca>
Date: Sat, 20 Feb 2021 01:29:49 -0500
Subject: [PATCH] removed LTC/HCW stuff, contact generation needs fixing, use
bytearrays instead of bitarrays for weights
---
.../src/CovidAlertVaccinationModel.jl | 8 +++--
CovidAlertVaccinationModel/src/agents.jl | 16 +++------
CovidAlertVaccinationModel/src/graphs.jl | 35 ++++---------------
.../src/mixing_distributions.jl | 20 +++++------
CovidAlertVaccinationModel/test/runtests.jl | 4 +--
5 files changed, 28 insertions(+), 55 deletions(-)
diff --git a/CovidAlertVaccinationModel/src/CovidAlertVaccinationModel.jl b/CovidAlertVaccinationModel/src/CovidAlertVaccinationModel.jl
index c9b198f..ca0f71b 100644
--- a/CovidAlertVaccinationModel/src/CovidAlertVaccinationModel.jl
+++ b/CovidAlertVaccinationModel/src/CovidAlertVaccinationModel.jl
@@ -41,11 +41,15 @@ default(dpi = 300)
default(framestyle = :box)
using BenchmarkTools
function main()
- agent_model = AgentModel(5000,5)
# display(size.(agent_model.demographic_index_vectors))
+
+ agent_model = AgentModel(5000)
u_0 = get_u_0(length(agent_model.demographics))
steps = 300
- # Random.seed!(RNG,1)
+ Random.seed!(RNG,1)
+
+ @btime solve!($u_0,$get_parameters(),$steps,$agent_model,$vaccinate_uniformly!);
+ # agent_model = AgentModel(10_000)
@btime solve!($u_0,$get_parameters(),$steps,$agent_model,$vaccinate_uniformly!);
# sol,graphs = solve!(u_0,get_parameters(),steps,agent_model,vaccinate_uniformly!);
# return aggregate_timeseries(sol)
diff --git a/CovidAlertVaccinationModel/src/agents.jl b/CovidAlertVaccinationModel/src/agents.jl
index 1632925..09dc805 100644
--- a/CovidAlertVaccinationModel/src/agents.jl
+++ b/CovidAlertVaccinationModel/src/agents.jl
@@ -1,11 +1,8 @@
-#Agent type definitions
-
+#agent type definitions
@enum AgentDemographic begin
Young = 1
Adult
Elderly
- AdultHCW
- ElderlyLTC
end
@enum AgentStatus begin
Susceptible = 1
@@ -22,8 +19,8 @@ struct AgentModel{GraphType,DistMatrixList1,DistMatrixList2,DistMatrixList3}
workschool_contacts_mean_adjusted_list::DistMatrixList1
rest_contacts_mean_adjusted_list::DistMatrixList2
home_contacts_duration_list::DistMatrixList3
- function AgentModel(num_households,num_ltc)
- pop_list,base_network,index_vectors = generate_population(num_households,num_ltc)
+ function AgentModel(num_households)
+ pop_list,base_network,index_vectors = generate_population(num_households)
pop_sizes = length.(index_vectors)
workschool_contacts_mean_adjusted_list = [symmetrize_means(pop_sizes,initial_workschool_mixing_matrix)]
@@ -65,6 +62,8 @@ function Base.show(io::IO, status::AgentStatus)
print(io, "I")
elseif status == Recovered
print(io, "R")
+ elseif status == Immunized
+ print(io, "Vac")
end
end
@@ -73,12 +72,7 @@ function Base.show(io::IO, status::AgentDemographic)
print(io, "<20")
elseif status == Adult
print(io, "20-60")
- elseif status == AdultHCW
- print(io, "20-60 HCW")
elseif status == Elderly
print(io, ">60")
- elseif status == ElderlyLTC
- print(io, ">60 LTC")
end
end
-
diff --git a/CovidAlertVaccinationModel/src/graphs.jl b/CovidAlertVaccinationModel/src/graphs.jl
index df195d6..e7b95cd 100644
--- a/CovidAlertVaccinationModel/src/graphs.jl
+++ b/CovidAlertVaccinationModel/src/graphs.jl
@@ -86,7 +86,7 @@ function generate_contact_vectors!(ij_dist,ji_dist,i_to_j_contacts, j_to_i_conta
end
#add a bipartite graph derived from mixing matrices onto g
function generate_mixing_graph!(g,index_vectors,contacts)
- for i in 1:5, j in 1:5
+ for i in 1:length(index_vectors), j in 1:length(index_vectors)
random_bipartite_graph_fast_CL!(g,index_vectors[i],index_vectors[j],contacts.contact_array[i,j],contacts.contact_array[j,i])
end
return g
@@ -109,7 +109,7 @@ function random_bipartite_graph_fast_CL!(g::SimpleGraph,anodes,bnodes,aseq,bseq)
return g
end
-#same algorithm but with a bitarray
+#same algorithm but with a bitarray, this cannot produce parallel edges
function random_bipartite_graph_fast_CL!(g::T,anodes,bnodes,aseq,bseq) where T<:AbstractArray
lena = length(aseq)
lenb = length(bseq)
@@ -124,41 +124,18 @@ function random_bipartite_graph_fast_CL!(g::T,anodes,bnodes,aseq,bseq) where T<:
end
#generate population with households distributed according to household_size_distribution
-#Assume 5% of elderly are in LTC roughly, these are divide evenly among a specific number of LTC ltc_homes
-#LTC homes and households are assumed to be complete static graphs.
-function generate_population(num_households,num_LTC)
- # LTC_sizes = rand(LTC_distribution,num_LTC)
- prob_adult_is_HCW = 0.01#assume 1% of adults are HCW
- prob_elderly_in_LTC = 0.05 #5% of elderly are in LTC
-
+function generate_population(num_households)
households_composition = sample_household_data(num_households)
households = map( l -> [fill(AgentDemographic(i),n) for (i,n) in enumerate(l)],households_composition) |>
l -> reduce(vcat,l) |>
l -> reduce(vcat,l)
total_household_pop = sum(sum.(households_composition))
- ltc_pop = Int(div(total_household_pop,(1 - prob_elderly_in_LTC))) - total_household_pop
- hcw_pop = Int(ceil(prob_adult_is_HCW*total_household_pop))
- ltc_pop_list = fill(ElderlyLTC,ltc_pop)
-
- adults_indices = findall(x -> x == Adult, households)
- households[sample(RNG,adults_indices,hcw_pop)] .= AdultHCW
-
- ltc_home_size = Int(ceil(ltc_pop/num_LTC))
- @assert ltc_home_size != 0 error("must have at least one LTC resident,make problem larger")
- if ltc_home_size != 0
- ltc_homes = [ltc_pop_list[i:min(i + ltc_home_size - 1,ltc_pop)] for i = 1:ltc_home_size:ltc_pop]
- else
- ltc_homes = Vector{Vector{AgentDemographic}}()
- end
- ltc_networks = map(LightGraphs.complete_graph,length.(ltc_homes))
household_networks = map(LightGraphs.complete_graph, sum.(households_composition))
- population_list = vcat(reduce(vcat,households), reduce(vcat,ltc_homes))
- network = blockdiag(reduce(blockdiag,household_networks; init = SimpleGraph()), reduce(blockdiag,ltc_networks; init = SimpleGraph()))
-
- index_vectors = [findall(x -> x == AgentDemographic(i), population_list) for i in 1:AgentDemographic.size+1]
-
+ population_list = reduce(vcat,households)
+ network = reduce(blockdiag,household_networks; init = SimpleGraph())
+ index_vectors = [findall(x -> x == AgentDemographic(i), population_list) for i in 1:(AgentDemographic.size-1)]
return (;
population_list,
network,
diff --git a/CovidAlertVaccinationModel/src/mixing_distributions.jl b/CovidAlertVaccinationModel/src/mixing_distributions.jl
index f3f8df6..8da1851 100644
--- a/CovidAlertVaccinationModel/src/mixing_distributions.jl
+++ b/CovidAlertVaccinationModel/src/mixing_distributions.jl
@@ -56,26 +56,24 @@ function from_mean(::Type{ZWDist{DistType,T}},α,μ) where {DistType <: Distribu
end
StatsBase.mean(d::ZWDist{Dist,T}) where {Dist,T} = (1 - d.α)*StatsBase.mean(d.base_dist)
+
const initial_workschool_type = Union{ZWDist{Geometric{Float64},Discrete},ZWDist{Poisson{Float64},Discrete}}
const initial_workschool_mixing_matrix = convert(Array{initial_workschool_type},map(t->from_mean(t...),[
- (ZWDist{Geometric{Float64},Discrete}, 0.433835,4.104848) (ZWDist{Geometric{Float64},Discrete},0.406326,2.568782) (ZWDist{Poisson{Float64},Discrete},0.888015,0.017729) (ZWDist{Geometric{Float64},Discrete},0.406326,2.568782) (ZWDist{Poisson{Float64},Discrete},0.888015,0.017729)
- (ZWDist{Geometric{Float64},Discrete}, 0.600966,0.975688) (ZWDist{Geometric{Float64},Discrete},0.4306,5.057572) (ZWDist{Poisson{Float64},Discrete}, 0.84513,0.021307) (ZWDist{Geometric{Float64},Discrete},0.4306,5.057572) (ZWDist{Poisson{Float64},Discrete}, 0.84513,0.021307)
- (ZWDist{Poisson{Float64},Discrete},0.887392,0.001937) (ZWDist{Poisson{Float64},Discrete},0.793004,0.00722) (ZWDist{Poisson{Float64},Discrete},0.940473, 0.022134) (ZWDist{Poisson{Float64},Discrete},0.793004,0.00722) (ZWDist{Poisson{Float64},Discrete},0.940473, 0.022134)
- (ZWDist{Geometric{Float64},Discrete}, 0.600966,0.975688) (ZWDist{Geometric{Float64},Discrete},0.4306,5.057572) (ZWDist{Poisson{Float64},Discrete}, 0.84513,0.021307) (ZWDist{Geometric{Float64},Discrete},0.4306,5.057572) (ZWDist{Poisson{Float64},Discrete}, 0.84513,0.021307)
- (ZWDist{Poisson{Float64},Discrete},0.887392,0.001937) (ZWDist{Poisson{Float64},Discrete},0.793004,0.00722) (ZWDist{Poisson{Float64},Discrete},0.940473, 0.022134) (ZWDist{Poisson{Float64},Discrete},0.793004,0.00722) (ZWDist{Poisson{Float64},Discrete},0.940473, 0.022134)
+ (ZWDist{Geometric{Float64},Discrete}, 0.433835,4.104848) (ZWDist{Geometric{Float64},Discrete},0.406326,2.568782) (ZWDist{Poisson{Float64},Discrete},0.888015,0.017729)
+ (ZWDist{Geometric{Float64},Discrete}, 0.600966,0.975688) (ZWDist{Geometric{Float64},Discrete},0.4306,5.057572) (ZWDist{Poisson{Float64},Discrete}, 0.84513,0.021307)
+ (ZWDist{Poisson{Float64},Discrete},0.887392,0.001937) (ZWDist{Poisson{Float64},Discrete},0.793004,0.00722) (ZWDist{Poisson{Float64},Discrete},0.940473, 0.022134)
]))
const initial_rest_type = Union{Geometric{Float64},Poisson{Float64}}
const initial_rest_mixing_matrix = convert(Array{initial_rest_type},map(t->from_mean(t...),[
- (Geometric{Float64},2.728177) (Geometric{Float64},1.382557) (Poisson{Float64},0.206362) (Geometric{Float64},1.382557) (Poisson{Float64},0.206362)
- (Poisson{Float64},1.139072) (Geometric{Float64},3.245594) (Poisson{Float64},0.785297) (Geometric{Float64},3.245594) (Poisson{Float64},0.785297)
- (Poisson{Float64},0.264822) (Poisson{Float64},0.734856) (Poisson{Float64},0.667099) (Poisson{Float64},0.734856) (Poisson{Float64},0.667099)
- (Poisson{Float64},1.139072) (Geometric{Float64},3.245594) (Poisson{Float64},0.785297) (Geometric{Float64},3.245594) (Poisson{Float64},0.785297)
- (Poisson{Float64},0.264822) (Poisson{Float64},0.734856) (Poisson{Float64},0.667099) (Poisson{Float64},0.734856) (Poisson{Float64},0.667099)
+ (Geometric{Float64},2.728177) (Geometric{Float64},1.382557) (Poisson{Float64},0.206362)
+ (Poisson{Float64},1.139072) (Geometric{Float64},3.245594) (Poisson{Float64},0.785297)
+ (Poisson{Float64},0.264822) (Poisson{Float64},0.734856) (Poisson{Float64},0.667099)
]))
-const contact_time_distribution_matrix = convert(Array{initial_rest_type},[Geometric() for i in 1:AgentDemographic.size + 1, j in 1:AgentDemographic.size + 1])
+const contact_time_distribution_matrix = convert(Array{initial_rest_type},[Geometric() for i in 1:(AgentDemographic.size-1), j in 1:(AgentDemographic.size-1)])
+
diff --git a/CovidAlertVaccinationModel/test/runtests.jl b/CovidAlertVaccinationModel/test/runtests.jl
index cb4d959..56772e1 100644
--- a/CovidAlertVaccinationModel/test/runtests.jl
+++ b/CovidAlertVaccinationModel/test/runtests.jl
@@ -3,12 +3,12 @@ using RandomNumbers.Xorshifts
using Test
using ThreadsX
import StatsBase.mean
-model_sizes = [(100,1),(1000,1),(1000,3),(5000,3)]
+model_sizes = [100,1000,1000,5000]
vaccination_strategies = [vaccinate_uniformly!]
vaccination_rates = [0.000,0.005,0.01,0.05]
infection_rates = [0.01,0.05,0.1]
agent_models = ThreadsX.map(model_size -> AgentModel(model_size...), model_sizes)
-dem_cat = AgentDemographic.size + 1
+dem_cat = AgentDemographic.size -1
--
GitLab