# Epidemiological Models for Influenza and COVID-19

Posted 3 years ago
48012 Views
|
37 Replies
|
32 Total Likes
|
37 Replies
Sort By:
Posted 3 years ago
 Thanks, your work has helped me to understand this better. Instead of relatively closed compartments, have you done any modeling of spread of disease through contact networks? Quarantine measures work backwards through networks from index cases to exposed cases, with the ultimate goal of limiting the disease space to the size of infected plus exposed network and hopefully eliminating contact with this space to the healthy network space. Information about the disease can spread much more quickly than the disease. In the internet age there might be a role for social media networks to offer messaging to close contacts to limit the spread of even the common cold.
Posted 3 years ago
 No, no modeling on networks yet.You are quite right about the standard practice of tracing contact networks backwards in order to quarantine potentially exposed people. Unfortunately, those kind of data are not available for this pandemic. Martcheva gives a good description and model on pp. 230-234 of her book An Introduction to Mathematica Epidemiology.Social media can play a role, and many people are posting constantly this week about steps to take to "flatten the curve".
Posted 3 years ago
 -- you have earned Featured Contributor Badge Your exceptional post has been selected for our editorial column Staff Picks http://wolfr.am/StaffPicks and Your Profile is now distinguished by a Featured Contributor Badge and is displayed on the Featured Contributor Board. Thank you, keep it coming, and consider contributing to the The Notebook Archive!
Posted 3 years ago
 A more manageable set of background notebooks is now available, starting with:EpidemiologicalModelsForInfluenzaAndCOVID-19--part_1.nbThere is a TOC in each one to assist navigation among them.
Posted 3 years ago
 I greatly enjoyed your video presentation. It seemed that the size of the susceptible population made a big difference, which you explained as mixing being inefficient, and limiting by quarantine. I have been playing with a much simpler logistic model mainly to capture the day to day dynamics with the assumption that the epidemic can only be stopped with effective quarantine: Logistic Quarantine ModelThe assumption is that the quarantine process is competing with the infection rate, but that the quarantine set can grow faster than the set of infected and exposed, such that the infected and exposed population is eventually contained in the quarantine set and that set becomes the limit to the susceptible set. I think if you incorporated this into your model, it could give a better estimate of the size of the susceptible set. The model has worked with China. The graph below uses a very simple measure on the reported numbers. The daily differences of the log of the daily accumulated case numbers, i.e. the continuous growth rate, which should be horizontal if the number of cases grows exponentially.The blue dots are the daily differences of log case number predicted by the logistic model on the last date of the data. Italy is shown below.Once the actual data closes in on the predicted, the quarantine has successfully limited the susceptible population. I suspect with refinement, your model might eventually be sensitive enough to discover other dynamics. For instance did the Chinese method of quarantine to put infected people together in close quarters increase the death rate because, while the people were infected, they may not have developed much of an immune response, and the quarantine method permitted reinfection causing higher viral titers. Household quarantine also may guarantee higher infection rates of household members than if the infected person were removed from the household. In the community where I live in rural central Arizona, the population has decided to remove itself from the general population. The nearest reported infection is a hundred miles away. The Walmart and two supermarkets have been emptied and people are staying at home. Fortunately we are surrounded by three million acres of National Forest, so there is plenty of room to go hiking.
Posted 3 years ago
 Congratulations! This is a very interesting and useful introduction into a highly topical issue. Following the yesterday's fascinating presentation of the models, I have tried unsuccessfully to run the notebook (Mathematica 12.1, Win 10 64). I am getting many error messages. Please advise! Attachments:
Posted 3 years ago
 Thanks for your interest. Did you also download the package that's mentioned at the beginning of the initialization section: https://www.wolframcloud.com/obj/rnachbar/Published/CompartmentalModeling.wl?The code for KineticCompartmentalModel is in the package and that's how the transitions are interpreted and the ODEs generated.I'll copy the notice about the package to the top of the notebook so it's more obvious.
Posted 3 years ago
Posted 3 years ago
 Very good modeling. I was wondering if you could take into account herd immunity— as the number of recovered individuals goes up, the transmission rate goes down. Once 70% of the population is exposed the rate of infection begins to drastically decrease. Hence it may take a long time for the last 30% to be become infectedAnother factor to consider with the data is the sensitivity of the test— at best it is 60% and in real world testing 40%.
Posted 3 years ago
 Hi Robert, This was an excellent introduction and great modeling. My question is this: It seems like the modeling is trying to determine a few parameters across the whole population, when in fact we know that the parameters vary by age group. Why not consider breaking down the model in to population groups of different ages, for instance 0-9 years, 10-19, 20-29, etc. In this way, factors such as total population, mortality, effective community size, etc, can be independently optimized. Then, the total number of infections, deaths, etc. can be summed over age groups. Cheers, -Kent
Posted 3 years ago
 Indeed, that is on my To Do list. I can easily get the age distribution for a country from Wolfram|Alpha, so that will help with the initial conditions. However, to estimate the values for an age-stratified parameter we really need the number of confirmed cases, recovered cases, and deaths broken down by age too. There are some data on deaths for various ages.See the University of Basel model for an example using age groups. They use estimates based on an analysis of the recent data from China, and are not fitting them. Their interface allows the user to adjust them, also.
