Skip to content
Snippets Groups Projects
Commit 203b7a3d authored by Peter Jentsch's avatar Peter Jentsch
Browse files

fixed a breakage moving mixing distributions to IntervalsModel from CovidAlterVaccinationModel

parent 54c67019
No related branches found
No related tags found
No related merge requests found
......@@ -10,15 +10,6 @@ function get_canada_demographic_distribution()::Vector{Float64}
source_bins = map(parse_string_as_float_pair, df[:,:demographic_data_bins])
return [sum(binned_data[1:5]),sum(binned_data[6:13]),sum(binned_data[14:end])]
end
function get_distancing_data_sets()::NTuple{3,Vector{Float64}}
f = readdlm(joinpath(PACKAGE_FOLDER,"data/csv/distancing_data.csv"), ',')
df = filter(row -> row[:country_region_code] == "CA" && row[:sub_region_1] == "Ontario", DataFrame([f[1,i] => f[2:end, i] for i = 1:length(f[1,:])]))
distancing_data = df[:,:retail_and_recreation_percent_change_from_baseline] ./ 100
workplace_data = df[:,:workplaces_percent_change_from_baseline] ./ 100
residental_data = df[:,:residential_percent_change_from_baseline] ./ 100
return (workplace_data, distancing_data,residental_data)
end
function get_canada_case_fatality()::Tuple{Vector{Tuple{Float64,Float64}},Vector{Float64}}
f = readdlm(joinpath(PACKAGE_FOLDER,"data/csv/case_fatality_data.csv"), ',')
......@@ -62,62 +53,6 @@ function get_household_data_proportions()
end
function parse_string_as_float_pair(s)::Tuple{Float64,Float64}
parsed = strip(s, ['(',')']) |> s -> split(s, ",")
if length(parsed) > 2
error("input tuple too long")
end
return (parse(Float64, parsed[begin]), parse(Float64, parsed[end]))
end
bin_str_to_int(bin_str) = occursin('<',bin_str) ? 1 : parse(Int, strip(bin_str, 's')) #this function parses the bin formatting in the input csv
bin_to_index(bin) = findfirst(((lower,upper),) -> lower <= bin < upper, age_bins)
function parse_cases_data()::Tuple{Vector{Date},Matrix{Float64}} #this type declaration probably not needed
f = readdlm(joinpath(PACKAGE_FOLDER,"data/csv/COVID_ontario_data.csv"), ',')
df = DataFrame([f[1,i] => f[2:end, i] for i = 1:length(f[1,:])])
dropmissing!(df,:Accurate_Episode_Date)
map!(x -> Date(x),df[!,:Accurate_Episode_Date],df[!,:Accurate_Episode_Date])
sort!(df,:Accurate_Episode_Date)
#there are some cases before jan 1 2020 but they are neglible so we assume that cases start here
begin_date = Date(2020,1,1)
cases_by_date = groupby(df,:Accurate_Episode_Date)
#Ontario health says that case data with 2 weeks of current date are subject to change
end_date = Date(maximum(df[!,:Accurate_Episode_Date])) - Day(14)
dates = collect(begin_date:Day(1):end_date)
infection_data = zeros(Float64,length(dates), length(age_bins))
for (i,day) in enumerate(dates)
cases_on_date = get(cases_by_date, (Accurate_Episode_Date = day,),nothing)
if !isnothing(cases_on_date)
cases_by_age_group = groupby(cases_on_date,:Age_Group)
for age_group_key in keys(cases_by_age_group)
cases = cases_by_age_group[age_group_key]
bin_str = age_group_key[:Age_Group]
num_cases = nrow(cases)
if length(bin_str)>0 && !occursin("UNKNOWN",bin_str)
bin = bin_str_to_int(bin_str)
bin_index = bin_to_index(bin)
infection_data[i,bin_index] = num_cases
else
#account for missing age data?
end
end
end
end
return dates, infection_data
end
using CSV
function load_mixing_matrices()
......
......@@ -13,24 +13,6 @@ function adjust_distributions_mean(distribution_matrix::AbstractArray{T},mean_sh
end
end
function ZeroPoisson(α,λ)
return ZWDist(α,Poisson(λ))
end
function ZeroGeometric(α,p)
return ZWDist(α,Geometric(p))
end
function symmetrize_means(population_sizes,mixing_matrix::Matrix{<:ZWDist})
mixing_means = mean.(mixing_matrix)
symmetrized_means = zeros((length(population_sizes),length(population_sizes)))
for i in 1:length(population_sizes), j in 1:length(population_sizes)
symmetrized_means[i,j] = 0.5*(mixing_means[i,j] + population_sizes[j]/population_sizes[i] * mixing_means[j,i])
end
return map(((dst,sym_mean),)->from_mean(typeof(dst),dst.α,sym_mean), zip(mixing_matrix,symmetrized_means))
end
function symmetrize_means(population_sizes,mixing_matrix::Matrix{<:Distribution})
mixing_means = mean.(mixing_matrix)
symmetrized_means = zeros((length(population_sizes),length(population_sizes)))
......@@ -72,31 +54,6 @@ const unemployment_matrix = alpha_matrix(
)
)
function make_workschool_mixing_matrix()
#all geometric with means given
ws_mixing = map(t->from_mean(t...),[
(Geometric{Float64}, 4.104848) (Geometric{Float64},2.568782) (Geometric{Float64},0.017729)
(Geometric{Float64}, 0.975688) (Geometric{Float64},5.057572) (Geometric{Float64},0.021307)
(Geometric{Float64},0.001937) (Geometric{Float64},0.00722) (Geometric{Float64}, 0.022134)
])
#symmetrize WS mixing with respect to the population proportions in the data
ws_mixing_w_unemployment_symmetrized = symmetrize_means(get_household_data_proportions(),ws_mixing)
#define a function that adjusts the means of W according to the unemployment_matrix
ws_adjust_mean(W) = (5/7) .* (1 .- unemployment_matrix) ./ ( mean.(W) + (5/7) .* (1 .- unemployment_matrix))
#create a zero weighted distribution where the zero-weight is given by the unemployment_matrix, and the non-zero weight is given by ws_mixing, symmetrized, with means adjust upwards
ws_mixing_w_unemployment_symmetrized_weekday_adjusted = ZWDist.(unemployment_matrix,Geometric.(ws_adjust_mean(ws_mixing_w_unemployment_symmetrized)))
return ws_mixing_w_unemployment_symmetrized_weekday_adjusted
end
@views function ws_sample(age)
return rand.(initial_workschool_mixing_matrix[age,:]) * (rand(RNG) < (5/7))
end
const workschool_mixing, rest_mixing = load_mixing_matrices()
......
......@@ -10,7 +10,6 @@ using Random
using KernelDensity
using StatsBase
using Distributions
using CovidAlertVaccinationModel
using ZeroWeightedDistributions
const PACKAGE_FOLDER = dirname(dirname(pathof(IntervalsModel)))
using DataStructures:OrderedDict
......@@ -58,7 +57,7 @@ end
priors_list::Vector{Distribution},
dists::Vector{Distribution},
particles::Integer;
alpha = 0.95
alpha = 0.98
)
# Arguments
......@@ -85,7 +84,7 @@ Number of particles to obtain representing posterior distribution of the distrib
# Optional Arguments
- alpha = 0.95
- alpha = 0.98
Threshold adaptive rate, see `KissABC.smc` for more details.
"""
......
const rest_distributions = CovidAlertVaccinationModel.initial_rest_mixing_matrix
const rest_data = get_rest_data()
const kerneldensity_data_rest = (
Y = InterpKDE(kde(rest_data.Y.durs,weights = rest_data.Y.weights)),
......@@ -11,6 +9,14 @@ const data_rest = hcat(
map(i->pdf(kerneldensity_data_rest.M,i), 0:143),
map(i->pdf(kerneldensity_data_rest.O,i), 0:143),
)
using CovidAlertVaccinationModel: symmetrize_means, get_household_data_proportions, from_mean
const rest_distributions = symmetrize_means(get_household_data_proportions(),map(t->from_mean(t...),[
(Geometric{Float64},2.728177) (Geometric{Float64},1.382557) (Geometric{Float64},0.206362)
(Geometric{Float64},1.139072) (Geometric{Float64},3.245594) (Geometric{Float64},0.785297)
(Geometric{Float64},0.264822) (Geometric{Float64},0.734856) (Geometric{Float64},0.667099)
]))
"""
err_rest(p::Vector{Number},dist::Distribution)
......
......@@ -33,7 +33,7 @@ function err_ws(p,dist)
for _ in 1:sample_repeat
for age_sample in YOUNG:OLD
#Sample neighbourhood sizes. Note that this includes the weekend/weekday coinflip, implemented in CovidAlertVaccinationModel.
neighourhoods = CovidAlertVaccinationModel.ws_sample(age_sample)
neighourhoods = ws_sample(age_sample)
for age_j in YOUNG:OLD #for a given age_sample loop over possible contact ages
if neighourhoods[age_j] > 0
......@@ -57,3 +57,30 @@ function err_ws(p,dist)
end
return errsum./sample_repeat
end
function make_workschool_mixing_matrix()
#all geometric with means given
ws_mixing = map(t->from_mean(t...),[
(Geometric{Float64}, 4.104848) (Geometric{Float64},2.568782) (Geometric{Float64},0.017729)
(Geometric{Float64}, 0.975688) (Geometric{Float64},5.057572) (Geometric{Float64},0.021307)
(Geometric{Float64},0.001937) (Geometric{Float64},0.00722) (Geometric{Float64}, 0.022134)
])
#symmetrize WS mixing with respect to the population proportions in the data
ws_mixing_w_unemployment_symmetrized = symmetrize_means(get_household_data_proportions(),ws_mixing)
#define a function that adjusts the means of W according to the unemployment_matrix
ws_adjust_mean(W) = (5/7) .* (1 .- unemployment_matrix) ./ ( mean.(W) + (5/7) .* (1 .- unemployment_matrix))
#create a zero weighted distribution where the zero-weight is given by the unemployment_matrix, and the non-zero weight is given by ws_mixing, symmetrized, with means adjust upwards
ws_mixing_w_unemployment_symmetrized_weekday_adjusted = ZWDist.(unemployment_matrix,Geometric.(ws_adjust_mean(ws_mixing_w_unemployment_symmetrized)))
return ws_mixing_w_unemployment_symmetrized_weekday_adjusted
end
@views function ws_sample(age)
return rand.(initial_workschool_mixing_matrix[age,:]) * (rand(RNG) < (5/7))
end
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment