Commit 5779702d authored by Peter Jentsch's avatar Peter Jentsch
Browse files

manually fitting lol

parent 59744455
...@@ -4,16 +4,16 @@ function approximate_mixing_matricies() ...@@ -4,16 +4,16 @@ function approximate_mixing_matricies()
sol = abm(p,nothing) sol = abm(p,nothing)
mean_mixing = zeros(3,3) mean_mixing = zeros(3,3)
display(p.O_distribution_shift) display(p.O_distribution_shift)
for (node_i,demo_i) in enumerate(sol.demographics)
for g_list in sol.inf_network.graph_list for g_list in sol.inf_network.graph_list
for g in g_list for g in g_list
# display(g.mixing_edges.weights_dict) # display(g.mixing_edges.weights_dict)
for node_j in neighbors(g,node_i) for e in keys(g.mixing_edges.weights_dict)
demo_j = sol.demographics[node_j]; demo_i = sol.demographics[e.a]
demo_j = sol.demographics[e.b]
mean_mixing[Int(demo_i), Int(demo_j)] += 1#get_weight(g,GraphEdge(node_i,node_j)) /2 mean_mixing[Int(demo_i), Int(demo_j)] += 1#get_weight(g,GraphEdge(node_i,node_j)) /2
mean_mixing[Int(demo_j), Int(demo_i)] += 1#get_weight(g,GraphEdge(node_i,node_j)) /2
end end
end end
end
end end
mean_mixing ./= (500 .* length.(sol.index_vectors) * 2) mean_mixing ./= (500 .* length.(sol.index_vectors) * 2)
display(mean_mixing) display(mean_mixing)
......
...@@ -102,6 +102,7 @@ function create_mixing_edges(demographic_index_vectors,mixing_matrix,weights_dis ...@@ -102,6 +102,7 @@ function create_mixing_edges(demographic_index_vectors,mixing_matrix,weights_dis
contact_array[j,i] = GraphEdge.(stubs_i,stubs_j) contact_array[j,i] = GraphEdge.(stubs_i,stubs_j)
end end
return MixingEdges(tot,contact_array,weights_distribution_matrix) return MixingEdges(tot,contact_array,weights_distribution_matrix)
end end
......
...@@ -27,7 +27,7 @@ function get_parameters() ...@@ -27,7 +27,7 @@ function get_parameters()
immunization_delay = 14, immunization_delay = 14,
immunization_begin_day =60, immunization_begin_day =60,
infection_introduction_day = 180, infection_introduction_day = 180,
O_distribution_shift = 0.2, O_distribution_shift = 1.0,
) )
return params return params
end end
...@@ -94,11 +94,11 @@ mutable struct ModelSolution{T,InfNet,SocNet,WSMixingDist,RestMixingDist} ...@@ -94,11 +94,11 @@ mutable struct ModelSolution{T,InfNet,SocNet,WSMixingDist,RestMixingDist}
adjust_distributions_mean!( md[1:3,3],params.O_distribution_shift) adjust_distributions_mean!( md[1:3,3],params.O_distribution_shift)
adjust_distributions_mean!( md[3,1:2],params.O_distribution_shift) #dont shift OO twice adjust_distributions_mean!( md[3,1:2],params.O_distribution_shift) #dont shift OO twice
end end
map_symmetrize(m_tuple) = map(md -> symmetrize_means(pop_sizes,md), m_tuple) map_symmetrize(m_tuple) = map(md -> symmetrize_means(pop_sizes,md), m_tuple)
ws_matrix_tuple = map_symmetrize(ws_mixing_tuple_preshift) ws_matrix_tuple = map_symmetrize(ws_mixing_tuple_preshift)
rest_matrix_tuple = map_symmetrize(rest_mixing_tuple_preshift) rest_matrix_tuple = map_symmetrize(rest_mixing_tuple_preshift)
# display(mean.(rest_matrix_tuple.daily))
is_app_user = app_users(demographics,params.app_user_fraction) is_app_user = app_users(demographics,params.app_user_fraction)
app_user_list = zeros(nodes) app_user_list = zeros(nodes)
......
...@@ -97,33 +97,29 @@ function fit_distribution_parameters(p_tuple) ...@@ -97,33 +97,29 @@ function fit_distribution_parameters(p_tuple)
) )
) )
function cost(p) function cost(p)
p_mergedbeta = merge(p_tuple_adjust, (β_m = p[2],β_o = p[2], )) #assume one beta output,model = solve_w_parameters(p_tuple_adjust, p_names,p)
output,model = solve_w_parameters(p_mergedbeta, p_names,p)
target_ymo_vac = ymo_vac .* sum(vaccination_data[1:end]) .* length.(model.index_vectors)
ymo_vaccination_ts = output.daily_immunized_by_age ymo_vaccination_ts = output.daily_immunized_by_age
normalize_by_pop(v) = v./length.(model.index_vectors) total_postinf_vaccination = sum.(eachrow(ymo_vaccination_ts[:,181:end]))
total_postinf_vaccination = sum.(eachrow(ymo_vaccination_ts[:,180:end]))
final_size = sum.(eachrow(output.daily_unvac_cases_by_age)) final_size = sum.(eachrow(output.daily_unvac_cases_by_age))
total_preinf_vaccination = sum.(eachrow(ymo_vaccination_ts[:,1:180]))
target_final_size = ymo_attack_rate .* length.(model.index_vectors) target_final_size = ymo_attack_rate .* length.(model.index_vectors)
target_ymo_vac = ymo_vac .* sum(vaccination_data[1:4]) .* length.(model.index_vectors) target_preinf_vac = ymo_vac .* sum(vaccination_data[1:4]) .* length.(model.index_vectors)
target_postinf_vac = ymo_vac .* sum(vaccination_data[5:end]) .* length.(model.index_vectors)
ymo_vaccination_ts = output.daily_immunized_by_age
total_preinfection_vaccination = sum.(eachrow(ymo_vaccination_ts))
return sum((normalize_by_pop(total_preinfection_vaccination .- target_ymo_vac)).^2) display((final_size,target_final_size))
+ sum(normalize_by_pop((total_postinf_vaccination .- target_ymo_vac)).^2) display((total_preinf_vaccination,target_preinf_vac))
+ 5*sum(normalize_by_pop((final_size .- target_final_size)).^2) display((total_postinf_vaccination,target_postinf_vac))
# display(final_size) display(sum.(eachrow(ymo_vaccination_ts)) ./length.(model.index_vectors))
# display( sum.(eachrow(ymo_vaccination_ts[:,1:end])))
# display( length.(model.index_vectors)) return sum((total_preinf_vaccination .- target_preinf_vac).^2)
# return sum(normalize_by_pop((final_size .- target_final_size)).^2) + sum((total_postinf_vaccination .- total_postinf_vaccination).^2)
end + 5*sum((final_size .- target_final_size).^2)
# display(cost((0.0,0.0005,-3.0,-3.0,-0.2,0.25,0.25,0.0,0.1))) end
display(cost((0.0000,0.00048,0.0005,0.16,-1.30,-1.24,-0.8,0.35,0.35,0.35,0.2)))
# #
out = smc(priors,cost; verbose = true, nparticles = 800, parallel = true) # out = smc(priors,cost; verbose = true, nparticles = 600, parallel = true)
return NamedTuple{p_names}(ntuple(i -> out.P[i].particles,length(p_names))) return NamedTuple{p_names}(ntuple(i -> out.P[i].particles,length(p_names)))
end end
function fit_parameters(default_parameters) function fit_parameters(default_parameters)
...@@ -133,8 +129,8 @@ function fit_parameters(default_parameters) ...@@ -133,8 +129,8 @@ function fit_parameters(default_parameters)
fit_dist_parameters_path = joinpath(PACKAGE_FOLDER,"abm_parameter_fits","fit_dist_parameters.dat") fit_dist_parameters_path = joinpath(PACKAGE_FOLDER,"abm_parameter_fits","fit_dist_parameters.dat")
# output = fit_distribution_parameters(default_parameters) output = fit_distribution_parameters(default_parameters)
# serialize(fit_all_parameters_path,output) serialize(fit_all_parameters_path,output)
fitted_parameter_tuple = deserialize(fit_all_parameters_path) fitted_parameter_tuple = deserialize(fit_all_parameters_path)
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment