From 662c083cf4561da98225fa02d9c318aec6d53390 Mon Sep 17 00:00:00 2001
From: Mark Penney <mark.penney@uwaterloo.ca>
Date: Thu, 4 Feb 2021 17:28:24 -0500
Subject: [PATCH]

---
 network-data/Prem/Canadian-Avgs.ipynb | 204 ++++++++++++++++++++++++++
 1 file changed, 204 insertions(+)
 create mode 100644 network-data/Prem/Canadian-Avgs.ipynb

diff --git a/network-data/Prem/Canadian-Avgs.ipynb b/network-data/Prem/Canadian-Avgs.ipynb
new file mode 100644
index 0000000..5492fe0
--- /dev/null
+++ b/network-data/Prem/Canadian-Avgs.ipynb
@@ -0,0 +1,204 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "%matplotlib notebook\n",
+    "import numpy as np\n",
+    "import matplotlib.pyplot as plt\n",
+    "import pandas as pd\n",
+    "from scipy import optimize\n",
+    "import math\n",
+    "import random\n",
+    "import csv\n",
+    "import itertools\n",
+    "from itertools import product\n",
+    "def Flatten(A):\n",
+    "    return list(itertools.chain.from_iterable(A))"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Average mixing numbers for Canadian population\n",
+    "\n",
+    "To adjust the POLYMOD results to fit Canadian demographic structures, we take the distributions extracted from POLYMOD and rescale them so that their averages match those determined by Prem."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Population statistics from Statcan 2019\n",
+    "N = [0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]\n",
+    "N[0] = 1943175\n",
+    "N[1] = 2039352\n",
+    "N[2] = 2031762\n",
+    "N[3] = 2114635\n",
+    "N[4] = 2476698\n",
+    "N[5] = 2625474\n",
+    "N[6] = 2603938\n",
+    "N[7] = 2580021\n",
+    "N[8] = 2421009\n",
+    "N[9] = 2396406\n",
+    "N[10] = 2502667\n",
+    "N[11] = 2749626\n",
+    "N[12] = 2511888\n",
+    "N[13] = 2096607\n",
+    "N[14] = 1706760\n",
+    "N[15] = 1164277+786704+510828+242554+74086+10795 #75+\n",
+    "Ntot = np.sum(N)\n",
+    "\n",
+    "# Load the Prem files\n",
+    "Malldf = pd.read_csv(\"Prem_Canada.csv\")\n",
+    "Mall = Malldf.values.tolist()\n",
+    "Mhomedf = pd.read_csv(\"Prem_home.csv\")\n",
+    "Mhome = Mhomedf.values.tolist()\n",
+    "Mschooldf = pd.read_csv(\"Prem_school.csv\")\n",
+    "Mschool = Mschooldf.values.tolist()\n",
+    "Mworkdf = pd.read_csv(\"Prem_work.csv\")\n",
+    "Mwork = Mworkdf.values.tolist()\n",
+    "Mrestdf = pd.read_csv(\"Prem_rest.csv\")\n",
+    "Mrest = Mrestdf.values.tolist()\n",
+    "\n",
+    "# Age brackets\n",
+    "Brac = ('Y','M','O')\n",
+    "Bracind = {'Y': (0,1,2,3,4), 'M': (5,6,7,8,9,10,11,12), 'O':(13,14,15)}\n",
+    "\n",
+    "# Determine contact matrices\n",
+    "allCont = {}\n",
+    "homeCont = {}\n",
+    "wsCont = {}\n",
+    "restCont = {}\n",
+    "for x, y in list(product(Brac,Brac)):\n",
+    "    tar = Bracind[y]\n",
+    "    src = Bracind[x]\n",
+    "    allSum = np.sum(Flatten([[N[j]*Mall[i][j] for i in tar] for j in src]))\n",
+    "    homeSum = np.sum(Flatten([[N[j]*Mhome[i][j] for i in tar] for j in src]))\n",
+    "    sSum = np.sum(Flatten([[N[j]*Mschool[i][j] for i in tar] for j in src]))\n",
+    "    wSum = np.sum(Flatten([[N[j]*Mwork[i][j] for i in tar] for j in src]))\n",
+    "    rSum = np.sum(Flatten([[N[j]*Mrest[i][j] for i in tar] for j in src]))\n",
+    "    Num = np.sum([N[i] for i in src])\n",
+    "    allCont[(x,y)] = allSum/Num\n",
+    "    homeCont[(x,y)] = homeSum/Num\n",
+    "    wsCont[(x,y)] = (wSum+sSum)/Num\n",
+    "    restCont[(x,y)] = rSum/Num"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Saving a copy of the mixing matrices\n",
+    "swap = ['Y', 'M', 'O']\n",
+    "\n",
+    "wsMat = np.full([3,3], 0.)\n",
+    "restMat = np.full([3,3], 0.)\n",
+    "for i, j in list(product(range(3), range(3))):\n",
+    "    wsMat[i][j] = wsCont[(swap[i], swap[j])]\n",
+    "    restMat[i][j] = restCont[(swap[i], swap[j])]\n",
+    "\n",
+    "pd.DataFrame(wsMat).to_csv(\"wsMixing.csv\", index=False)\n",
+    "pd.DataFrame(restMat).to_csv(\"restMixing.csv\", index=False)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 7,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "[[2.22220876 1.27618792 0.15504874]\n",
+      " [0.84971392 3.72295479 0.74083542]\n",
+      " [0.21683524 0.782842   0.66709872]]\n"
+     ]
+    }
+   ],
+   "source": [
+    "print(restMat)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Rescaling workschool distributions from POLYMOD to Prem\n",
+    "\n",
+    "Since the workschool distributions from POLYMOD are fit to a 2-variable distribution we have to make a choice as to how to adjust those parameters. We first modify only p, keeping alpha from POLYMOD fixed."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Input the POLYMOD parameters\n",
+    "POLY = {('Y', 'Y'): [0.4378787156156721, 0.093112367309826], ('Y', 'M'): [0.4144210645686216, 0.2988493341114929], ('Y', 'O'): [0.9071801708997953, 0.8256463657373431], ('M', 'Y'): [0.76732665579341, 0.41260047034256747], ('M', 'M'): [0.421073517853328, 0.12039762005212719], ('M', 'O'): [0.8424225575069703, 0.5573421127782885], ('O', 'Y'): [0.8976481193486225, 0.40369064730244425], ('O', 'M'): [0.7855029449186991, 0.33072950060358974], ('O', 'O'): [0.9407010943356611, 0.4756835100017134]}\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "Premparams = {}\n",
+    "for x, y in list(product(Brac,Brac)):\n",
+    "    paramval = (1-POLY[(x,y)][0])/wsCont[(x,y)]\n",
+    "    if paramval <1:\n",
+    "        Premparams[(x,y)] = ['geo', POLY[(x,y)][0], paramval ]\n",
+    "    else:\n",
+    "        Premparams[(x,y)] = ['pois', POLY[(x,y)][0], 1/paramval ]"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "print(Premparams)"
+   ]
+  },
+  {
+   "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.3"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
-- 
GitLab