notebooks/sircn/r.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## SIR model on a configuration network in R using simecol"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Author*: Simon Frost\n",
    "\n",
    "*Date*: 2018-07-12"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Description"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Equations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### References"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Implementation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "library(simecol)\n",
    "library(reshape2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "sir.cn.ode <- new(\"odeModel\",\n",
    "  main = function(time, init, parms, ...){\n",
    "    with(as.list(c(init,parms)),{\n",
    "      dtheta <- -beta*theta+beta*(dpsi(theta,k)/dpsi(1,k))+gamma*(1-theta)\n",
    "      S <- psi(theta,k)\n",
    "      I <- 1-S-R\n",
    "      dR <- gamma*I\n",
    "      list(c(dtheta,dR))\n",
    "    })},\n",
    "  equations = list(),\n",
    "  parms = c(beta=0.1,gamma=0.05,k=5),\n",
    "  times = c(from=0,to=125,by=0.01),\n",
    "  init = c(theta=0.999,R=0),\n",
    "  solver = \"lsoda\"\n",
    ")\n",
    "poisgn <- list(\n",
    "  psi = function(theta,k){theta^k},\n",
    "  dpsi = function(theta,k){k*theta^(k-1)},\n",
    "  dpsi2 = function(theta,k){k*(k-1)*theta^(k-2)}\n",
    ")\n",
    "equations(sir.cn.ode) <- poisgn\n",
    "sir.cn.ode <- sim(sir.cn.ode)\n",
    "sir.cn.out <- out(sir.cn.ode)\n",
    "sir.cn.out$S <- sir.cn.out$theta^parms(sir.cn.ode)[[\"k\"]]\n",
    "sir.cn.out$I <- 1-sir.cn.out$S-sir.cn.out$R"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "sir.cn.out.long <- melt(as.data.frame(sir.cn.out),\"time\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Visualisation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "library(ggplot2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {},
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0gAAANICAMAAADKOT/pAAAAHlBMVEUAAAAAv8QzMzNNTU18\nrgDHfP/r6+vy8vL4dm3///8RzeqeAAAACXBIWXMAABJ0AAASdAHeZh94AAAgAElEQVR4nO3d\njVrjuLKFYdgJOcP93/Dp/OI4skolVUkq+1vP3gwwsIRkvSMn0PTXLyGkOV+jvwBC9hAgEWIQ\nIBFiECARYhAgEWIQIBFiECARYhAgEWKQZkiXbKR/3xz3fu8BvGcQfgK1M7DQoQiQhP7o+zD8\nBIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8BIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8\nBIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8BIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8\nBIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8BIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8\nBIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8BIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8\nBIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8BIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8\nBIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8BIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8\nBIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8BIDUsAiK5fLuj74Pw08ASA2LoFgu7/7o+zD8\nBHYE6fR65XRa/OMel0VQLJd3f/R9GH4C+4H0YnO6mzotaAGpfQTn+ugT2A2k1/lzP42e/3hG\nmM1/76lakuwA5o2r/uj7MPwEdgPppeYT0v/+Rfjc/7ZS/RUTMmGGQcIT2VOaIF0jnK/ftwia\nqo7uxwANn1vUH/3OKPwEjnBrd40wm+9XfDBNehU1IzjXR58AkG6z+V5Gus+rWa6qRVb0R9+H\n4SewN0inuqe/v1extjTpVdSM4FwffQI7hFTzDdk1JOt7vEmvomYE5/roE7jO4JyM9Gld4/0j\nQglI398XO0pAkuoDTiDtBkgpSrcYUAKSVB9hAoVyjgzp37/PUbpsaSq9BECS6iedQJ2dQ0O6\nJifpksZUdj2AJNXPNAEDPQeHdElZev/QOklAkurHT8CWz+EhXT4trT+6ghKQpPphE/DyA6Rr\n8pLWlPLltwHkD2kKkOQR3t8coib1ZfVN/z8ha0sJSFJ9rwm4udnPN2SzkWaTeGde0kUlCUhS\nvfsEzOXYzMBChyJjfmdDXtJFIQlIUr1jv5sdixlY6FBk0C8/ESRdiiUBSar36PfXswiQsotg\nJAlIUr1xfy8+fwFSfhEESZcySUCS6s36+/L5C5CkRbCQBCSp3qJ/CKBngCQugiDpUiAJSFJ9\nY/84QM8ASV6EdklAkurr+0cLegRIBYvQLAlIUr3vD2dfpr0EFjoUGf27v1slTXoVNSM41+v7\nlYfQpJfAQocioyFdBEiSpEmvomYE53pVf82N3KSXwEKHIsMhlUtKDyD2t+VIkNSEHiNYfakb\nAVLhIjRJmvQqakZwri/rrzN0H8HoK93sB5ISUo2kSa+iZgTnerm/3tB9BJsvdLsfSIWL0PIw\nadKrqBnBuV7obzJ0H8Hiy8z1A6l0ERqOpEmvomYE5/pMfzui2wgWX2auH0jFi1AvadKrqBnB\nuX6r38LQfYT2LzLfDyQzSNuSJr2KmhGc65P9VohuIzR/jUI/kMoXofph0qRXUTOCc/1nvyGi\n2witX6LUDyTFItTe3E16FTUjONev+o0VXaa9BBY6FIkCaUvSpFdRM4Jz/bLfHNFthMavUOwH\nkmYRiiWtBijtr8yOILkoukx7CSx0KDINpEvdkTTpVdSM4Fx/73dCdBvBoCPbD6QqSDpJk15F\nzQjO9dd+P0WXaS+BhQ5F5oFUdyRNehU1IzjX/7oqukx7CSx0KDIhJJWkSa+iZgTfemdFl2kv\ngYUORSaCJB5JqecbJr2KmhEcu/0VXaa9BBY6FJkRkuZImvQqakZwa+6h6DLtJbDQochMkGqO\npEmvomYEn9pOii7TXgILHYpMCUlxJE16FTUjeJR2U3SZ9hJY6FBkKkgVR9KkV1Ezgn1lT0bT\nXgILHYrMCan8SJr0KmpGsC7squgy7SWw0KHIXJD0R9KkV1Ezgm1db0bTXgILHYpEhfSUNOlV\n1IxgWfapKNgEUv1AqlgE7ZE06VXUjGBXlTqMQk0g3Q8kT0gPSZNeRc0IVkXpW7pAE9jqB1LN\nIiiPpEmvomYEm5qtR0ZhJrDdDyRXSHdJk15FzQgWJdtPMASZQK4fSC6QLkBaJ/c8XYgJ5PuB\nVLUIunu7Sa+iZoTWgvzT3QEmIPUDyQfSm6RJr6JmhLZPl75pNP0E5H4g1S2C6kia9CpqRmj5\nZInR9BMo6QeSE6SlpEmvomaEhs8VGc0+gaJ+IAGpZITqz5SPo8vcEyjsB1LlImieAZ/0KmpG\nqPy8IkYzT6C4H0hekC5AKrmru9fPOoHyfiABqWSEmk8qZTTtBDT9QKpdBIWkSa+iZgT9p5Qz\nmnQCun4gAalkBPVnKBjNOQFlP5D8If036VXUjKD8eM1xdJlxAup+IFUvQrmkSa+iZgTdh+sY\nTTgBfT+QgFQyguaDtYymm0BNP5AcIW38LS/mmWsf6h1NNoGqfiDVX0blH6bwy0z7sILRXBOo\n7AcSkEpGKP3AGkZTTaC2H0iekDpJmmYfVh1Hl4kmUN8PpIbLOMuRNMs+rGQ0zwQa+oEEpJIR\nSj6o2tEsE2jpB5IrpD6SptiH9YwmmUBbP5BaLiOQXmlxNMUEGvuBZABp9L3dBPuwhdEUE2jt\nB5IvpC6Shu/DpuPoMsEE2vuB1HQZgXRNI6PxEzDoB5IFpMH3doP3YbOj0ROw6AeSM6Qeksbu\nw3ZHQOoUIOUzch8aMAJSr8wLaY4HSQP3oYkjIHXK/JDGShq3D20cAalTgJTPsH1owwhIvQKk\nfEbtQytHQOqUiSEVPEjylzRmHxrd1t3qgdQlASANPZKG7ENDR0DqFCDlM2IfGjICUq8Eh+Qu\nacA+NHUEpE5phuSZB6Tchzwh9fqS/LNwNPpLIeWZ+USa4d6u+3/Qbc8jTqReiQ7pd2eQrB0B\nqVP2AslLUud9aM0ISL0yNaSC7yTtC5K9IyB1SghIGUm/lx1BcnAEpE7ZDSQnST33oYcjIHUK\nkPLptw/Nn2Z41AOpS+aGJD9I2g0kJ0dA6pQYkLYl/Tp/K6nXPvRyBKROAVI+nfahmyMgdQqQ\n8umzD/0cAalTJockPki69oeH5McISL2yI0guknrsQ09HQOqUIJA2Je0Kkks9kLoESPl02Ieu\njoDUKXuA5PkgyX8f+joCUqfMDkl6kLSE5CHJfR86OwJSpwApH+996PpEwzVA6pMokLYk3fvD\nQnJ3BKRO2RUkB0m++9DfEZA6ZXpIwr1daEgdHAGpU4CUj+c+7OEISJ0SBtKGpEd/REhdHAGp\nU/YFyV6S3z7s4whInQKkfNz2YSdHQOqU+SHlHyTFh+S8D4HUJzuB5PYgyWsf/p1HQJL6gdQf\nkrkkp324uK8DktQPJFNIaUkxIZ2BpOgHEpDSeXuiAUhSP5CMLmMJJK8HSR778P0JOyBJ/UAa\nAMlaksM+XD3xDSSpH0hASmX1DSQgSf1AsoWUlBQP0vobsUCS+oHUE5LTgyTzffjxAw1AkvqB\nNAKSsSTrffj5k0FAkvqBZHUZ9wMp8RN2QJL6gQSkdRI/qQokqR9IxpBSkhb9ASClfuIbSFI/\nkIZAspVkug+Tf3ICSFI/kID0lvQfQQKS1A8ks8u4M0jvI9gNkAqQ+mRHkFweJBnuw7QjIIn9\nQBoDyVSS3T7ccAQksR9I1pASksJASj9AugBJ7gcSkP6y5QhIYj+QOkPyeJBktQ83HQFJ7AeS\n3WVUQbKUZDSBbUdAEvuBBKRHNh8gXYAk9wMJSI9kHAFJ7AeSOaRPSe/9k0LKOQKS2A+kUZAM\nJVlMIOsISGI/kIB0C5Da+oFkeBkDQ8o7ApLYD6TukOwfJLVPQHAEJLEfSEAq+OtbgCT1A8ke\n0oekNCQ7SXaQNkdoHEAIkPoESPm0TkB0BCSxH0hAkh0BSewHkuVlLINk/iAJSPII3v1AGgfJ\nTFLbBAocAUnsB9LRIZU4ApLYDyQHSGtJQBICpD7ZGyTrB0ktEyhyBCSxH0gDIVlJsoAkjFA/\nQEmA1CdhIG08SJoYUtmBBCSxH0iHhlToCEhiP5BGQDJ+kNQOSRyhdoCyAKlPwkFaSZoWUumB\nBCSxH0gjIRlJqp1AsSMgif1AAhKQDPqBZHsZI0EqdwQksR9IQyDZPkhqhFQyQtUAxQFSn+wX\nko2kugkoDiQgif1AOiokjSMgif1AcoH0LmlqSGUjVAygCJD6ZIeQTB8k1UxAdSABSewH0sEh\nFY6gH0ATIPVJIEjJB0kZSCaSKiagO5CAJPYD6ZCQlI6AJPYD6dCQikfQDqALkPpkj5AsHySp\nJ6A9kIAk9gPJB9KbpBwkC0nVkMpHUA6gDJD6pADS6XR6/vP26uMfj7gsQjoRIKkPJCCJ/TuB\ndLr9L/3Wb1dIqXu7ySDpHQFJ7N8HpNPrxe2VEJAMHyQBSR7Bu3+XkH7vt3b3t/73L05fVzIP\nSPIHPiD5f0WrvBx1H5kMjg7S842HpNkhdZeEo8OmAtLbK1Pf2hnc26kmUHFjx62d3L/DW7s5\nIC0lbfSPhaQbQfcFaQOkPlFBOq3f8QukRaoOJCCJ/fuA9HzC+/n/5WOka1wWYSNqSO2SaiAp\nR9B9uDZA6pPib8guIC2/H9sVUuJB0kSQ6g4kIIn9e4GUj8sibGVmSJWOgCT2A2kcJLMHSXpI\n6hG0n6CsB1KXACmf4gnUHkhAEvuB5AVpIUmA1CxJDUk/gvozdPVA6hIg5VM6geoDCUhiP5AG\nQrK6t9NCqhhB/ymqeiB1SSxInw+SJoFUfyABSewH0gSQWiUpIdWMUPE5mnogdQmQ8imbQMOB\nBCSxH0iHg1Q1Qs0nKeqB1CUxIf1J2u7vCKnlQAKS2A+ko0GqG6Hqs8rrgdQlu4fUKKlkAk0H\nEpDEfiA5XMaZIVWOUPdpxfVA6hIg5VMwgbYDCUhiP5CGQrJ5kKSAVDtC5eeV1gOpS/YPqU2S\nPIHGAwlIYj+QDgWpeoTaTyysB1KXBIX0kjQcUuuBBCSxH0hjIZk8SCqGVD9C9WeW1QOpS4CU\njzSB5gMJSGI/kDwuox5Sk6RSSA0j1H9qUT2QugRI+QgTaD+QgCT2A+k4kFpGaPjcknogdcme\nIVk8SMpPwOBAApLYD6RJILVIKoPUMACQxH4gOUJ6ShoKyeJAApLYD6SjQGroB5LcD6TRkAwe\nJOUmYHIgAUnsB5LLZZwPUkP9BUhyP5BmgdQgKTMBmwMJSGI/kI4Bqb78PkLj50v1QOoSIOWz\nPQGjAwlIYj+QhkNqf5AkQ6rvfozQWiDUA6lLwkJ6SCqDVC8JSPII3v1A2jMkK0dAEvuB5HMZ\ngaSqB1KX7BxS84OkrQlYPdUAJLkfSAeAVF38N0J7RbYeSF1yEEjVkjYmYHcgAUnsB9L+IdXW\nLkcw6MjVA6lLgJRPegKGBxKQxH4guUK6SxL7XSHVtr6NYFGSqQdSlxwFUq0kIMkjePcDaa+Q\nLB0BSewHktNlBJKmHkhdsntIjQ+SUhMwdQQksR9IQCoZwaZmsx5IXXIYSJWSEhOwfO77AiS5\nH0j7hlRX+DmCUc9WPZC6JDCkm6SCfmNIxgcSkMR+IO0aUl1fYgSroo16IHXJcSDVSQKSPIJ3\nP5C8LuNISNaOgCT2AwlIJSOYNaXrgdQlB4DU9CBpPQHrpxqAJPcDaS5IVZK2INV0bYxgV5Ws\nB1KXACkfIMkjePcDaX+Q7B0BSewHkjOkq6SifiC5juDdD6TdQbJ/qgFIcj+Q3C5jHaQaSWlI\nFUXbI1iWJeqB1CVAygdI8gje/UDaGyQPR0AS+4E0C6SGB0lAkkfw7gfSbJAqJC0n4OIISGI/\nkIBUMoJt3Uc9kLokNKRvIMn1QOqSY0Cqf5C0mICPIyCJ/UACUskIxn3reiB1SUxIl0pIekl/\nE3ByBCSxH0hAKhnBunBVD6QuOQik6ns7IMkjePcDaU+QvBwBSewH0nyQ1JKAJI/g3Q+kHUFy\ncwQksR9I7pC+y68ikNxG8O4HEpBKRrCvfKsHUpccDZJW0mMCfo6AJPYDyfEyAqm0HkhdAqR8\ngCSP4N0PpJkgVT5Iuk/A0RGQxP5jQBqUByTFZzyPpJrRnpBqPpccIoc7kZRH0m0CngcSJ5LY\nf4wTyWUR5ACptB5IXRIc0reiH0hOI3j3A2lKSDpJ1wm4OgKS2A8kIJWM4FP7qgdSl0SFdAFS\nYT2QuuRAkKoeJP2bgK8jIIn9QJoTkkoSkOQRvPuBtAtIzo6AJPYDCUglI3gVP+qB1CVHglTz\nIAlI8gje/UDqAOlb8zkVR5K7IyCJ/UACUskIbs33eiB1CZDyAZI8gnc/kFwvYwWkigdJ7o6A\nJPYDaVZICklAkkfw7gcSkAoCJKkfSPEh+TsCktgPpOkgqR8kAalgBO/+npC+vtJvfElQjgmp\nVFIHR0AS+4EEpIIASeofdmt3KEiuD5KAVDKCdz+Q5oOke5DUwxGQxH4HSHcXt5df19zfuP9v\n9a7Fh77e2hWkDvd2QCoawbvfFdLXXcsS0uJdXwtTz7eBdA2Q7Efw7neGtBC1+a6/92/d5AFp\nO10cAUns93iM9Afj9/d1Ii2YLN71+lggLQMk8xG8+30hPe/X1g+bUpAy93ZHhVQiCUhlI3j3\nu0L6+jxvVu9anUhbhYeA9PMvj1fLIfVxBCSx3+Xp76/H3k88Alq9C0iP/DxzewtI5iN49/tA\n+voqgLR6suH15h4hSZJ+Frm+DSTrEbz7fb4h+zpu/p7+XoDJPP19UEg/PytJpUdSJ0dAEvtd\nIT3YrJ9seL7rCN+QLbq3+/lZSwKS9Qje/fwxivGQ1o4WkqR2IJWO4N0PpOGQPh39k1R2JD0d\nRd+HQOqTXUNKOfr5AZLxCN79QJoG0tsbP0CyHcG7H0iDIb1/A+migvRyFH0fhp8AkBoWoSwC\npJUjnSQglY/g3Q+kPpC2JK0dfUjKdQOpfATvfiANhfRxIGkgvZ77Dr8Pw08ASA2LUJYySKn3\n/QDJcATvfiB5X0Y1pEvpkQQkxQje/UAaCSnpqPRI+vuphvD7MPwEgNSwCIXRQ3qXtFkMJM0I\n3v1AGghpw9Hq5m6rGEiaEbz7HSD932ZqHRwWUu5IWvy8avh9GH4CQGpYhMLUQLoAyXoE734g\njYO07ejtSNroBZJqBO9+IHWClJCUgSQfScs/iRR+H4afAJAaFqEwbZA2jyQg6Ubw7geS+2Xc\ngpRz9HYkJf89kHQjePcDaX5IKUlvf8Y8/D4MPwEgNSxCaeogCUcSkJQjePePhFTsY5eQBEdA\nsh3Bu38YpK+1j4yWQ0LKS1o6ir8Pw0/ACdJ/6RhBEv+Cv8+4LEJpPCCdgaQcwbt/FKTb74n8\n+9WQr1/BmkSyR0iio+Uz4B//DkjaEbz7h55IX6+D6faPr68PNClIFbd6LotQmvQ3kmRIFxnS\n463w+zD8BMbf2r1B2kKyPpGyf5nS/iB9SjoDSTuCd/8MkB5/89jv5u//3jGk7GcCyW4E7/4Z\nID3fUXxrp47LIhQnBangQAKS5Qje/ZNA0tzaHQfSlqSVo/j7MPwExn8f6WvxtMPmXzX28WSD\n+Hf8Ack2QJL6x0H6ekL6+9v6nn93kgDpaz+QhM/8e7rh7d1A0o/g3R/wZ+3+npkoLnBZhOIk\nIJUdSBtH0tpR/H0YfgJAaliE4gBJrI8+ASA1LEJxGiD9SVq8D0gVI3j3x4UU5vtIqW8klTpK\nHklrR/H3YfgJhIT0q/5+7HSQig+k1NMNHwdS/H0YfgIxIf1qHQ2G9HlvVw4pcW8HpJoRvPtD\nQlLHZRHKYwLpJQlINSN49wNpIKSiT14fSZ+O4u/D8BNwgWSfD0jKOzsgtQZIUn9ESNGebPiA\npLiz+7y3A1LVCN79AW/tvhZ/rrYwLotQniZIqyMp4Sj+Pgw/gaCQ3v9ZEJdFKM8GpMLPBpLF\nCN79QOoH6SVJB+ldEpDqRvDuB1J/SLo7uzSktw8Ivw/DTyAkpNfvHiovcFmE8hhCSh1I8fdh\n+AmEg/T1liiQftsgLSUBqXIE734gDYNU/vlAah/Bu3/gn5AtpxD9JxvsICUdxd+H4Scw9nc2\nFD/OKfio0+n0euX26usd00FS39ktJAGpdgTv/rGQCs+ajZ9sWLz7dPvf47X3dwDJIECS+l0g\nndMxgpR4jHT6I3Rav2MPkF6fk3QUfx+Gn8BoSFW3dolPWkK63dMt3vG/fykaxTEPSPc3HiZU\nBStIDl8imS6FkOqfbMhDulqaFNJdEpBIUbxPpISk1Z3cadJbu9u9Xc2d3ere7qM/+p1R+AmM\nvrWzerIBSM4BktQfEFL2yYbTlJAuQMrXR59AUEifH/F4fuH0fIw02dPfCUjajqWkj/7o+zD8\nBIJ+HynxIY/vwj6ftZvsG7K2kD77o+/D8BMI97N225KycVkExXIBSaiPPoGQkAL+0OoCUt1D\npMtS0md/9H0YfgJAalgExXKZQkr0R9+H4ScQElJFXBZBsVyfkPQlL0iff8d5+H0YfgJAalgE\nxXIt/oxsNaTkX0zx6I++D8NPICSkwLd2QNqojz4BIDUsgmK5/n4jV/VDpMvlvCkp/D4MP4GQ\nkHLv3IjLIiiWywLSGUgtI3j3O0CyT9pMnN8iZA1pLSn8Pgw/gdiQAt/aVZSct4+k8Psw/AQi\n39odDNL5vH0khd+H4ScQEtIenmzQdwCpbQTv/siQygtcFkGxXC9ILQ+Rtu/twu/D8BMICaki\nLougWC5zSO+Swu/D8BMAUsMiKJZrDUlfcQZS2wje/cEgBf2VxVaQNiSF34fhJwCkhkVQLNfl\nKQlI6froEwgHafG+WE82NEI6ryC9Swq/D8NPYCSktl+ir3I0DaTq5xqekD7+auZ7f/R9GH4C\nI39nw9dv8fdUPz5KxwhIzQGS1D8OUppIMonfa6dQ9DsdJH3D09ElKSn8Pgw/ASdI3+kkIJXl\n4xuyis+9xWURFMt1aYR0BlLrCN794yApOCw/soLRTiH9SQq/D8NPYCCkyicbqhztCVLqSAq/\nD8NPYCSkO4siB3v5PlLrcw231xNHUvh9GH4Cwx8jAakkZyA1j+DdHwxSZVwWQbFc1xdGd3Yp\nSeH3YfgJDH2y4bf++0jauCyCYrmuL4CUqY8+gbg/2aCLyyIoluv6oh7Su6MlpIek8Psw/AQC\n/6ydKi6LoFiu6wszSJ9HUvh9GH4CQGpYBMVyXV8YPddwSRxJ4fdh+AkAqWERFMt1fdEM6fX2\nz/pICr8Pw08ASA2LoFiu6wuzO7vPIyn8Pgw/ASA1LIJiua4vgJSpjz4BIDUsgmK5ri8MIa0l\nhd+H4ScApIZFUCzX7WUjpOW7lpD+28E+DD8BIDUsgmK5bi8rn2tIHEhA0o/g3e8AyT5AWkFa\nSQq/D8NPAEgNi6BYrttLIG3XR58At3YNi6BYrttLu4dIK0j/hd+H4ScApIZFUCzX7aXhgbSS\nFH4fhp8AkBoWQbFct5d+kP6Lvg+BlAqQPga4vXxs/W/d5xZB+vz7mY0DJKkfSL0g/TRB+nj/\nuySTLzMTIEn9QJob0saB1PtIApLUD6TOkL6NIfWRBCSpH0ghICX+DZBUI3j3D4NU/5tWK+Ky\nCIrlur5whuQrCUhSP5CmhrR5Z3fpKwlIUj+QOkH6AVKmHkiJ/N///aQDpFpIyX+3huQpCUhS\nP5BmhpQ7kLpKApLUD6S+kK6/KLP8E4FkNoJ3P5DiQ7r4SwKS1A+kAJC2/vVzYX/dJQFJ6gfS\nxJDyB1IKkpckIEn9QOoD6bkGLpD8jyQgSf18HykwpL8jyfthEpCkfn7Wriuk7xpI2x/wB8lZ\nEpCkfiDNC0k8kN5+oYqrJCBJ/UDaHSQPSUCS+oHUH1KpJB0kV0lAkvqB1AXSTwOk7Md0kgQk\nqR9I00IqOJBWv3USSNsjePcDaR+QnCUBSeoHUk9I3+aQVr8I2UsSkKR+B0j2OTQk4aM2IBlL\nApLUD6RZIZUdSJ0kAUnqB9K+IDlJApLUD6QekH76QVpIMqQEJKkfSJNDkj+whyQgSf1AGgCp\nRFLpgZT4C8wcJAFJ6gdSR0gXF0g9JAFJ6gdST0iXPpDsJQFJ6gfS3JBK6j8hLSWZUAKS1A+k\nOSEpDqT0X/NsKwlIUj+QOlzGnwGQbCUBSeoH0o4gbUpqpgQkqR9IU0Mq6/9NQrKUBCSpH0j9\nIF3KIakOpC1Ib5LaKAFJ6gfSCEiiJCWkIkktlIAk9QNpz5DMDiUgSf3HgDQ4j13+77UnJOkz\nHo60I9zGeM/7oaT6ssneEvxEWh4XZSeS7kC6TmDzSLK5v+NEkvqPcSK5LEJxBkNaSaqiBCSp\nH0j7gJSV1E4JSFI/kLpBur2uglQ6gAhpLUltCUhSP5AmhKQ8kAogtVICktQPpJ1AEiV9UlJY\nApLUD6QhkPKSfCAlJBVTApLUDyT3y/hTC6l4hPsEZEkpSmWWgCT1A2k+SNoDqRxSmlIBJiBJ\n/UDqBenxlj+kvKQNSoIlIEn9QNoNpLIj6ZoNSxlMQJL6gTQtpPIRtJC2KW1pApLUD6T9QNJI\nylr61AQkqR9I00FS39nVQRIovXMCktQPJO/L+NMfUqmkAksPT0CS+oE0BlJGktrR3wTUkK4p\nwpS637MLkPpkD5Ceb84H6aKw5IQKSH1yKEj6O7vFBGolXWow2ckCUp/sCZL4IGkUpFsaNdVK\nA1KfHBGSZohPSNWSrnHUdLxIV65rgJTPYgIWkO4ZvQN3EunKdU1kSB87W4BUcWfnA+mR0Rsx\neqQr1zVAymc5AXNJ9xEAVRnpynXNDiD9vSMopGVGb85Ika5c1+wR0pakCkfvE3CRJKzQ6N06\ncYR17ZsDQao5kMZDEjN6O4+LsK59A6R8kpBMJfGzdlI/3wtn4AYAABK+SURBVEfqDSn/IKkd\nksuRBCSpH0i+lzGxq0sg6QYBkjyCdz+QpoJUdSCtJ+AgCUhSP5C6QFq+C0if9UDqEiDls56A\nvSQgSf1AmhGSchAgySN49wNpFKSUpLoDaROSnSQgSf1Acr2MqR3tD8n+SAKS1A+k7pAy93bW\nkMwkAUnqB9KEkLSjfE4ASOsRvPuB1APS+zuB9FEPpC45CqTKO7vUBIwlAUnqBxKQSkYw6tmq\nB1KXHAySepTUBGwlAUnqB5LnZUzv5s3nv4HkN4J3P5DmgVR7Z5eFZCMJSFI/kHYKyfZIApLU\nD6QOkNb9XSGZSAKS1A+k6SDph0lPAEjLEbz7gTQNpOoDSYBkIQlIUj+QdgvJ8kgCktQPJMfL\n+AOkwnogdcmhIFWMszUBO0lAkvqBNA7Su6T6AwlIBSN49wPJH9LHVewFyU4SkKR+IA2AlLy3\nA5LrCN79QJoMUs042xOwkgQkqR9IQCoZobVAqAdSlwSF9KOB1HBnVwCpVRKQpH4g7RqS1ZEE\nJKkfSO6QEldxG1LVQAWQGiUBSeoH0kBIS0lOkIyOJCBJ/UCaA1LLnV0RpDZJQJL6gTQC0ue9\nnRskmyMJSFI/kKaCVDdQEaQmSUCS+oHkdhl/5oBkciQBSeoH0hSQmu7sCiG1SAKS1A8kb0iX\n0ZAsjiQgSf1AOg6kBklAkvqBNBLSS1KTI3ECQALSLbuH1HYgFUOqlwQkqR9IXpfxJwfp0hVS\n+5EEJKkfSEeCVC0JSFI/kCaCVDuUPAEgefcDyRnSbYDPf/sGqfFAUkCqlQQkqR9Ih4DUeiQB\nSeoH0rEgVUoCktQPpKGQ7pIaHRVNAEi+/UByuow/k0KqkwQkqR9IYyAt7+1a7+zKJtAkCUhS\nP5B8Id0HSPx7IC3rgdQlx4BUP1bZBFokAUnqB9JhILU83wAkqR9IwyE139lpIVVIApLUDySf\ny/gzIaSGIwlIUj+QxkL6HgFJLwlIUv9eIJ1Op7fXTqe/98SA1DBY8QSA5Na/E0in2/+er/29\n9YjLImTztmOT/U9I7QeSHpJaEpCk/n1AOr1e3F4CSfiygGTevzdIpyek153d//7F8WvbyGPD\nZj7iAen3Cannl5X9wshuozmRfl+3dg9JQyAV7NcVpM5fGJKOGB2kv6cdXv/e5VjO5f0WKndr\nZ3Bnp5lA5aMkbu2k/r3d2i35hID03RdSpSQgSf27g3Rav+Mal0XI5X2zzgSp7vkGIEn9+4D0\nfJ7u9Lu4rRv5faQCSO/3dk2jVUFSSQKS1L8TSM9vw96/Efv2rN01LouQy8SQqo4kIEn9e4GU\nj8siZPITAZJGEpCkfiAdDlKNJCBJ/UCaAlLbcECSR/DuB5IjpNcAyQ9aPm3XNpx2AnpJQJL6\ngXRASPrnG4Ak9QPpyJCKJQFJ6gfSMEiLB0mNw+knoJUEJKkfSA6Xcb1LgSTVA6lLgJRPxQSU\nkoAk9QPJD9LfAOkPGwlJ+XwDkKR+IA2H9N3qqAlSmSQgSf1AOigknSQgSf1AAlLJCBUDKAKk\nPgkG6WOHbvUPhaSSBCSpH0jjITUPWDkBIBn2A8kN0mKA9MedJ4FUIAlIUj+QjgtJIQlIUj+Q\nxkNqllQ9ASCZ9QPJ/DJ+7s7pIYmSgCT1A2kgJKt7u/oJlEoCktQPpGGQzkB61QOpS0JCWg6Q\n/LgZIJVKApLUD6RjQyqUBCSpH0gTQGqVBCR5BO9+IFlfxsTGTPaf54BUJglIUj+Qjg6pSBKQ\npH4gjYRk9CAJSPII3v1AcoL0NkDiw87TQCqRBCSpH0hAKpAEJKkfSED6g7QpCUhSP5CML2Nq\nS4qQGiU1TwBIzf1AGgPpPBUkURKQpH4g+UB6H+Dzw85nw3s7gwkIkoAk9QMJSNcAqbEfSEMh\n2fzaBosJ5CUBSeoHku1lTG7Hz/7zdJDykoAk9QNpCkhtkmwhpSQBSeoHEpAeAVJLP5BcIK0G\n+PiwpyOb321nNIGMJCBJ/UAaAek8JaSMJCBJ/UAyvYzpnRgO0ockIEn9QALSX4BU3Q+kgZCu\nr1o822A3gS1JQJL6geQBaT3A6u3ztJC2JAFJ6gfSYEgW93YekMSfzbAMkPokDKTCG6OJIW0c\nSUCS+oE0DtL9jckgpSUBSeoHUn9I5xSkFkk+v71F/s18ZgFSn0SD9DHA+5tRIMl/DYBVgNQn\nQMrHeB8mJAFJ6geSw/czPwZ4f/PNkcWDJOt9CCR9P5C6QzpPD+lTEpCkfiABKZG1JCBJ/UCy\nh/Q5wNtbK0cGD5IcIT0mAySpH0i9Ia0PpCkhrSUBSeoHkv2fQvgcYPnGB6T2ezuPffguCUhS\nP5CAlA6QVP1AGgRp8Z45Ib1JApLUDyRzSIkBFq9/HkjtD5Kc9uFCEpCkfiABaSuLh0lAkvqB\nZHYZN+/sokJaSAKS1A+kvpASjtofJLntw/W3k7wCpD4BUj5++7CTJCD1SShIqQH+Xs1Bqpbk\nuA/7SAJSn4SAlNlwn5De//3EkPpIAlKf7AdS8kCaGlIXSUDqk71Dan2Q5LsPgVTSDyRjSMkB\nXq+lHc0NqYMkIPXJbiBtHEit93be+9BdEpD6JAKk3F6LDsldEpD6BEj59IPkJAlIfRIIUnqA\n5ysbjlofJPnvQ2dJQOqTvUDaOpBaj6QO+9BXEpD6JACk7DbbASTfMwlIfbIzSImPmB+SqyQg\n9clOIG0fSI0Pkvrswx8/SkDqkziQNga4/yM4JEdJQOqT+SHlN1g5pDpJvfahmyQg9ck+IOUc\nxYDkJglIfQKkfPrtQydJQOqTMJC2Bri9zEJqepDUcR/6SAJSn0wPSdhbt/68o6Yjqec+dJEE\npD5phuSd59bKfcwT0sa/fkJy+OqM8ydp9FdClOFEyqfzf9DtzyROpD6JAmlzgIvoqOlBUu99\naH57B6Q+mR2StKs0kGokdd+H1pKA1CdAyqf/PjSWBKQ+2QEk0VEsSMaSgNQnk0MSd1QRpIYH\nSSP2oakkIPXJoSBVSBqyD38MKQGpT+JDKnAUDZKlJCD1SQxImQF2Ccnw9g5IfTI3JHkvlUGq\nf5A0bB9aSQJSn4SHVOSo/kgatw+Nbu+A1CdAymfgPrSRBKQ+mRpSwTbaMSQbSUDqk+iQyhzV\nSxq7Dw0oAalPgJTP4H3YLglIfRIBUu5Ddg7pTVIVpdETMOgHUuNlNDyQqp8AH74PWyUNn0B7\nP5AmhKSVNME+bJM0wQRa+4HkDul8BEhth9IME2jsB1LbZTQ9kCJDapI0xQTa+oE0E6RKSZPs\nw3pKk0ygpR9I3pA0jmJDqpc0ywQa+oHUdBmND6TgkKopzTOB6n4gTQWpTtJE+7BO0kQTqO0H\nkjMknaPwkN4llVKaaQKV/UBquYzmB1J8SFWU5ppAVT+QfCEpHdVJmmwf6iVNNoGafiA1XMaC\nzXJESHpK001A3w8kA0jbH6H4qYZHdgFpJUmkNN8E1P1A6gJJUVohacZ9qJI04wSU/UCqv4z2\nTzVcsxNIqkNpygno+oHkCanC0QuSQtKk+7Cc0qQT0PQDaTZIFUfStPuwlNK0EyjvB1L1Ing8\n1XDNjiCtJW3+TWzTTqC4H0gdICn71ZJm3odFlGaeQGE/kGoXQb5hOQPpmgJKc0+gqB9IrZC2\nP+J1Y1cJqVjS7PtQpDT7BAr6geQG6VwLSX0kTb8P15LWizb9BOR+IFUuguYpu1pIpZIC7MM8\npQATkPqB5AXpXA9JeySF2Ic5SiEmkO8HUt0iaA4kIN2zTSnIBHL9QHKCdG6BpJQUZh9uUQoz\nge1+IFUtwvYTUM8svxcLpFc+KN1WMNAEtvqB1ARp8wPObZB0kkLtwxSlUBNI9wOpZhF0BxKQ\n3vNJKdgEUv1AaoG0+QFnK0hFksLtw09KrX8vuhQgXTMbJOWBVHUVNZLCQepPCUjXzApp8wPO\n7ZA0N3cBISUpOVoC0jWTQdIeSG2QCiSFhJSm5GUJSNdMCmnzA9Z/DqlukXcP6dLxWALSNXNB\nEq/4x5/na4MkS4oLqduxBKRrpoKkvrGrvorFkiJD2qBkbQlI18SC9PkHzIEkxN8SkK6ZCVKF\no+qrWCopPKTf3zQlO0tAumZGSJsfYAipVNIOIF22jiUjS0C6ZiJINQcSkMT6R/+GJQNMQLpm\nHkhVjhquYpmk3UDaptRqCUjXHBdSmaQdQbpkLLVgAtI100Cqc2QCKSdpX5AuHpaAdM0skCod\nNV3FEkm7g3TJWarCBKRrJoEkXsoNR21XsUDSHiFdspbUmIB0zRyQyh2ZQiqQtFNIl7wlFSYg\nXTMZpK0P2HLUehVFSfuFdE0WU/Hfm27xZeb6gVS6CPWOmq+iJGnfkC6SpRJMQLpmBkgNjtqv\noiBp95CuETAJmoB0zQSQxAuWcWRwFb+zlA4B6SJbyv1VZlZf6lY/kIoWodyRD6S8pKNAukbG\nlLxEQLpmOCTJ0TnryOQqfmcoHQnSNQWY1lcKSNcMhiTeOgiObK7i9zalo0G6pgjT4reLW365\niQBJXoRmR0ZX8XuT0hEh3VKm6fa7XG2/4I8ASVwEydFZdGR1Fb+3KB0W0jWlmHx/cR6QhEUQ\nL0SBI7Nt+L1B6dCQbinX5AQKSNlFENf/XOLIcBuuKd0tAekWlSZjUEDKLIK87EWMTLfhh6Sr\nJSD9RavJCBSQthfB6Di6GG/DBKVvIK1SwalNFJA2FkFe4nOxI+ttmKKU+O6SZaJBuqeOU5Up\nICUXQVzW81nhyH4bblDy0xQT0j0NnMpNAelzEeSV1DFy2YbblFw8RYZ0H6HZU14VkFaLIK/d\nWpHsyGcbCpRsRe0A0iM2nj72BpD+FqHkP0GfimRGftuwzJKFq/1AesYW1NZ+Eb+svnGHVLIs\nCUOFjjy3od5SlbT9QXrGC1QZLAsdinhDkhZhw1AhI+9taG3poAFSQfKzySzBNqFyRl1+YnL0\nNtxVjgzpdDq9v/b3jt8qSHlCCkWXfj96PHoD7jBHg3S6/W/x2t87rhFmo1ekYdT9Z/hHb77d\nZv+QTusXf6/fIsxGiUin6DLwD8OM3nn7zsEg/e9fhE/XIGqZx+iM3nm7TwLS6Gv+lh6QCgyF\nRlSS0ftwjxl9Td/ifWvncDf3PkDTZxf089Pf4gje/fv4yQZPSFULtBrAoCPbH30fhp8AkG6z\ncdGzHMC27rM/+j4MP4G9QHp/0lv99Lf7Knv3R9+H4SewG0iP77+eXq9pviELJHkE5/roE9gP\npHxcFkGxXN790fdh+AkAqWERFMvl3R99H4afAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4af\nAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4afAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4af\nAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4afAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4af\nAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4afAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4af\nAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4afAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4af\nAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4afAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4af\nAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4afAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4af\nAJAaFkGxXN790fdh+AkAqWERFMvl3R99H4afwFEg5SP+HbOzhwkMT4wZACkfJjA8MWYApHyY\nwPDEmAGQ8mECwxNjBs6QCDlGgESIQYBEiEGARIhBgESIQXwhnU4n137v3L76uJN4fOVRJ3D/\nqoNMwhXS6bkYQXO7dnEncbrPIOoEHoCCTMIT0un1ImZOr8sXcxL3/wxEncDpBSnEJIC0nR1A\nOkXYg1u5GQJS6D14zeKGIugkTqfQE3h+0SH+awCkrZzCQ4ryH/OtPL7oGP81ANJWTrdEnsQ+\nIJ3eX8waIOUS4r+Fm9kFpLWmWcPT37kEf/r7FOOZ4608Vv/5+tyT4BuyuYT4VuB2gnwvcyv3\nZ+1O9y9/9knwI0KEGARIhBgESIQYBEiEGARIhBgESIQYBEiEGARIhBgESIQYBEjD8vXK9fXR\nXw1pCxdwWIC0p3ABhwZAewkXcmiAtJdwIYfmCel5e3d7+/HO+00fiREu1dCsIN3sPAA9Hz6R\nEOFKDc07pL8Xz3+DpCjhQg3N+tbu7wWQYoULNTQ5SNzbRQoXamikE4lECZdraIC0l3C5hiZ7\na/cLpzjhQg3NNqRfHiOFChdqaDKQ+IZsqHCpCDEIkAgxCJAIMQiQCDEIkAgxCJAIMQiQCDEI\nkAgxCJAIMQiQCDEIkAgxCJAIMQiQCDEIkAgxCJAIMcj/AzkNBQ3QJQCdAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "plot without title"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "ggplot(sir.cn.out.long,aes(x=time,y=value,colour=variable,group=variable))+\n",
    "  # Add line\n",
    "  geom_line(lwd=2)+\n",
    "  #Add labels\n",
    "  xlab(\"Time\")+ylab(\"Number\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "R",
   "language": "R",
   "name": "ir"
  },
  "language_info": {
   "codemirror_mode": "r",
   "file_extension": ".r",
   "mimetype": "text/x-r-source",
   "name": "R",
   "pygments_lexer": "r",
   "version": "3.4.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}