Posted 3 years ago
 Running the EpidemiologicalModelsForInfluenzaAndCOVID-19--part_1.nb gives many errors. I have downloaded the CompartmentalModeling.wl file and put it in the notebook folder and load it inside the notebook also.Example of the first error: Running the code: {susceptibleSEIR, infectedSEIR, exposedSEIR, recoveredSEIR} = {\[ScriptCapitalS], \[ScriptCapitalI], \ \[ScriptCapitalE], \[ScriptCapitalR]} /. ParametricNDSolve[ Join[odesSEIR /. forceOfInfectionSEIR, {\[ScriptCapitalS][ 0] == \[ScriptCapitalN] - I0, \[ScriptCapitalE][0] == 0, \[ScriptCapitalI][0] == I0, \[ScriptCapitalR][0] == 0}], Head /@ varsSEIR, {t, 0, 100}, {\[ScriptCapitalN], I0, \[Beta], \[Zeta], \[Gamma]}]; I get the errors: ParametricNDSolve::dspar: 1.5 cannot be used as a parameter. ReplaceAll::reps: {ParametricNDSolve[{(\[ScriptCapitalE]^\[Prime])[t]==0. -\[Zeta] \[ScriptCapitalE][t]+1.5 \[ScriptCapitalI][t] \[ScriptCapitalS][t],(\[ScriptCapitalI]^\[Prime])[t]==0. +\[Zeta] \[ScriptCapitalE][t]-0.0478 \[ScriptCapitalI][t],(\[ScriptCapitalR]^\[Prime])[t]==0. +0.0478 \[ScriptCapitalI][t],(\[ScriptCapitalS]^\[Prime])[t]==0. -1.5 \[ScriptCapitalI][t] \[ScriptCapitalS][t],\[ScriptCapitalS][0]==-I0+\[ScriptCapitalN],\[ScriptCapitalE][0]==0,\[ScriptCapitalI][0]==I0,\[ScriptCapitalR][0]==0},{\[ScriptCapitalE],\[ScriptCapitalI],\[ScriptCapitalR],\[ScriptCapitalS]},{t,0,100},{\[ScriptCapitalN],I0,1.5,\[Zeta],0.0478}]} is neither a list of replacement rules nor a valid dispatch table, and so cannot be used for replacing. Set::shape: Lists {susceptibleSEIR,infectedSEIR,exposedSEIR,recoveredSEIR} and {\[ScriptCapitalS],\[ScriptCapitalI],\[ScriptCapitalE],\[ScriptCapitalR]}/. ParametricNDSolve[{(\[ScriptCapitalE]^\[Prime])[t]==0. -\[Zeta] \[ScriptCapitalE][t]+1.5 \[ScriptCapitalI][t] \[ScriptCapitalS][t],(\[ScriptCapitalI]^\[Prime])[t]==0. +\[Zeta] \[ScriptCapitalE][t]-0.0478 \[ScriptCapitalI][t],(\[ScriptCapitalR]^\[Prime])[t]==0. +0.0478 \[ScriptCapitalI][t],(\[ScriptCapitalS]^\[Prime])[t]==0. -1.5 \[ScriptCapitalI][t] \[ScriptCapitalS][t],\[ScriptCapitalS][0]==-I0+\[ScriptCapitalN],\[ScriptCapitalE][0]==0,\[ScriptCapitalI][0]==I0,\[ScriptCapitalR][0]==0},{\[ScriptCapitalE],\[ScriptCapitalI],\[ScriptCapitalR],\[ScriptCapitalS]},{t,0,100},{\[ScriptCapitalN],I0,1.5,\[Zeta],0.0478}] are not the same shape. 
Posted 3 years ago
 The package should be loaded automatically when the initialization cells are evaluated.The Cloud notebook was accidentally edited yesterday, and I thought I had restored the changes to the original state, but I might have missed something.I will test the whole notebook now. Thanks for reporting the problem.
Posted 3 years ago
 Apparently you assigned the value 1.5 to the symbol [Beta] and 0.0478 to [Gamma]. I can reproduce your error messages when I do that. Also, if you look at the second message, you can see that the last argument to ParametricNDSolve is {\[ScriptCapitalN], I0, 1.5, \[Zeta], 0.0478}, which shows numeric values instead of symbols for [Beta] and [Gamma].It would probably be best to quit and restart the Kernel to fix the problem.Please let me know if you did not assign values to [Beta] and [Gamma], because that would suggest that the Manipulates are "leaking" dynamic assignments, which I made every effort to code against.
Posted 3 years ago
 Don't worry everything are ok, my bad. :) After including Clear["Global*"] at the start, everything runs fine. Probably some leftovers from previous notebook.Thanks for all of these very very good project!
Posted 3 years ago
 thank you! - where can we find the notebook part_2 and other parts if any? thank you so much Gerhard
Posted 3 years ago
Posted 3 years ago
 thank you- Gerhard
Posted 3 years ago
 Hi, I hope you are okay. I want to applied your code to my data, How can I upload my data for it could have the same format as your data?
