InstituteforDiseaseModeling / covasim

COVID-19 Agent-based Simulator (Covasim): a model for exploring coronavirus dynamics and interventions
https://covasim.org
MIT License
250 stars 223 forks source link

Interaction between the pop_scale, pop_infected, and rescale parameters #397

Open AndrewC19 opened 1 year ago

AndrewC19 commented 1 year ago

Hi,

I am slightly confused about how population scaling is supposed to interact with the initally infectious population. This seems to change depending on whether dynamic rescaling is enabled or disabled.

For example, if I have:

  1. pop_size=5000
  2. pop_infected=100
  3. pop_scale=10 (i.e. I want the population to represent 50000 individuals)

Should the scaling factor also apply to pop_infected i.e. should this start the simulation with 100 infected individuals or 1000 infected individuals?

I couldn't find the answer in the documentation so I performed some tests and found that the answer appears to change depending on whether rescale is set to True or False.

The plots below highlight the discrepancy, which show two identical simulations: one with rescale=True and another with rescale=False. The first plot has 100 day-zero cumulative infections (scaling has not been applied to the initial infectious population) and the second plot has 200 day-zero cumulative infections (scaling has been applied to the initial infectious population):

image

image

Here's my code:

import covasim as cv

def no_pop_scale():
    sim_pars = {"pop_size": 50000,
                "pop_infected": 100,
                "pop_type": "hybrid",
                "n_days": 200}

    sim = cv.Sim(sim_pars)
    sim.run()
    print(sim.results["cum_infections"][0])
    sim.plot_result("cum_infections")

def with_pop_scale():
    scaling_factor = 10
    sim_pars = {"pop_size": 50000 / scaling_factor,
                "pop_infected": 100,
                "pop_type": "hybrid",
                "n_days": 200,
                "pop_scale": scaling_factor}
    sim = cv.Sim(sim_pars)
    sim.run()
    print(sim.results["cum_infections"][0])
    sim.plot_result("cum_infections")

def no_rescaling():
    sim_pars = {"pop_size": 50000,
                "pop_infected": 100,
                "pop_type": "hybrid",
                "n_days": 200,
                "rescale": False}
    sim = cv.Sim(sim_pars)
    sim.run()
    print(sim.results["cum_infections"][0])
    sim.plot_result("cum_infections")

def no_rescaling_with_pop_scale():
    scaling_factor = 20
    sim_pars = {"pop_size": 50000 / scaling_factor,
                "pop_infected": 100,
                "pop_type": "hybrid",
                "n_days": 200,
                "pop_scale": scaling_factor,
                "rescale": False}
    sim = cv.Sim(sim_pars)
    sim.run()
    print(sim.results["cum_infections"][0])
    sim.plot_result("cum_infections")

if __name__ == "__main__":
    no_pop_scale()
    with_pop_scale()
    no_rescaling()
    no_rescaling_with_pop_scale()

Regardless of whether dynamic rescaling is enabled or not, it seems odd that this would affect the reported initial conditions of the model (i.e. the number of cumulative infections, exposed individuals etc. on day zero).

My intuition is that dynamic rescaling should only provide an efficiency gain and should not interefere with the interpretation of the parameters/reporting of results, i.e. If I set initial infections to 100, pop size to 5000, and population scaling to a factor of 2, I anticipate the simulation to start from 100 exposed individual and 10000 total individuals whether scaling is used or not (or 200 exposed individuals and 10000 total if scaling should also be applied to the initally infectious population).

Can anyone help me to understand why pop_infected is only affected by pop_scale when rescale is set to False?

Thanks.