springer-math / Mathematics-of-Epidemics-on-Networks

Source code accompanying 'Mathematics of Epidemics on Networks' by Kiss, Miller, and Simon http://www.springer.com/us/book/9783319508047 . Documentation for the software package is at https://epidemicsonnetworks.readthedocs.io/en/latest/
MIT License
151 stars 61 forks source link

Question about gamma distributed incubation period #70

Closed jsheen closed 4 years ago

jsheen commented 4 years ago

Hi Joel,

I was wondering if it was possible to add a gamma distributed incubation period in an SEIR model using this library. For example, if I wanted to use the incubation period from: https://www.acpjournals.org/doi/10.7326/M20-0504 where the parameters of the gamma distributed incubation period are: rate = 1/0.948 and shape = 5.807 (mean ~5.5 days).

The idea I had for an approximation of this was when creating the "H" graph from https://epidemicsonnetworks.readthedocs.io/en/latest/functions/EoN.Gillespie_simple_contagion.html, set the global E->I rate as 1: H.add_edge('E', 'I', rate = 1, weight_label='expose2infect_weight')

Then draw the weights from the gamma distribution to create the incubation time dictionary: node_attribute_inc_dict = {node: 1 / np.random.gamma(shape=incperiod_shape, scale=1 / incperiod_rate, size=1)[0] for node in network.nodes()}

I don't think this should work exactly, but wanted to double check if there might be a better way to get at the gamma distributed incubation time.

Thanks in advance for any help, Justin

joelmiller commented 4 years ago

Hey Justin,

There is a way to do it exactly in a way that the E class doesn't transmit and the I class does, but it comes with the downside that the code doesn't distinguish the E and I classes in the data it records. It may be possible to get it to record the data for you.

This would use the fast_nonMarkov_SIR algorithm. You'd need to define a trans_and_rec_time_fxn function which is called when a node becomes infected. This function needs to calculate both the eventual recovery time and the eventual transmission times to all neighbors. So the way you would do this is to first calculate the time until the E->I transition, and then calculate the time until recovery and the transmissions. It would then return the combined time from initial S->E transition and the I->R transition as well as the times from S->E transition to transmitting to neighbors. If any transmission time is after recovery time, you don't do anything with those.

I think it is possible to edit the trans_and_rec_time_fxn function so that it records some data to a global variable which would allow you to include the E to I transition.

More details here: https://epidemicsonnetworks.readthedocs.io/en/latest/functions/EoN.fast_nonMarkov_SIR.html#EoN.fast_nonMarkov_SIR

An alternative would be to generate an equivalent fast_nonMarkov_SEIR function. I would definitely welcome having such a function added, but it would be some effort.

jsheen commented 4 years ago

Hi Joel,

Thanks a lot, that makes sense after thinking through it.

Let me take a look into potentially making a pull request for a fast_nonMarkov_SEIR function if the global variable to record "E" doesn't work out. It might be worth pursuing in the first place actually.

Thanks again, Justin