Posted 3 years ago
 fitDataWDR returns a list of 3 data series. Each data series is a list of pairs of values of the form {time, number}. The first data series in the number of confirmed cases (not the cumulative number), the second is the number of recovered cases (which happens to be cumulative because it's a sink), and the third series is the number of deaths (also cumulative because it's a sink). For example the first 5 elements of each series for Beijing is In[18]:= fitData = fitDataWDR[Entity["AdministrativeDivision", {"Beijing", "China"}], "21 Jan 2020", "dateRange" -> All]; In[19]:= Take[#, 5] & /@ fitData Out[19]= {{{1., 14}, {2., 22}, {3., 35}, {4., 39}, {5., 66}}, {{1., 0}, {2., 0}, {3., 1}, {4., 2}, {5., 2}}, {{1., 0}, {2., 0}, {3., 0}, {4., 0}, {5., 0}}} If you don't have data for one of those series, then leave it out of the list of series, and also leave the corresponding model out of the list of models given to sumSquaredError.Hope this helps!
Posted 3 years ago
 When I ran the code on your last reply I got the following error. Can I work with the data from Argentina? Attachments:
Posted 3 years ago
 Sorry for the delay in responding.In your notebook, you have In[19] := Take[#, 5] & /@ fitData which looks like a copy & paste error. Remove the "In[19]:=" from the code.
Posted 3 years ago
 Once again: A very impresive project! Statistical fit does not work though (see attached notebook). Please advise! Attachments:
Posted 3 years ago
 I had updated the package last week, but unfortunately missed two necessary changes needed in the full version notebook.There is an updated version of that notebook in the Cloud now. It has today's date for the most recent revision.
Posted 3 years ago
 Any pointer to that notebook?
Posted 3 years ago
Posted 3 years ago
 Thank you. I am testing them.
Posted 3 years ago
 This would be a great example of Dynamic.https://public.flourish.studio/visualisation/1712761/?fbclid=IwAR1YOkUVgzFG3QkaflyZZNoalYXxl59Kv-OPLoUUWwpYjxbM9fDinveDChk
Posted 3 years ago
 Many thanks for your contribution. I Excellent work. But what would happen if used with QuantileRegression.m?
Posted 3 years ago
 Thank you for your suggestion. At the moment, the bigger issue is finding the right structural form for the model, that is, the compartments and their connections, and adequately mapping them to the available data. Once that is done, then perhaps using quantile regression might provide more robust estimates of the parameters for more reliable predictions.There is also an issue with parameter identifiability, regardless of the estimation method. This will take a certain amount of domain knowledge to properly choose which parameters to hold fixed and which to estimate.
Posted 3 years ago
 Hi Robert,I am trying access the covid19 data for Iran by using the following fitData = fitDataWDR[Entity["Country", {"Iran"}], "1 March 2020", "dateRange" -> All]; but getting following error : fitDataWDR::locnf: Locale Entity[Country,{Iran}] not found. Could you please tell me what wrong am I doing here. Regards
Posted 3 years ago
 There is a syntax error in the Entity specification. Use fitData = fitDataWDR[Entity["Country", "Iran"], "1 March 2020", "dateRange" -> All]; 
Posted 3 years ago
 Have you considered adding a quarantine? components to the SEIsaRD model?It seems that would be the next logical model improvement step, especially for Chinese data which had a very aggressive Quarantine implemented. Then add multiple different population cohorts (child, young adult, middle age adult, old adult, compromised health-for all age groups).Of course that is a lot more work!
Posted 3 years ago
 Thanks for a really nice post and nice work! I've made my own model that includes quarantine but instead of creating a compartment (state variable) for quarantined I simply made beta a time varying function that drops a certain level once quarantine is started. I guess that's a good enough approximation even though it might not catch all dynamics.Now my question. Looking at e.g. Hubai outbreak you get beta = 0.79 and gamma = 0.04, which would correspond to an R0 of =19.8 (beta/gamma). But I've seen work estimating R0 to between 2 and 2.5 for Corona. Why the discrepancy?Regards Peter Aronsson
Posted 2 years ago
 Hi, I have a problem with notebook 1, I already downloaded the package and put it in the same directory
Posted 2 years ago
 Indeed you do! How did you "download" it?Using the "Download" button in the Wolfram Cloud, using the "Open from Cloud..." button on the Desktop Mathematica Welcome screen, or using File > Open from Wolfram Cloud... menu work OK.The only time I've seen this behavior is with NotebookOpen["https://www.wolframcloud.com/obj/rnachbar/Published/\ EpidemiologicalModelsForInfluenzaAndCOVID-19--part_1.nb"] from the desktop.Try one of the three methods I first mentioned.PS. I see that you're using Mathematica 11.3. Some things may not render correctly when the notebook is first opened (e.g., some cell styles are not defined), but reevaluating inputs should take care of that.
Posted 2 years ago
 Congratulations! This post was featured in the Wolfram Blog Wolfram Community Takes on Current Events: COVID-19 Data, Green Spaces, Homeschool Puzzles and More. We are looking forward to your future contributions.