Skip to content

Commit

Permalink
Merge pull request #35 from fccoelho/master
Browse files Browse the repository at this point in the history 
First round of inference using PyMC3 with Nuts
  • Loading branch information
@maxbiostat
maxbiostat committed Mar 18, 2021
2 parents 50fc539 + 7de369d commit ef1149d
Show file tree
Hide file tree
Showing 5 changed files with 1,689 additions and 1 deletion.
730 changes: 730 additions & 0 deletions code/New_notebooks/Data Simulation.ipynb
View file

Large diffs are not rendered by default.

596 changes: 596 additions & 0 deletions code/New_notebooks/Inference PyMC3.ipynb
View file

Large diffs are not rendered by default.

239 changes: 239 additions & 0 deletions code/New_notebooks/Inference cmdstan.ipynb
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,239 @@
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Inference comparing different models and simulated data"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"ExecuteTime": {
"end_time": "2021-03-17T15:05:50.552831Z",
"start_time": "2021-03-17T15:05:50.011627Z"
}
},
"outputs": [],
"source": [
"import numpy as np\n",
"from matplotlib import pyplot as plt\n",
"from scipy.integrate import odeint\n",
"import cmdstanpy"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Full SIR model\n",
"Let's start with the SIR model (S(t),I(t)), with LogNormal data. See Data Simulation Notebook for details on how the data was generated. "
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"ExecuteTime": {
"end_time": "2021-03-17T15:10:41.576000Z",
"start_time": "2021-03-17T15:10:41.565578Z"
}
},
"outputs": [],
"source": [
"data = np.load('logNSIRdata.npy')"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"ExecuteTime": {
"end_time": "2021-03-17T15:10:45.139430Z",
"start_time": "2021-03-17T15:10:45.118708Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"array([[0.65719681, 0.01228651],\n",
" [0.80807178, 0.01356656],\n",
" [0.89714913, 0.02439541],\n",
" [0.81356525, 0.01492339],\n",
" [0.80193999, 0.02341071],\n",
" [0.843973 , 0.02940851],\n",
" [0.88415507, 0.06102952],\n",
" [0.78496068, 0.06830704],\n",
" [0.95299406, 0.03634341],\n",
" [0.82480258, 0.05400302],\n",
" [0.87924861, 0.07700426],\n",
" [0.69566286, 0.09291214],\n",
" [0.84106409, 0.2931043 ],\n",
" [0.84037437, 0.10659906],\n",
" [0.78747029, 0.15183002],\n",
" [0.8413808 , 0.17485457],\n",
" [0.69513859, 0.27376009],\n",
" [0.8957417 , 0.204253 ],\n",
" [0.54047104, 0.26318554],\n",
" [0.70447339, 0.29281707],\n",
" [0.55916381, 0.30169365],\n",
" [0.74555253, 0.2288304 ],\n",
" [0.42922743, 0.29885803],\n",
" [0.42492554, 0.4964244 ],\n",
" [0.36974463, 0.3714413 ],\n",
" [0.33873785, 0.51803187],\n",
" [0.25375758, 0.63782453],\n",
" [0.30105992, 0.45043982],\n",
" [0.16212611, 0.3671034 ],\n",
" [0.20123444, 0.61551334],\n",
" [0.17485733, 0.83546652],\n",
" [0.12368927, 0.58877006],\n",
" [0.13517663, 0.6402482 ],\n",
" [0.12336577, 0.47835116],\n",
" [0.06223319, 0.85452369],\n",
" [0.08016665, 0.47882747],\n",
" [0.05191069, 0.77020783],\n",
" [0.09964947, 0.86436132],\n",
" [0.06974944, 0.51175903],\n",
" [0.05422395, 0.58962803],\n",
" [0.03459974, 0.63428953],\n",
" [0.03512599, 0.24933931],\n",
" [0.0523298 , 0.90175561],\n",
" [0.02541269, 0.97546729],\n",
" [0.02417837, 0.49537834],\n",
" [0.02280534, 0.88672294],\n",
" [0.01652361, 0.39648975],\n",
" [0.02000476, 0.3898582 ],\n",
" [0.01568821, 0.50271297],\n",
" [0.01683631, 0.72305858],\n",
" [0.01073431, 0.32665781],\n",
" [0.01146334, 0.5589468 ],\n",
" [0.01299654, 0.72767495],\n",
" [0.00930209, 0.4327625 ],\n",
" [0.00619795, 0.60399564],\n",
" [0.00982066, 0.63619781],\n",
" [0.00880199, 0.34508012],\n",
" [0.00637482, 0.30628235],\n",
" [0.00453392, 0.25234106],\n",
" [0.00440076, 0.41386706],\n",
" [0.00406945, 0.48632646],\n",
" [0.0035637 , 0.47785857],\n",
" [0.00515154, 0.50758805],\n",
" [0.00546414, 0.48795409],\n",
" [0.00372785, 0.22840271],\n",
" [0.00309657, 0.36058463],\n",
" [0.00248 , 0.27360166],\n",
" [0.00218348, 0.15912421],\n",
" [0.0027339 , 0.2814214 ],\n",
" [0.00343316, 0.33153157],\n",
" [0.00308386, 0.29411404],\n",
" [0.00135047, 0.21676085],\n",
" [0.00184293, 0.29936522],\n",
" [0.00163758, 0.33561504],\n",
" [0.00165209, 0.44955844],\n",
" [0.00229176, 0.28955305],\n",
" [0.00123067, 0.25554528],\n",
" [0.00118941, 0.23016772],\n",
" [0.00099515, 0.31081312]])"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.6"
},
"latex_envs": {
"LaTeX_envs_menu_present": true,
"autoclose": false,
"autocomplete": true,
"bibliofile": "biblio.bib",
"cite_by": "apalike",
"current_citInitial": 1,
"eqLabelWithNumbers": true,
"eqNumInitial": 1,
"hotkeys": {
"equation": "Ctrl-E",
"itemize": "Ctrl-I"
},
"labels_anchors": false,
"latex_user_defs": false,
"report_style_numbering": false,
"user_envs_cfg": false
},
"toc": {
"base_numbering": 1,
"nav_menu": {},
"number_sections": true,
"sideBar": true,
"skip_h1_title": false,
"title_cell": "Table of Contents",
"title_sidebar": "Contents",
"toc_cell": false,
"toc_position": {},
"toc_section_display": true,
"toc_window_display": false
},
"varInspector": {
"cols": {
"lenName": 16,
"lenType": 16,
"lenVar": 40
},
"kernels_config": {
"python": {
"delete_cmd_postfix": "",
"delete_cmd_prefix": "del ",
"library": "var_list.py",
"varRefreshCmd": "print(var_dic_list())"
},
"r": {
"delete_cmd_postfix": ") ",
"delete_cmd_prefix": "rm(",
"library": "var_list.r",
"varRefreshCmd": "cat(var_dic_list()) "
}
},
"types_to_exclude": [
"module",
"function",
"builtin_function_or_method",
"instance",
"_Feature"
],
"window_display": false
}
},
"nbformat": 4,
"nbformat_minor": 4
}
63 changes: 63 additions & 0 deletions code/New_notebooks/sir.d
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,63 @@
/+dub.json:
{
"name": "sir_stoch",
"authors": [
"Flávio Codeço Coelho"
],
"description": "Example usage of SIR model from epistochmodels",
"copyright": "Copyright © 2020, Flávio Codeço Coelho",
"license": "GPL-3",
"dependencies": {
"epistochmodels": "*"
}
}
+/
import std.datetime.stopwatch;
import std.stdio;
import std.file;
import std.algorithm;
import std.range;
import core.exception : RangeError;
import std.typecons : tuple, Tuple;
import models;

/**
Save the simulation as a CSV file
*/
void save(string filename, string[] varnames, Tuple!(double[], int[][]) data)
{
import std.array: join;
import std.conv: text;
File outf = File(filename, "w");
auto head = join(varnames, ",");
outf.writeln("t,", head);
foreach (i, int[] row; data[1])
{
try
{// writes t on column 0
outf.writeln(text(data[0][i]), ",", row.map!text.join(","));
}
catch (RangeError e)
{
writefln("Some error %s: %s", e, row);
}
}
outf.close();
}

void main(){
double beta = 0.9;
double gam = 0.1;
int N = 100;
int I0 = 10;
double tf = 1000;
auto sw = StopWatch(AutoStart.no);
auto model = new SIR(N, beta, gam);
model.initialize(N-I0, I0, 0);
sw.start();
auto sim = model.run(0, tf);
sw.stop();
save("sir_stoch.csv",["S","I","R"],sim);
writefln("Time of the SIR run with N=%s: %s", N, sw.peek());
writefln("Number of steps: %s", sim[0].length);
}
62 changes: 61 additions & 1 deletion code/old_notebooks/Inference.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -569,7 +569,67 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.6.7"
"version": "3.8.6"
},
"latex_envs": {
"LaTeX_envs_menu_present": true,
"autoclose": false,
"autocomplete": true,
"bibliofile": "biblio.bib",
"cite_by": "apalike",
"current_citInitial": 1,
"eqLabelWithNumbers": true,
"eqNumInitial": 1,
"hotkeys": {
"equation": "Ctrl-E",
"itemize": "Ctrl-I"
},
"labels_anchors": false,
"latex_user_defs": false,
"report_style_numbering": false,
"user_envs_cfg": false
},
"toc": {
"base_numbering": 1,
"nav_menu": {},
"number_sections": true,
"sideBar": true,
"skip_h1_title": false,
"title_cell": "Table of Contents",
"title_sidebar": "Contents",
"toc_cell": false,
"toc_position": {},
"toc_section_display": true,
"toc_window_display": false
},
"varInspector": {
"cols": {
"lenName": 16,
"lenType": 16,
"lenVar": 40
},
"kernels_config": {
"python": {
"delete_cmd_postfix": "",
"delete_cmd_prefix": "del ",
"library": "var_list.py",
"varRefreshCmd": "print(var_dic_list())"
},
"r": {
"delete_cmd_postfix": ") ",
"delete_cmd_prefix": "rm(",
"library": "var_list.r",
"varRefreshCmd": "cat(var_dic_list()) "
}
},
"types_to_exclude": [
"module",
"function",
"builtin_function_or_method",
"instance",
"_Feature"
],
"window_display": false
}
},
"nbformat": 4,
Expand Down

0 comments on commit ef1149d

Please sign in to comment.