abm_sensitivity_analysis.jl

using CovidAlertVaccinationModel
using OnlineStats
using ThreadsX
using Plots
const samples = 5
##Univariate tests
const len = 4 #number of points to evaluate
gr()

const default_parameters = (
    sim_length = 600,
    num_households = 5000,
    I_0_fraction = 0.005,
    base_transmission_probability = 0.001,
    recovery_rate = 1/7,
    immunization_loss_prob = 0.0055, #mean time of 6 months
    π_base = -4.0,
    η = 0.0,
    κ = 1.0,
    ω = 0.001,
    ρ = [0.0,0.0,0.0],
    ω_en = 0.00,
    ρ_en =  [0.0,0.0,0.0],
    γ = 0.0,
    β = 5.0,
    notification_parameter = 0.001,
    vaccinator_prob = 0.2,
    app_user_fraction = 0.4,
    notification_threshold = 2,
    immunizing = true,
    immunization_delay = 14,
    immunization_begin_delay = 50,
    infection_introduction_delay = 100
)

#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)),
    (:recovery_rate, range(0.1, 0.5; length = len)),
    (:immunization_loss_prob, range(0.00, 0.05; length = len)),
    (:π_base, range(-4.5, -3.5;  length = len)),
    (:η, range(0.0, 0.01; length = len)),
    (:κ, range(0.5, 1.5; length = len)),
    (:ω, range(0.0, 0.01; length = len)),
    (:ω_en, range(0.0, 0.5; length = len)),
    (:γ, range(0.0, 0.5; length = len)),
    (:β, range(2, 30; length = len)),
    (:notification_parameter, range(0.00, 0.05; length = len)),
    (:app_user_fraction, range(0.05, 0.25; length = len)),
    (:notification_threshold, (1:len)),
    (:immunization_delay, [7,10,14,20]),
)

const univariate_path = "CovidAlertVaccinationModel/plots/univariate/" 
function univarate_test(variable, variable_range)
    parameter_range_list = [merge(default_parameters,NamedTuple{(variable,)}((value,))) for value in variable_range]
    solve_fn(p) = mean_solve(samples, p,DebugRecorder)

    univariate_outlist = ThreadsX.map(solve_fn, parameter_range_list)
    
    p = plot_model(variable,parameter_range_list,univariate_outlist,default_parameters.infection_introduction_delay,default_parameters.immunization_begin_delay)
    return p
end

if !ispath(univariate_path)
    mkdir(univariate_path)
end
function univariate_simulations()
    plt_list = ThreadsX.map(univarate_test_list) do ur
        out = univarate_test(ur...)
        display("done $(ur[1])")
        return out
    end
    for ((varname,_),p) in zip(univarate_test_list,plt_list)
        savefig(p,"$univariate_path/$varname.pdf")
    end
end
univariate_simulations()