{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Thalamo-cortical neuron\n",
"# Analysis of Ca2+ and voltage-dependent kinetics of the hyperpolarization-activated current"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Imports"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"import matplotlib.cm as cm\n",
"from matplotlib import ticker\n",
"from matplotlib.colors import LogNorm\n",
"\n",
"from PySONIC.neurons import getPointNeuron"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Functions"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"def plotIhKinetics(Vm, CCa, gatings, ylabel, cmap='viridis', fs=18, lw=2):\n",
" \n",
" mymap = cm.get_cmap(cmap)\n",
" sm = plt.cm.ScalarMappable(cmap=mymap, norm=LogNorm(CCa.min(), CCa.max()))\n",
" sm._A = []\n",
" fig, ax = plt.subplots(figsize=(5, 3))\n",
" for key in ['top', 'right']:\n",
" ax.spines[key].set_visible(False)\n",
" ax.set_xlabel('$V_m$ (mV)', fontsize=fs)\n",
" ax.set_ylabel(ylabel, fontsize=fs)\n",
" for c, gating, in zip(CCa, gatings):\n",
" ax.plot(Vm, gating, linewidth=lw, c=sm.to_rgba(c))\n",
" ax.xaxis.set_major_locator(ticker.MaxNLocator(2))\n",
" ax.yaxis.set_major_locator(ticker.MaxNLocator(2))\n",
" fig.subplots_adjust(right=0.85)\n",
" cbar_ax = fig.add_axes([0.87, 0.1, 0.02, 0.8])\n",
" fig.add_axes()\n",
" fig.colorbar(sm, cax=cbar_ax)\n",
" for item in ax.get_xticklabels() + ax.get_yticklabels() + cbar_ax.get_yticklabels():\n",
" item.set_fontsize(fs)\n",
" cbar_ax.set_ylabel('$[Ca^{2+}_i]\\ (uM)$', fontsize=fs)\n",
" return fig"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Parameters"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"pneuron = getPointNeuron('TC')\n",
"Vm = np.linspace(-100, 50, 100) # mV\n",
"CCa = np.logspace(np.log10(0.01), np.log10(10.0), 10) # uM\n",
"\n",
"# rate constants\n",
"alpha = pneuron.alphao(Vm)\n",
"beta = pneuron.betao(Vm)\n",
"\n",
"# proportion of regulating factor in unbound state (-)\n",
"P0 = pneuron.k2 / (pneuron.k2 + pneuron.k1 * (CCa * 1e-6)**4)\n",
"\n",
"# Extend to match dimensions (nCa, nV)\n",
"alpha = np.tile(alpha, (CCa.size, 1))\n",
"beta = np.tile(beta, (CCa.size, 1))\n",
"P0 = np.tile(P0, (Vm.size, 1)).T"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Open form"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAZ8AAADzCAYAAABZoxsBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3XecVOW5wPHfM7O9sgtL7wiCIgqoWJEoir3FGk001xKvsZeYm6hRE29ijzdeTSA3sWuMvYXYewUUxIIoIIjSWerWOc/945zFdZm2087OzPP1cz6ze857zvvOsu4zbxdVxRhjjMmkgN8FMMYYk38s+BhjjMk4Cz7GGGMyzoKPMcaYjLPgY4wxJuMs+BhjjMk4Cz7GGJOlRORYEXnd73IkwoKPMcZkGREJisjFwL2A+F2eRFjwMcaY7PMH4AjvNStZ8DHGmOxzs6pOBBb6XZBEWfAxxpgso6rf+l2GZFnwMcYYk3EWfIwxxmScBR9jjDEZZ8HHGGMyTESmicgrEa4NEZFHRWSNd9wtInUZLmLaie3nY4wxmSMipwF/BV5V1UkdrnUHZgJFwK1AAXApsAjYVVWbM1rYNCrwuwDGGNOVTflBua5eE4p4feacpo+Bxnanpqrq1I7pRCQI/Bq4Kkp2FwH9gR1U9VPvvneB54FTgGmdLX9XZcHHGGOiWL0mxHv/HhjxerDP/EZV3TnaM0SkBHgXGAPcDewXIekJwCttgQdAVV8QkXnetZwJPtbnY4wxUShKi7ZGPOJUAlQBx6vqKcBWN4pIDTAUt9mto1lA1ACXbazmY4wxUSjgELVvvFpEpgJPqepTEdKsB4arRo1W/bzXpWGufQtUiUi1qq6LVeZsYMHHGGOicGs+kft8gHWqembUZ6g6gBMjq0rvdXOYaw3eazlgwSdbHXjggTp9+nS/i2GM6Tqirgwdo+aTKm3dINEyixXAskZe9vmsWrXK7yIYY7KEAi04EQ+8ZjcROSzJrDZ4r6VhrpV2SJP18rLmY4wx8VIgFH0+ZMxmtzgt9l77hLnWF6hX1U0pyKdLsOBjjDFRKEpLBprdVLVeRBYC48JcHgvMSHshMigvm92MMSZuCqEoB6lrdgN4BJgsIiPbTojIZGBb4MEUPL/LsJpPDM3NLRQVFfpdDGOMTxShJfp4hFQ1uwFcD/wEeFFEbsKdH/QL3Lk/96Yojy7Bgk8Me51+KwClTcqAijJOPGxXDjlyF59LZYzJFAWcDC2BqaorRWQicAtwDe6w68eBS1W1KTOlyAwLPlEsWbySlrIAIsKGcviEZq6Y/jpXPfEalZsdDtt5BOdddCjBYNDvohpj0kSB5hT3UKjq4CjX5gEHpzTDLsj6fKLoP6AHNdtWQe8gVEKoyCFUCK0lAeprC7hnwQJ2O+2PHPOTW/l6yUq/i2uMSRNHJeJBavt88obVfKJoag2xurGBxkDInXtc6cbqQoTAegeaQQmwAIcjrr6Hfg3Cnf9zOrXdK6M/2BiTNRyEZqK2bqSyzydvWPCJoqSwgBm/OJsvVq5m7rcr+PDrb3n3qyUsWbsOqtwOyKADgQ1KQYPwdSFMuXgq+/buxXV/ONnn0htjUsWr4ZgUsuATQ2EwyKjePRnVuyfHjh0NwOI19fz7s/k8+/HnfLJsBaFqobUKgpuUokCA5+tXMPOEG3n0z/9JVbdyn9+BMSYZitCs1q+batbnk4CBtd04Y49deOyMk/jnaSdy6PbbEggIrRVCQx20VAira4Lsf+4dvPnax34X1xiTBHdV60DEwyTGfnJJGtO3NzcdfTBPnfUT9tlmCBqE5mpoqlGaq4Kc//fpvPj8h34X0xiTIFW35hPpwAYcJMSCT4oM61HL1BOP5PbjDqemrJRQqdDQQ2mpDHDZAy/yzFM5tTKGMXnFQSIeeAMOouzlY8Kw4JNi+207jCfOOJldB/VHC4TG7m4A+s1jr/LR7EV+F88Y00lun09BxMMkxoJPGvSqquBvJx3NYaNHQkBorFWaqwKcccPDbN7YEPsBxpguw/p80sN+cmlSGAxy/ZEHcvy4HSAgNNUoDbVBDjv9f/0umjGmE9pGu0Xp8zEJsOCTRgERrj54P47ecTs3ANUqq7sHOe+iO/0umjGmExwNRDywAQcJseCTZuIFoPED+qJBoakW3li/mnmfLvG7aMaYOCgQIhDxwAYcJMSCTwYUFRRw27GH0beqEqcImmsDnHFNTm3NYUzOUoQWDUY8TGIs+GRIbXkZN//QXai2pUJZ1yPIlVf9w+dSGWNiUYWQBiIeJjH2k8ugsf378h+7jQcRmmqEZxctoamp2e9iGWOisJpPeljwybDzJ+1B/6pKtBAaewY556K7/C6SMSaGGH0+JgH2k8uwksICrjvqQABaymFm83o2bbC5P8Z0VVbzSQ8LPj7YeWB/9h4yEALQ1D3A6Rf+3e8iGWMicLfRjjrU2iTAfnI++eWUSaDQWgbzgg2sq9/od5GMMWHEUfOxeT4JsODjk23qunP46JEg0Ng9wIW/us/vIhljIgghEQ9snk9CLPj46KL99kIcCJXC7OYNfhfHGBOGqtDiFEQ8TGIs+PioT3UlBwwfCkBjjyA33WwfnIzpatyFRaNuqWASYMHHZ+dO3hOA1lJ4eNZnPpfGGNORIrQ4wYhHrhKRESJylIj8TETO9L4enqrnW53RZ8N79mBYWRVfbl7Pht5BXn91LnvvM9rvYhlj2smX+TwiMgo4CzgW6NV22ntVL81y4CHgL6r6aaJ55cdPtIv75RH7AdBaIfz3tOd9Lo0xpj1FaNVgxCMXiMgwEXkYmAucBswGrgZ+AhwMHOJ9fY137XRgroj8U0SGJpKn1Xy6gL2GDaIsFGBz0OGbbg6hUIhgMDd+qY3Jdu7abjnft/MJ8BFwKvCoqm6KllhEyoFjgPO8e0s6m6HVfLqAgAhnTZoAQEt1gNv+d7rPJTLGtFGEVicY8cgRx6nqzqp6T6zAA6Cqm1T1LlUdDxyfSIYWfLqI43bdEVQJFcMjMxJuRjXGpEGMeT5ZT1WfyPS9Fny6iJqyUgYEykBgbY+ArfdmTBcRR83HVjhIgAWfLuSSo34AQEuFcPnv/ulzaYwxbWLM88mJFQ5EZK2IvCQiN4jIiSIyIp352YCDLmTfkdsQaFWcQuGtJcv9Lo4xBnfAQS7P52mnGtgHmOR9ryKyEfgQmAXM9F4/VVVNNjMLPjFMeexvOCj9yqsZUFnNsOpahnfrwajaOmpLylKaV1EwyLjansxYv5KNdUFWrVpPjx5VKc3DGNM5iuDk/mg3gFVABfA48AqwPTAeGAvs7aVRoEFEZgMzVfW8RDOz4BNFyHH4ct0aWtVhfv3qra4PrOzG2Lo+7NV3MBP7DaZXWWXSeV527GSO/b8HCJUJ197wOLdc95Okn2mMSZwCrfmxdcI2uHN7zgYmABer6vkiIsAIYBxuMBoH7ATshjvUOiEWfKIIBgLMPPEclm5az9cb1rF4Yz1f1q/m8/rVfLJmBYs31LN4Qz1PLHBHp42r68tRw7bn0CEjqSkpTSjPHfr0ItiihAqFt1dZ05sxXUE+7NujquuBC0XkL8AtwCMi8hJwvqp+AswDHmhLn+jk0jYWfGKoLi6huriE7Wp7fu98q+Pwef0q3lu2hNeWLuTtZUuYtfIbZq38ht+9/zLHDd+BM0fvwoDKbp3KT0QYWdqNj1vXsbEmSFNTM8XFRal8S8aYTlCVLlfzEZEJwJ9xaySzgVNUdX4qnq2qnwEHicihwE3AhyLyZ+A3qrq2XboFyeTTtX6iWaQgEGC72p6cut14/rb/Mcw84efcOvFQ9u47mKZQK/d89gGTHpnGVe+8yLqmxk49+5LjvOV2yoVb/+fZdBTfGNMJjkrEI9NEpAR4DLgBqAGmA3emOh9VfRq33+dXuEvrfC4i/+k1wyXNgk+KlBUWccSw7bhnynH8+8ifcvSw7VHgzk9nst+jf+WpBfFPHJ0wZADSqmgBPP3Jl+krtDEmJgVanUDEwwc/wB3efb+qNgPXAqO9RUFTSlVbVfVGYAdgPnAbcFkqnm3BJw22ranj5omH8PThp7BLr/6satzMua8+xeVvP0dTqDXm/cFAgH7qLpW0tnOtdsaYFGsb7dZVaj7ASGDL/iuqGgIWAikJPiIyREQOF5FficgDIjIX+AJ3gIEAsf+IxcGCTxptV9uThw46kd/uvj9FgSD3fvYhxz57PysbYi6dxLmHuyMbWysCPPSPN9NdVGNMFF1sM7lyoOMSKJuBpOZ+iMjbIrIeN9A8DlyJG+jeB36BW+Oq8WpCSbPgk2Yiwo9HjuWRQ06if0U1c1Yt44R/PcDyzdG3zT5op1HgKE4R3P38jAyV1hjTkWqXa3bbDHQcTlsGbEzyuRNwY8L9wH5AhaqOVdWfquqtqvqqqq5LMo8tLPhkyA49evP4oSczsqaOL9et4fhnH+Cbjesjpi8uKKC6wf3nWVrckqliGmPC6GLNbp/hjnIDQESCwFDaNcUloQz4EfACMF9EHvaa3w4SkV4x7u0UCz4Z1KO0nAcOPIHR3XuxaEM9P33+YTa2NEVMv//IYQA0VQWoX5vshxpjTCIUIeQEIh4+eBnoLiKnikgR8GtgvjdEOhnVuE1rl+DO52kAjgR+BzwNfCMiX4vIkyJylYgcnkxmFnwyrKaklPumHM+w6lrm1a/ioteexYmwTNI5R00CIFQi3Pg/z2SwlMaY9lLd5yMi00TklQjXhojIoyKyxjvuFpG6tuuq2oC7s+jPgdXA/sBxCRWkHVXd4DWt3aKqJ6vqdrgBaW/gQuAeYC1wEG5/0GPJ5GeTTH1QXVzCtP2O5sin7+G5xfO55YM3uHjc3lul61NdSUGT0losvPLV1z6U1BijSkprOCJyGu421K+GudYdt2ZTBFyH+zf6UmCMiOzqDa1GVWcCu6SsUBF4G8u96R1tZSzBXe9tbDLPtpqPT4ZW13LbpMMJiPCn2W/zzrLFYdP18/oV11XlxcKGxnRBqRlqLSJBEbkSmBYl2UVAf2A/Vb1OVa/F3a56R+CUZN5FZ4jIFSJyqIj063hNVRtV9W1VvT2ZPCz4+GhivyGcu+PuAFz2xnQ2tzRvleaMg3YDoLUcZryfktUzjDGdoBCrz6eHiMxod5zZ8RlebWEW7sKd9wBLI2R3AvCKqm6Zla6qL+Cuq3ZCit9aNFcDTwCLRWSFiPxbRK4TkRNEZNtUrHJgwcdnPx+zOyNr6vhqQz03zHp9q+uHT9jBG3It3HLXiz6U0Jg8p27TW6QDWKWqO7c7poZ5SglQBRyvqqcQZqKmiNTgjlqbGeb+WcDOKXtPsQ0Cfgj8NzADdxXrS4H7gE+A9SKS1ARECz4+KwoGuXHvgwmKcOcnM5mx/Pt9O8UFBZQ2uB8y5rXaiDdjMk2BkAYiHnFaDwxX1YeipGlr4gpXK/oWqBKR6vhLnjhVXaKqj6vqFap6sKr2wh3efQ3QBCzADaYJs+DTBYzu3ouzdpiAAr9972U6bhI4plt3ABqq7Z/LmMyL2edTLSJTReSwSE9QVUdVYy1L07Yh2OYw19pWNCjvfPlTQ1W/UNWrcSeg9gQOTOZ59tesizh7zG70KC1n9qpveXbRvO9du/C4HwAQKhWefvp9P4pnTF5zHIl44C7yeaaqPpVkNm1/j6NtUe0kmUfSVPVt4FncfqGEWfDpIsoLi7hwpz0BuH7mazSHQluu7TSkv7vKdRDu/Ne7fhXRmLzk9u1IxCOF2tbcCrcTZWmHNH77EHeuUcJ8Dz6xJlRFua/Om6i1XETWi8grIrJbJsqcLsePGMPQ6lq+2lDP/fM+3HJeRCj3+n0Wauf2BjLGJC/kSMSDOJrd4tQ236JPmGt9gXpv3k3aiciidkvrTBGRHh2S7EiSTYC+Bp92E6p2w51QdRNwOPC8t2xEpPsqgddwZ/XeAVyB21n3koiMTne506UgEOCX4/cB4LY579DY+l0T8Y7d3X/7xirfPy8Yk3di1HxS0uymqvW4WyOMC3N5LO6os0xZAOyLu7TOs8ByEflKRJ4VkdeAnwJvJJOB33/JEp1Q9UtgW+AwVb1KVW/FXQIC3KW/s9b+A7dhVG1PVjVs4vEFH285f+GJkwFoLRGefvI9v4pnTN7J8H4+jwCTRWRk2wkRmYz79+7BVGcWiaruq6q1wDDgeNzKwafAaNzdTacDZyWTh9/Bp9MTqrzJTacAz6jqa+3uW4a7IN7Wk2WyiIhw1uhdAZg69/0t676NHtAbaVUIwt+mW/AxJmNi9/mkqtkN4HpgDfCiiFwkIr8CHsad+3NvCp7fKaq6UFUfVtVfqeqBqjpQVbur6iGqGn5Zljj5FnySmFA1GLeJ7XnvOSIiFQCqeruqRlu6IiscPGRb+pVXsWDdGl5Y/AXg9ftsdj9lfSXW72NMJqkjEQ9SN9oNVV0JTARm486puQB3Y7eDVDXyEvhZyM+aT6ITqoZ7rytE5AagHtggIl+k6JOH7woDQU7b3o29f5n7XS1nbA93HEZjpd8VVmPyS4wVDhJ4ng5W1UkRrs3zJnZWqGpPVT3VC0ppIyL7JXHv5ETu8/OvWKITqrp5r7/FHep3PvAT7zmPR/pBiMiZbWsvrVyZ1n/HlDh+xBiqi0qYuWIpH6z8BoCLT9ofgNZS4dlnwlUYjTGppgrqBCIeOWK6iLzkLSYajJVYRApF5CgReRV3QEKn+fmTS3RCVbH32g3YU1XvVNV7cKuq9cDvwz1IVae2rb1UVxdzJLfvyguLOH7EDgA8MG82ACP79URaFALw9+k238eYTIlR80lln49fxuKuN/ck8K2I3Cci53vBaA8R2VNEDvP6oR7CbZ16GPdD/06JZOjnfj6JTqhqG+f+qKqubTupqvUi8iRwiohUqGrWL4R2wogdmTr3fZ5a+BlX7LovlUXFlDUImwphodMQ+wHGmBTY0rcTyTpV3Wol62yiqnOBA0Rkd+Bs4AjgRLauHAjuOnWPAneoasJLrvgZfBKdUNXWR7QizLUVuD+cCiDrg8/Q6lom9B7Au8uW8MSCTzh55FhGVnRjJvU0VNj+PsZkTIJ9O9nGWzrnba/pbTywHVCH+xNYCcwFPlDVpJf58a3ZLYkJVXNxV1XdPsy1IUAj7g8pJ5w4YgwAD34+B4Dzjp0IQKgUPpj1pW/lMiZvaMzRbjlHVUOq+p7XrXGDqt6oqnep6sxUBB7wf55PpydUebWhJ4FDRWT7dvcNwV0d4QlVDYW7NxsdOGhbqotKmLt6OR+tWsaEUUORkKJB4bb7X/G7eMbkCYly5ESfT8alNPh4c246rgEUTcwJVSIyVEROFpGh7e77Be7ggpdF5Ncicinu5NIG4FepeC9dRUlBAUdv48bYf8yfg4hQ7I0PnLup3seSGZNHnChHCuf55JOUBR8ROQa3JnOJiDwlInvGuifOCVUTcbedndjuvkW468G9iru73uW4q6zuoaoLUvWeuopjh7uj3p5ZOI8WJ8TAgjIANlX4WSpj8oQCKpEPk5BUDjj4L1UdDyAi3XCbxiZGv8WdUAUcHOX6ncCdYc4vAI5NsKxZZVRNHcO7dWd+/WreWLqI06bsxmWvvUxrmbDsm9X07tvd7yIak9MSnUxqIktls1tIRNrGex+NO5jApICIcMTQ7QB4YsGnHLb3GHAULRBumPpvn0tnTB5wJPJhEpLK4HMc8H8i8jpwEnBOCp+d9w4fOgqA5xbPp9kJUeT1+7y7fJmPpTImP4hGPkxiUhl8jgWW4+6t0wQclMJn572Bld0YW9eXza0tvLD4C3o5hQCsK7PffmPSSqPUelK7mVxeSWWfz38Ao1U1JCJzgKeBh1L4/Lx3+NBRfLDyG55Y8CmHjh/FHfM/orVMaGxspqQk4t57xphkRf+Ml/UrHIjIlUk+4m5vIFjcUhl8FgJ7isg7uEszvJPCZxvgkMHb8tv3XuLVpQv4/bFncceNc3AKhTv+9gIXnh1xzIYxJlkpmVbZpV2FG2IT6cRS3F1NF3XmplQGn+OBM3H7fj4gy3cU7Yp6llUwodcA3l62mDe//YqCRmgthWfnzudCvwtnTK5qG2qd+y4EnujkPbWE35MtppQFH1XdANyUqueZ8A4cNIK3ly1m+lfzqW4MsLpUWVnY6nexjMlpkvs1H4BVqvpVZ24QkYTX0Ix7wIG3esFe3r44PxORA0Uk3H47Jo0OGOTupffK1wvYuX9vAJpskVFjTHJ2BxKZt1Hv3dvp1a3jCj4iMh74FHdFgTuA24FngNUi8lcRGdDZjE1i+pRXslNdHxpDrUzYfwQATrHw9DOR1mE1xiRLHIl45AJVfVdVVydwX8i7N9z2N1HFDD4iMgh4HhiBu2PdtcCvgb/jbmHwH8BsETmks5mbxBw0yA0676/5hkCTgsDfn3svxl3GmIRojMOGWicknj6f/8Ld8nqKqj7f/oKICO78nj8Cj4jIkao6PfXFNO1NGTSC3894lReXfEFZo7KxWFhkm8sZkzYx+nyyfqi1H+JpdjsAuKtj4AFQ10PAjsAnwN0iUpPiMpoOBlfVMLKmjg0tzfTs4S4y2lCeG9V/Y7qk6DUfk4B4gk9fIm/sBmxZnfoIoBxbVicjDvSa3vqPrgPczeU+/niJn0UyJieJujWfSIdJTDzBZwNQGiuRqi7B3QDuyGQLZWKbPHAbAOZsWgEtDgSEm+95zudSGZOj8mhLBRHZIxP5xBN8PgPiLcwHwNCYqUzStq/tSZ+ySpZv3kih49b9526wzeWMSYc8q/k8nYlM4gk+jwBHiMiYOJ9XnFyRTDxEZEvtp6zaHTeyyWZdGZMe+dXnk5HqXDzB58/AEuAJERkVI+1EoFMzZE3i2oJPQR93UdHWMmH5cqv9GJNS+dfnk5GQGjP4qGojcAzQDXhfRK4Rkb4d04nI2cBRdH5tIJOg3XoPoLygkG9aNgIOGhRunGYj3Y1JuRyu+YjIQhFZ0HYAVd7XC73v0yKutd1U9QMR2Rd3i4TLgV+KyPvAl0AZMAYYBiwGrktTWU0HxcEC9uk/lGcXzUOCiobg7W+/9btYxuScHN80blK7rwWYDfwg3ZnGvbabqn6AO5/nMtwgsztwMu6W2cNwO6n2UtW1aSiniWDyALfpTSrd79eV5mY7gDG+yuEVDlT1q3bHIiDU/ly68u3Uqtaquhm4AbhBRAYCg3B3uvhUVdekoXwmhh/0H0pAhKYKpWCd0mKbyxmTWpp3Kxx0mQEHYanqYlV9XVXftMDjn5qSUsb37IeD4hQ5aKFw67REFqc1xkSUw30+YdyfiUwSDj6m69hvwDD3iyL3/4Tpn33pY2mMyS1Cfo12U9WfZyIfCz45oK3fxylXFGV1ScjnEhmTQ9RbYifCYRJjwScHDKuuZWBlN5wC0CKludzt9zHGpIgT5TAJseCTA0RkS9ObU+z1+/yfrfNmTKpYzSf1LPjkiLamNy3x+n0++cLP4hiTW7JwwIGIHCsir8eZ7o8i8lMRKehw7Zl0lc+CT47YpVd/KguLcIpBg8oq6/cxJjWybHkdEQmKyMXAvcQYNi0i5wN/wt0w9JfAmyJS2y7J3ukqpwWfHFEUDDKx3xAAQqUOLRXW72NMymRXzecPuPur/SGOtGfj7lJ9GrAdMBN4qV0AStucHws+OeSAgcMBcEoULRBusfk+xqRENtV8gJtVdSKwMI60fVR1NoCqhlT1bOBF4GUR6U4aw6sFnxwyqf9QgiJosaKi/HuezfcxJmnRaj0+1HxE5BgR0TDHnQCq2pkFHleKyJD2J1T1YuBl79jSByQi+6eg+FtY8Mkh1cUlTOg9AMSt/awutX4fY5IldLnRbo/h9tF0PH6WwLNeBE7teFJVLwBeAUranf5d2xciMltE/i4i54nI3iJtq0vGz4JPjtky4bTUoaVcWLtmvc8lMib7daXg4zWPbQxzNCXwuHOIsBOBqp4HDG73/YR2l3+N26y3L+7AhrUiMl9EHoo3Yws+OWZ/b4M5p0RxgnDln57yuUTG5IAcnWSqqs3egtGRri+OcP5pVb1GVY9U1UFAL+DnwIx487bgk2MGVHZjZE0dBECLlXdXLfO7SMZktzQsryMi00TklQjXhojIoyKyxjvuFpG6JN5B2qnqalV9TlWvj/eeTm2pYLLD/gO34bO1KwmVKhsr7POFMclK5ag2ETkNOB14Ncy17rgd/UW4zWEFwKXAGBHZVVXjnj+hqncCdyZYRgEuBg4HioH5wAfALGCWqq5L5Lnt2V+mHLRlyHWpQ2up8vpbn/hcImOyXApGu3mTP68EpkVJdhHQH9hPVa9T1WuBY3A38jwlobIn5lLgeqAvEAJ+hLuX2wvAGhH5sjP9O+FY8MlBo7v3on9FNQRBi+HGf77id5GMyV6xVzjoISIz2h1bbSwnIiW4tYargXuApRFyOwF4RVU/3ZK96gvAPO9apvwUeAvYFrf2A3AW8HugGXcU3JRkMrDgk4NEhIMGjQDAKVW+osHnEhmT5aLXfFap6s7tjqlhnlACVAHHq+opQGvHBCJSAwzFXWWgo1nAzil4J/EaCPxDVUN8V7/7QlUvBw7ALf/IZDKw4JOjDhzsBp9QqUNDJbS2bvW7boyJQxybyVWLyFQROSzKY9YDw1U1WlNVP+81XK3oW6BKRKoTeQ8J2IRbw6HdawmAqr6OOz/oimQysOCTo8bW9aVnaTkUgFMGv7Mh18YkTFQjHsA6VT1TVSP+T6aqjqrG+gTYNlEz3NDntuaL8s6XPiFfAsMAVHWDl/+Adtff47vmuIRY8MlRAREOGrwtAE6Z8q8F8SzzZIzZSuZWtW77exxtGEOmZhY9DxzS7vv3gEPbfd8T6JZMBhZ8ctiBg75reltX1TWX3zUmK0Tv84mn2S0eG7zX0jDXSjukSbdbgGtEpC3fO4BDvDlHvwYuBOYkk4EFnxy2a6/+1BSVQgGEKuHuf2w1rcAYE4cYNZ+YzW5xaltNoE+Ya32BelXdlGQecVHVtar6D1Vt8L5/CLgROAn4Le7w6wuTycOCTw4LBgIcuc12gNv09n9vfOBziYzJQmlY4SBsNqr1uOuljQtzeSydWLomEeI6S0R2i1C+X+COgtsDGKyq7yaTnwWfHHfkUDf4hMoa7yE9AAAVXklEQVQcVpbbKtfGdFYco91S6RFgsohsGcYsIpNx59s8mPLcvu9HwP8Ce0ZKoKpLVfUdVd2YbGYWfHLcmB696V1UDkFoqVXefHee30UyJvuoRj5S1+cD7qoCa4AXReQiEfkV8DDu3J97U/D8aE7Cncx6c7REIvIzEfmNNy8pYRZ8cpyIcOL2OwEQKlOuud92NzWmU2KPdktVnw+quhKYCMwGrgEuAB4HDkpwy4TO2Bn4p6rGakx8HneOzzHJZGbBJw8c4TW9OaXKktIWn0tjTPZJdbObqg5W1UkRrs1T1YNVtUJVe6rqqV5QSrdqYFGsRKq6AHeNt0NjpY3Ggk8eGFxVw7CSbhCA5h7Ko8+873eRjMkqKVjhIBusAbrHmfY1YFQymVnwyRM/HbsLAKEK5abpb/hcGmOyiBKrzydlzW4+mwfsE2fabwk/JDxuFnzyxJHDtiMQcjeYW9U9y7dfNCbDMjjazU+PAweJyISYKSHpze0s+OSJisJiDh3kjt5srVWuve1Jn0tkTHYQMjPPpwv4C/AN8KiIjI2R9gfAF8lkZsEnj5w5zv1AEyp3ePSL+T6XxpgsEa3JLebAsOzhrWZwDO4Cp2+JyA0iMqxjOhE5D3dbhceTyc+CTx4Z3b0X1c2FEID1fZVlK9b6XSRjskKeDDhAVd/D7fdZhLuN9jwR+UBEHhSRR0Tkc9x13xYQYz5QLBZ88sxlEycB0Fqt/Md/3+dvYYzJEjGa3XJlwAEAqvoB7rbdZwMfAWOA44CjcDe7exrYx9tqIWEWfPLMUcNHE2gFLVK+7N7od3GM6foUCGnkIweparOq/llVx+Jun7AzMAGoU9UjVPWbZPOw4JNnSgsKOXnYGMCd83PR79O9XJQx2S9PBhyEpaqrVXWWqr6vqilrqy9I1YNy1QvLPqLVCSEIIiAIARECCCIBAiIIQlCEgPd9gQQokCAFgaD7dSBIUaCAwkABxYFCSoKFFAcKEBFf3tPFe07i7gVz0BLl35u+9qUMxmQTcXI7yojIHsA8VV3dyfuCuDWijzrbDGfBJ4Zr5z7KptbUL6kkCKXBQkoLiqkoKKayoJTqojK6FZZTU1xOXXEVPUuq6F3Sjf5l3akuKktZ3tXFJexR0pe3Wr6hoY/y4BNvc8IRu6fs+cbklO82jctlrwM/Bu7v5H3dvHv3B17qzI0WfGLYr/cONLY2o6j3O6g4qlte274OqYOj7mtIHVqdEK3q0KohWp0QLU6IZqeVJqeFxlALzU4rm0PNbA41s7op9geGqsJShlT0ZHhlb0ZU9mWHbgMZUlFHQBJrOb3l6COY8MAdOKXKte+8acHHmAgEkOh9O9UiMhV4KosHHQjQXUQGdvK+Wu/eTrPgE8Plo49Oy3ND6tDQ2symUBMbWxrZ0NrA+pYG1jZvYk3TRlY1rWd54zq+bahn6ebVrG9pYPbar5i99qstz6gsKGHn7sPYrcdw9qobSV1JVdz59yqrZAdq+UjWsHGActtdL3DOKZPT8VaNyXoSfT7POlU9M1NlSaM/ekdnJVQvtODjk6AEqCgsoaKwhF4l1VHTqiqrmzbw5cblzN+wjE/Wfc2c+sWsaFzHy8s/5uXlHyMIO9YM4oA+Y5jSZ0cqC8NtA/99d570I8bfextaqty+8EPOwYKPMVtRhRzv8wGuTvL+BZ29wYJPFhARepRU0aOkigk9hm85v3TzGt5dNZ+3Vn3OO6vm8+HaRXy4dhG3fvYvDugzhuMH7cGIqshr/3UvKWNy6UBeaF5MQ1+Hq/7nca4678hMvCVjskquj2pT1WSDT6dJ7H2Dcs/OO++sM2akdTv0jNvY2sgbKz7jqaUzeX/1l1vO795jBKcO3YextUPC3tcUamXUtJtxiqFoRYA5559LSUlxpoptTFcRsd+iqqKfTtjp7Ig3vvDm5TNVdee0lCqH2TyfHFFRUMKBfXfif3c5jYf3vogTBu1JSbCQt1d9zs/em8YFM+7k8/VbzwsrDhZw1iB3DcHmHg77X35HpotuTNeXB2u7ZZoFnxw0sLwHF406hCf3uYzTh+1LebCYt1Z9zslv3cbVcx5mVeP676X/xUH7U1UfhAB8PayZB594y6eSG9M1iaMRD5MYCz45rFtRGWcOn8xj+1zCjwbvRaEEeeabWRzz+s3cveA1Wp3QlrQvnnEm0gJaAld89hYtLbbdtjFbRK/55MzCoplkwScPdCsq54KRB/PgXhcwsecoNoeaue3z6Zz81p/4YM1CAOoqKjm1ZjsAWuocdv3NrX4W2ZguQ1SRUOSDHFtYNFMs+OSRAeXduXHcj7l1/Kn0L6tlwcYV/Oy9afxu7qOsa97Mb449lF6ri0Bg7TCHQy+/3e8iG9M1WJ9PylnwyUO7143g/j3P5/Rh+1IoQZ78egbHvXELzyydxZsXnktxvUAA5g7cyCV/fMjv4hrjrzxc1ToTLPjkqZJgIWcOn8x9e57LuJohrG3exNUfPcw5M/7GP398HIEGoBAeKVvEOTc+4HdxjfGVqEY8TGIs+OS5wRU9uWPX0/nNDsdQU1TOrLUL+dnMaZw4fhDBZgcthKerl3DK7+/0u6jG+ETBcSIfJiEWfAwiwiH9xvHQXhdy9IBdUVWeXzubwUNbKCtsggLl1Z4rGH/lTTYKzuQfxfp80sCCj9miuqiMX25/JHfv8XPG1QxhQ6iBih4NdO+xnpLyJlYPaWXUn/7IPU/bPCCTX2KMdjMJ8D34iMgQEXlURNZ4x90iUpeu+0xsI6r6cseup3PTuB8zrKIXwUKlqlsD3Xuup6h/A1eueI0dr76RRUtX+V1UYzLDaj4p5+vCoiLSHXgZKAKu88pzKTBGRHZV1eZU3mfiJyLs3XMUe9RtywvffsRdC1/liw3LqKhqpLyykaaaQg546S8Uf1XEjZOmMGWvMX4X2Zj0UIWQ9e2kmt+rWl8E9Ad2UNVPAUTkXeB54BRgWorv67RvGpagxPeLJ53eU2nr9OGeIAhI29Plu3O4QWLLfyJ4G3xvOR+QAAECBCRIQIIEJUCAYNxbeAclwJS+O3JAnzG8vepzHvrqbd5a+TklpS2UlLbgdN/MJcsf5cK/PUXtN2Vcse9kpuyxQyd/DsZ0cdFrOLmwmVzG+bqqtYh8CSxU1ckdzn8GLFXV/VJ5X5vOrGp9w9zjaNGGuNKmQrR/DQ0TmlS3Pud4D1LvHvdgSwh1NIAQRAgSCBQiUkRBoIgCKaYgUERhoJgi7ygOlFASLKUkUEpJsIzSYBmNIeXDNV/z3KLPaCpUWp0grRoAhNaWAC2NBbCugMr6IibWDOXiow6grjb+je6M8UHET2PVxb11j34nR7xx+sKbbFXrBPhW8xGRGmAo8HCYy7OAQ1J5X6LO6/4hhQltEpt9nC0B67vDARwV99WBUEhwEHaoghPHCCEVWlVoIUCzE6BJAzQ7QRqdIA1OAQ1OkGe/+TMbFxeyobmYDU3FbGgoYfOmEloaKgg0daNboA+Duw1gzID+7LztUCrLY2+EZ0zmKKg1u6Wan81u/bzXpWGufQtUiUi1qq5L0X0JaaUYR1tT8aiYOh/jtq4nSZjrstX1784L0NYCF4hUgO/tpJWemrKjXqBbD5vWecEOCHmBTxFCCg7i9vPy3Xn1zrc9p+3dtdX80O+XWju88r1rsvW1ON5yuFqp8df0z3fh0ilTU/MwG1iQcn4Gn0rvdXOYa23tXOVAxyCS0H0iciZwJsDAgQPjLmR5n4/iTpvNVB0gBIRQpwWlCdVG1GnEoQFHm1BnM6oNOLoZdTai2oDqJlQ34zgbCTkbaHU24DjuObQRoRENNYE2UyAhCkUpECXY9goEcQPflqGXGQh2JvfNqA73+TQBig04SAM/g0/b35pof13C/YsndJ+qTgWmgtvnE08B84lIAPdHW4gES/guxqefqqLaSEtoLas3fsvC5Yv5eu3XbG5ax+aWdYS0CZVGNNCMSCsirQTEIRBwCEgIRAkElAAOAVHEO9z3pQSkrRHxu1qehKsVhtkrOdy5rdIk9e5Numyo3y11D7OaT8r5GXw2eK/hGvhLO6RJxX2mixIRREopDpTSt6YvfWvG+10kY76jCqFQ7HRdhIgcAfw3MAD4HLhAVd/wt1Rb83OS6WLvtU+Ya32BelXdlML7jDEmMVkyyVREhgJ3A2cD3YBbgCdFpNrXgoXhW/BR1XpgITAuzOWxQNix0IneZ4wxifEmmUY6upZBwDRVfVVVHVW9D7cbYqTP5dqK38vrPAJMFpEtPxgRmQxsCzyYhvuMMaZz1B2QE+noSlT1ZVW9pO17EdkNqMBtfutS/A4+1wNrgBdF5CIR+RXu/J2ZwL3gViNF5GSvOhn3fcYYkzJdqOYjIseIiIY57uyQbhvcD+pXqOrajBc0Bl+X11HVlSIyEbdd8hrc4dOPA5eqapOXbCLwd+CnwIJO3GeMMclT7Wr79jxG+OGoW/Y7EZEJwJPA7ap6Q6YK1hm+Lq/jFxFZCXzViVt6ALaEs8k0+73LnFWqemC4C9XBHrpbaeSFU57bdHeXWl5HRA4GHgAuVtW/+l2eSPxeWNQXqtqprRdEZEZX+uUy+cF+77qKrjeqLRIRGQT8AzhVVR/xuzzR+N3nY4wxXZvizvOJdCRARKaJyCsRriWzV9mFuCu83CUiG9sdkxIqaBrlZc3HGGPipYA6qav5iMhpwOnAq2GuJbVXmapeAFyQssKmkQWf+KRodUJjOsV+77oCVTQFKxyISBD4NXBVlGQZ26vMb3k54MAYY+IlItNxB39EUgI0tvt+qreWZPtnlADvAmNwVyDYD/hCVSd1SJfUXmXZxGo+xhgTRaRRcJ1UAlQBx6vqQyKyqGOCTO9V5jcLPsYYk37rgeGqUTcHy+heZX6z0W4diMgwEWmINDpERP5TRD7z0nwkIickk86YjpIc7WS6IG+dtVi7Usa7V1lOsODTjoh0w509XBLh+iXA7cAc3BElS4EHROT4RNIZ01G70U674Y52ugk4HHheRIr8LJtJu0T3OMtK1uzmEZFRwKNEWP3VC0xXAfer6kneuWnAK8ANIvKwqobiTZfed2OyWN6MdjJbyau9yqzmA4jIKcCHQHcg0nIUh+FWee9oO6Hukra3427atEcn0xkTzgnAK22BB0BVXwDmeddM7sqrvcos+LjG4G7FsAPwZoQ0bcuczOpwflaH6/GmM+Z72o12mhnm8izsdyen5dteZRZ8XP+lqqeo6vIoafoBa1W1Y2fgt97rwE6mM6ajuEY7ZbA8JvPyZq+ynO3zEZHeMZJsVNWNALGWrPBUEt8olHjTGdNRvKOdcmKorQnreuAnuHuV3YQ7+OkX5OBeZTkbfPiuphHJtcDlnXhegPCjUNrOOZ1MZ0xHeTXayWwtn/Yqy+Xgc0aM6x908nkbCD8Kpazd9c6kM6ajvBrtlM9UdXCUa/OAgzNXGn/kbPBJwyZKi4FaESnu8Amkr/e6tJPpjOkor0Y7mfxmAw7iNwsQYKcO58d6r+93Mp0x35Nvo51MfrPgE79ncFeuPbfthIgEgLNxt+R+p5PpjAknb0Y7mfyWs81uqaaqq0XkD8BVXjB5CfghsDdwXNuqBfGmMyaCvBntZPKbBZ/OuQbYBPwcOAr4HDhWVTsugR5vOmO+J59GO5n8ZpvJGWOMyTjr8zHGGJNxFnyMMcZknAUfY4wxGWfBxxhjTMZZ8DHGGJNxFnyMMcZknAUfY4wxGWfBxxhjTMZZ8DHGGJNxFnxMXhKRfiKySkSGJPmcl0UkFG3nXBGpEZFmEXnE+/6v3rptxuQtCz7GFyJyjYioiOweJU2BiHwmIhtFpH+Ki/BH4EFVXZjkc+7F/f/o6ChpfggU8t3CoNcAZ4nImCTzNiZrWfAxfpnjvY6OkubnuFsJ/F5Vv05Vxt7CnUcC16XgcQ8DTcAxUdIcD6zF3W4DVV0MPIC7eKgxecmCj/HLR95r2OAjIrXAlbibq6W6iepC4HVVXZLsg1R1HfA0MFFEena87p37AfBPVW1ud+l+YF+r/Zh8ZcHH+OUL3E33ItV8rgJqgYtVtTFVmYrIAOAw3G0K2p9fJCK3icjpIvK5iDSIyPsisquI9BaRh0Rkg4gsFZFrvb2a2twHBHG3z+joGO9ax714XgPWAOek6r0Zk01sSwXjGxGZCfRX1V4dzm8LzAVeVdXJKc7zDGAqsI2qftnu/CLcIBHE7Q8KAL8G6oF1XnnaNgY8ADhVVe/y7i0ClgGzOpZXRF4DBgJDtMP/bCJyHzBRVQek8j0akw2s5mP8NAfoKSJ1Hc7f7L2en4Y898Ld6G9BmGv9gCmqer2q/gG4A+gPzFXVE1R1Ku7AgmbcAASA15z2MDBJRHq0nReRfl5+93UMPJ45QP9kR9wZk40s+Bg/bdXvIyL7AwcDt6vqx2nIcyiwKEIw+FJVP2r3/efe62NtJ1R1E7AC6NPh3ramtyPbnTsOECJvf90WAC34mLxjwcf46XvBR0SCuLWeVbh9PunQHVgf4dryDt+3eq8rOpwPsfX/O68Bi4Fj2507Abcp7tMI+bWVo0eE68bkLAs+xk8daz5nel9frqpr05SnQ+Tf+9YI52N2jHo1qQdwR7DVishgYFci13poV45QrOcbk2ss+BjfqOoyYCUwWkSqgauBD4FpHdOKyGEi8rGIXCEiX4vIWm9k2ngRedMbifaiiJTFyHY5bu0nHe4FCnBH0x2DG1QeiJK+rRwda1zG5LwCvwtg8t5HwHjcOT11wA9V1QmTbhxuf81CYDDwM+AG4HncEWjNuB34hwIPRcnvK2A3EQmqakprHKo6V0TmAIcAvYAXvQAbSduqDYtTWQ5jsoHVfIzfPgKqcUe2Paiqr0dINw6Yqqr3qmorMBsoAc5S1WWqugZYQuwPVC8BZURfWSEZ9wIHAnsSvckNYDfgC2/FA2PyigUf47e2fp8m4BdR0o0Fnmz3/Y7AG17QQUQE2A74JEZ+/8bt99k7odLG9gBQjvt+HouUyJukujvwrzSVw5guzSaZmi7PmzuzEqhtG4ggIn8FVqrqf3nfb4M7EbRSVVtiPO8xoE5V90pvyaOWYX/gOWAnVZ3tVzmM8YvVfEw2GIc7N2dth3MzO3z/UazA47kR2NMLWH75CfC8BR6Tryz4mGwwFpjV9o2IFALb8/3g87000ajqm8BTwGUpLGPcvBUNjsFdvseYvGTNbiYveQuMzgZ2ab/GW4by/juwTlUvyGS+xnQlFnyMMcZknDW7GWOMyTgLPsYYYzLOgo8xxpiMs+BjjDEm4yz4GGOMyTgLPsYYYzLOgo8xxpiM+3/uxZqsLqROkAAAAABJRU5ErkJggg==\n",
"text/plain": [
"