{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "## STAT5 Dimerization Model" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The goal of the notebook is to study practical identifiability of [STAT5 Dimerization](doi:10.1021/pr5006923) Model with *CICOBase*. STAT5 Dimerization is one of the Benchmark models for [dMod R package](https://github.com/dkaschek/dMod). We have translated the model to [Julia language](https://julialang.org/). [dMod BenchmarkModels repo](https://github.com/dkaschek/dMod/tree/master/BenchmarkModels/Boehm_JProteomeRes2014) contains the model files, experimental data and best-fit parameters.\n", "The model is defined by the following system of differential equations:" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "using DiffEqBase, CSV, DataFrames, CICOBase, LSODA\n", "using Plots\n", "\n", "# constants\n", "const cyt = 1.4\n", "const nuc = 0.45\n", "const ratio = 0.693\n", "const specC17 = 0.107\n", "\n", "\n", "# ode system\n", "function stat5_ode(du, u, p, time)\n", " # 8 states:\n", " (STAT5A, pApA, STAT5B, pApB, pBpB, nucpApA, nucpApB, nucpBpB) = u\n", " # 6 parameters\n", " (Epo_degradation_BaF3, k_exp_hetero, k_exp_homo, k_imp_hetero, k_imp_homo, k_phos) = p\n", " \n", " BaF3_Epo = 1.25e-7*exp(-1*Epo_degradation_BaF3*time)\n", "\n", " v1 = BaF3_Epo*(STAT5A^2)*k_phos\n", " v2 = BaF3_Epo*STAT5A*STAT5B*k_phos\n", " v3 = BaF3_Epo*(STAT5B^2)*k_phos\n", " v4 = k_imp_homo*pApA\n", " v5 = k_imp_hetero*pApB\n", " v6 = k_imp_homo*pBpB\n", " v7 = k_exp_homo*nucpApA\n", " v8 = k_exp_hetero*nucpApB\n", " v9 = k_exp_homo*nucpBpB\n", "\n", " du[1] = -2*v1 - v2 + 2*v7*(nuc/cyt) + v8*(nuc/cyt)\n", " du[2] = v1 - v4\n", " du[3] = -v2 -2*v3 + v8*(nuc/cyt) + 2*v9*(nuc/cyt)\n", " du[4] = v2 - v5\n", " du[5] = v3 - v6\n", " du[6] = v4*(cyt/nuc) - v7\n", " du[7] = v5*(cyt/nuc) - v8\n", " du[8] = v6*(cyt/nuc) - v9\n", " end;" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's load the experimental dataset " ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "data = CSV.read(\"data_stat5.csv\", DataFrame);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Observables, initial values and solver settings are set according to [model's files](https://github.com/dkaschek/dMod/tree/master/BenchmarkModels/Boehm_JProteomeRes2014) and [dMod settings](https://github.com/dkaschek/dMod/blob/master/R/PEtab2dMod.R)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "saveat = Float64.(data[!,:time])\n", "tspan = (0.,saveat[end])\n", "\n", "u0 = zeros(8)\n", "u0[1] = 207.6*ratio # STAT5A\n", "u0[3] = 207.6 - 207.6*ratio # STAT5B\n", "\n", "prob(p) = ODEProblem(stat5_ode, eltype(p).(u0), tspan, p)\n", "\n", "function solve_prob(p)\n", " _prob = prob(p)\n", "\n", " # solution\n", " sol = solve(_prob, lsoda(), saveat=saveat, reltol=1e-7,abstol=1e-7) #save_idxs=[1,2,3,4,5] \n", " STAT5A = sol[1,:]\n", " pApA = sol[2,:]\n", " STAT5B = sol[3,:]\n", " pApB = sol[4,:]\n", " pBpB = sol[5,:]\n", "\n", " # observables\n", " pSTAT5A_rel = (100 * pApB + 200 * pApA * specC17) ./ (pApB + STAT5A * specC17 + 2 * pApA * specC17)\n", " pSTAT5B_rel = -(100 * pApB - 200 * pBpB * (specC17 - 1)) ./ ((STAT5B * (specC17 - 1) - pApB) + 2 * pBpB * (specC17 - 1))\n", " rSTAT5A_rel = (100 * pApB + 100 * STAT5A * specC17 + 200 * pApA * specC17) ./ (2 * pApB + STAT5A * specC17 + 2 * pApA * specC17 - STAT5B * (specC17 - 1) - 2 * pBpB * (specC17 - 1))\n", "\n", " return [pSTAT5A_rel, pSTAT5B_rel, rSTAT5A_rel]\n", "end;" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The [best-fit](https://github.com/dkaschek/dMod/blob/master/BenchmarkModels/Boehm_JProteomeRes2014/parameters_Boehm_JProteomeRes2014.tsv) parameters values:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "p_best = [\n", " 0.026982514033029, # Epo_degradation_BaF3\n", " 0.0000100067973851508, # k_exp_hetero\n", " 0.006170228086381, # k_exp_homo\n", " 0.0163679184468, # k_imp_hetero\n", " 97749.3794024716, # k_imp_homo\n", " 15766.5070195731, # k_phos\n", " 3.85261197844677, # sd_pSTAT5A_rel\n", " 6.59147818673419, # sd_pSTAT5B_rel\n", " 3.15271275648527 # sd_rSTAT5A_rel\n", "];" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's run the simulations and plot results for best-fit parameters values" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdd3xT1fsH8HNH9k7apgMKLS0tFFr2KHsP2Utkg4qiTBFB9KsookwBRVFBpqAgsgVUNsjeo0DLLqV7pNnJHb8/rvZXS1sKJE1u8rxfvny1aZo8DW0/Peee8xyMZVkEAAAA+Cvc0wUAAAAAngRBCAAAwK9BEAIAAPBrEIQAAAD8GgQhAAAAvwZBCAAAwK9BEAIAAPBrEIQAAAD8GgQhAAAAvwZBCAAAwK9VUhAyDPPuu+9W/P4syzIM4756/ARN054uwRfAy/jiWJaFbo4vDr4V3QSrnO9Oq9Wq0+ksFksF709RlNPplEgkbq3K5xmNRoVC4ekqeA9exhfncDgQQkKh0NOF8Bt8K7oJTI0CAADwaxCEAAAA/BoEIQAAAL8GQQgAAMCvQRACAADwaxCEAAAA/BoEIQAAAL8GQQgAAMCvkZ4uwPPsNHpkZlPN6KGJfWhCqWY21cw+NKF0C6sQYMFSFChGgWIsRIoCxVigGOklmF6CAiVYoBgRmKerBwAA8GL8NAiv5LFLrjFJBexDE5trQ2EyrKoMVZNj4XJUX4f1roZXk6NgKWZwsJlWlG1ls20o3YLuGdkz2SjDwmTZULaNzbUhnfifdAyRYlxeBktRkBgFSrAgMQqTYWLC018qAMB1LBbLjRs3PPjsUqnUU8/uWRiGxcfHk6RbMsvvgvBEJvvFZfpCDppcBx8bi4fLUbAEw8sY2OlEWKQCIVT6h2kWZdtQjo3NtKJMK5ttRVk29mQmyrahbBuTYUVZVramCmsaiDULwpoEYbGqMp8IAMALGzdunDlzZnh4uEeenWVZDPPTXyIpKSnbt29v166dOx7cj4LweAY77wp9LR9NjsM3t8clL/ylExgKlqBgCVZHg0oNS4pBtwzs35nswcfsgivMIzNbV4s1DMBaBmOtg3E9NFIFgG8oiurfv//y5cs9XYjf6dSpk/t6jvt+EDIs+j2VmXOJMTnRe/H4to44WVkrhEgcxWmwOA02NhYhhAqd6Eou+3cmuy6FeetvWoCjhgFYSz3eQo81DMBePJgBAAA8B1/+7etk0M93mHmXGRGBJtfBh0bhnl3bohSglsFYy2BsOsIRQo8t7N+Z7PEMdsZZ5lIuW02OtQzGuFCsrfHX6Q8AAKh0PhuEW+8z75xiopToq0SiQ6g3xkqoFBsYgQ2MQAghK4XO57BnstndD9n/nWcsFNs0EGsShLcJxlrosUobwgIAgB/yzSDMtKI3jtPbOpItg70xAp8kIf8ZLHLvplvQ6WzmdBb77mnmrpHtXAXvURXrWhXXiTxbJgAA+CDfDMIPztFjauJ8ScEnhUhRn2p4n2roi8Yo24b2pjI7HrATTjojFViPcKxnON4gAOZOAQDANXwwCC/nsb8/ZG4OFHi6ENcIFKMR0fiIaEQxxKksdncqM/QQbWdQ5zCsRzjWOQwXwVZFAAB4AT4YhNNO0580JFRCT9fhaiTOTZ8Scxuju0Z21wP2q+vMsEN0kyCsR1V8QAQWJoNRIgDgPx4/fpyWlqbRaCIjI3Ecz8/PpyiqxH3kcrlEImEYxmAwKBQKkiQdDofBYChxNwzDAgICjEZj0SMQBKFUKovuQNN0YWEh9wgvUvPZs2cLCgo6der0Ig/yTHxtGcb2B0yaGY2p6WtfVwmRCmxSHfyvbuS9wYKxsfj5HLbuVipuCzXjLH08g2U9XR4AwOMcDkf//v0bNGgwY8aMnj176vV6mqZ79+4dHx8fHx8fFhYWFRXFvf3jjz8ihHbs2KHVapcuXYoQ+vPPP7kPxcbGBgcHc283bNgQIdS5c+fw8PAaNWrUqFGjV69exZ9x06ZNWq32+++/f8HKDx48+Ouvv77ggzwTnxoROhg0/QzzdSLhP8sstSI0MAIfGIFolriUy+56yEw+RT80sV2r4D2rYS3UmMLTFQIAKgfDMDabragH29q1a5OTk+/fvy8WixFCycnJOI4fPXqU+2hiYuLIkSPfeOONok9fvXp13759V69ePXXq1B49eqSnpyOE9u/fP3ToUO7tIps2berevfuTBRQ9wttvv/3Uas1ms1AoFAgECCGn04lhmJvap1WETyXGsutMTRXqHOaPM4QEhhoGYLMaEOf6kCd6kQ0CsO9vMLE7hAMP0DsfME7G0/UBANymS5cu06dPr1WrVs2aNTt27JiTk4MQSk9PDwkJ4VIQIVSzZs1y2rNlZmYePHjw+++/NxgMZ86cKf/pTCZTamoqw/zn10pqaurp06dXrlz58OHDy5cvl/W5q1atGjRoUJcuXaKjozds2JCWlta9e/eaNWtWr179jTfecDgcz/Blu47vjAhzbGjeFfpoD9/5ip5blBKbXAebXAd/lGfaly1ddJV5/Tj9ciQ+LApvEuiPfyUA4D5OBh3NYCvzgkRtDQqV/ucH2Wg0Hjx48MKFC1KpdOTIkTNmzFi5cuWAAQMWLlzYunXrTp06tW/fvnnz5jhe5shnzZo13bt3DwwMHDJkyOrVq5s0aVJOAZMnTyZJ0mQyLV68eOTIkdyNq1at6tu3r1arfeWVV9asWbN48eJSP9dut2/fvv3o0aPNmjVjGKZTp06dO3fes2cPN5G7fPnySZMmPftL8qJ8JzY+vkAPqYHHqOAX/f9TCdjXYvDXYvBUM7vxNjv8MI1j6OVIbEQ0HqmAFwoAF3hkZuddpivzwvyoaHxoVMmf39dff10mkyGEpkyZ0rVrV4RQ7dq1b9y4sXHjxqNHj86bN69evXp//PEHd58nrV27dv78+QihkSNHJiYmLlq0qKxjLrZu3RoSEoIQ2rt3b79+/Zo2bRobG8uy7Pr167mrgyNHjuzcufPcuXNFotI3Prds2bJZs2YIoezs7CNHjkydOnX//v0IoXr16v31118QhM/vRgG75R6TNMBHtky4XFUZNj0Bm56An89h16UwiTupKCU2IhofXANXwmsGwAuIUGB/dvP8L1KtVsu9odPp8vPzuXMqwsLCpk2bNm3atOzs7Pr1669bt27cuHFPfu6JEydu3LixdevWXbt2IYScTue2bduGDh1a6hNxKYgQ6tatW0JCwt9//x0bG3vo0KF79+798ssvmzdvRggZjcbdu3f379+/1EcIDAzk3sjKyiIIYseOHUUfatGixXN+/S/G8/9+LvHOKfrDegQ0XnmqhgFYwwBiQVPij0fMr3fZGWedLfTYiGi8TzVc4FPXiwHwL0lJSdwb165di4yMLHE5MDAwMCoqqrCwsNTPXb169YABA7p168a9KxQKV61aVVYQFqEoKiMjgwvg1atXv/zyy507d+Y+hGHYqlWrygrCotoiIyNZlp0xY0ZERESFvki38YUg3JPK3jehN2vBL/KKEuKoZzjeMxwVOIidD5gfbjLjjtP9I/DhUTxuxwOAP1uzZk1MTIxWq3333Xe52cV58+bdv3+/ffv2Wq324MGD586d+/bbb5/8RLPZvHnz5r/++qvoumDbtm2rVKly9+7dyMjIEndOTU1duHBhq1atcBxfs2aNSCTq3LmzwWDYunXrqVOn6taty92tRYsWERERjx8/Dg0NLadmmUz2/vvv9+/f/5NPPgkICLh+/bpIJBo+fLgLXo5nxPsgpBj03hl6UVMCBjTPQS3k2tbgD03sz3fY0UdpIY4GRmKjovHqcBERAP6YNWvWyZMnU1NTp02bNnr0aITQoEGDNm/evGPHjsLCwqioqNOnT9euXbvo/oMHD46Pj0cIPXr0aMqUKY0bNy76UGBg4KJFi7KysiIjI6tVq1Z8i4VKpVKr1T///DOGYY0bN167dq1MJrt+/fr7779flIIIodDQ0Llz55YahPXr11epVEXvfvLJJ/Hx8du3b8/JyYmNjR0xYgRCqHHjxpU8RsQqZ7WT1WrV6XQWi6WC96coyul0SiRPP7v26+vMjgfM/u68T3R3MBqNCsWz7STkLiL+cpeJVmIjovEhNXC5319EfI6XEZTALYsXCnnf8Om77767fPmytx3Mm5iYOHPmzB49eni6EDfq1KnT9OnTO3bs6I4H53d+5NvRnEs0pKALcRcR5zUhdj5g1t9m3j9Lj66JT66DV4H+bQCAZ3H48OETJ04UvyU+Pt4705rfEfLpRbp/BF5HA7+jXUxMoEGR+KBIPMOKvrtBN9hGda2CT0/A4+ClBsD7lMgbL6HRaKpXr178loCAAA/V8hQ8DsI7hexPt5lr/f1+5s6dgiVoVgPinbrE6ltM1310rBpNjMN7hsP1WAD44fbt26mpqTqdrnbt2iRJPn782Ol0lriPRqNRKpUURT18+DAkJEQikVgsluzs7BJ3wzAsPDw8LS3Nbrdzt4hEorCwMISQyWTKyspCCBEEERwczO0gTEhISEhIeI6aT58+nZmZWaKRqVvxOAinnGJmJBD6p19GBC9KKUCT6uDjauO/3GGmn2FmXWAmxeFDo3ACxocAeCubzdanT59bt27VqVMnIyMjOTk5Ly9vzJgxt2/fRgilpaUpFAru7Ihp06a98cYb27ZtGzRo0Ny5c6dPn37kyJEJEyYghKxWa2ZmJjewEwqFSUlJAwYMuH//Prfdvnbt2tzWwx07dowZM6ZKlSoIoaysrAEDBqxcuZIgnvOIuBMnTly+fBmC8OkOPmZvFLBbOsJZfJVHiKMR0fiwKHzbfWbBVebzS8y78fjwKDgQEQCv4HA4TCZT0c76NWvWZGZmJicnc42tU1NTcRzft28f99FSm26/8sorq1ateu+997p168blJdd0m3u7yI8//vhk0+3atWtfvHgRIZSTkxMXF7d79+7evXuXU212djZ3/BNCyGw2YxhWVi+bSsDXOa5Nd5kJcbiQr+XzGI6h/hH4qV7k9y2JbfeZyE3UvMuMwTOdcgEACCHUpUuXyZMnJyQkNGrUqEWLFhkZGQih7OxsrVbLpSBCqGrVquU03U5LSzt27NiyZctsNtvJkyfLf7q8vLykpKSiCdISdDqdWq0u66M//fRTnz59Onbs2LRp040bN6ampnbo0KF+/fpxcXGDBw+2Wq0V+oJdja8jwpsF7MuREIOe1CYEaxNCXstnv77ORG12DqmBvxuPV4XFpcDPsJTTlnQWMXSlPaOwWgyhCSp+i9FoPH369IULF8Ri8dixY6dPn7527dqXX3558eLFjRo16ty5c7t27dq2bVsUik9as2ZNz549tVrtsGHDVq9enZiYWNY9MQz75JNPCIJ4/PjxggULisaUjx8/njFjht1uv3DhQv369cua2LTb7Xv27Dl+/Di3f79Dhw7du3efOnUqTdPDhg1bvHjxzJkzn+dFeTF8DcJbBjZG9fS7AXero8G+b0l80pD47gbdEBaXAv/DFOZbLhyq1KfEkOS/QYgQevXVV7lpxvHjx3Nnu9esWfPWrVu//vrr4cOHv/vuu8jIyEOHDpW1I3b9+vVLlixBCI0aNapBgwaLFy+Wy+Wl3nP37t3c7Ovhw4e7devWokWLOnXqIIREIhHXiUYkEm3ZsiU5OZnbsP+kZs2acSmYl5d3+PDhd955h2u6XbVq1YMHD0IQVlSBA5kpFAqDD69RYnFpPR2aFEd09MuDIYG/IbRBulEfeLoKpFaruTc0Gk1BQQHXdDswMPCtt9566623CgoKEhIS1q1bV+qRuUePHr1169YPP/ywatUqhBBFUb/99lvR+UolFF2DbNu2bb169U6dOsUFoU6nGzt2LPchh8Px2WefcQ24n1S0iSI7O5sgiGPHjnHv4jjes2fP5/vyXxAvg/BWARujKnu2G3hI8cWlE07SMhJNhMWlAFSKK1euDBgwACF0+fLl6OhoDMO4LOQ+qlarq1WrZjabS/3c1atXDxo0iPt0hFBISMiqVavKCsIiDocjLS2t1K2BTqezrGuExUVERJAkOXz48Li4uKfe2a14GYQ3DWysGn65eqmixaW/pzJfXGI+vchMqI2PjcUlvPxeA4AfNmzYUL169cDAwClTpnCzi3PmzElOTm7fvn1AQMDBgwcvXbq0cuXKJz/RZDJt2bLlyJEjDRo04G7p0KFDaGjo7du3o6KiStz5wYMHn332Gdd0e+3atSqVijv7ECGUmZk5b948hND9+/fXrVu3ZcuWp9YsFApnz549YMCAjz76SK/XJyUlCYXComFlZeLlLyduROjpKkB5cIw74AI/nsHOu0LPvUy/UQufFEdo4KgsANxg9uzZ169fP3LkyOeffz5o0CCE0OjRo7dt23b8+HGbzVa9evWLFy/WqFGj6P6jRo2qX78+Qojr012UggghrVb75ZdfZmVlRUVFRURETJkypehDOp2uRo0aBw4cQAh169Zt7NixYrEYIRQXFzdixIj8/Hwcx+Pi4q5cuVL8uYqLj48v3nJ26tSpCQkJO3fuzMvLi4yM5M5+atq0abVq1Vz7+pSPl023+++nX47EBsGq0afxnm7RV/PYBVeY3anM8Ch8WjzPOpd6z8vIX9B0262g6fYL4ueI0MDGqCEF+aSuFlvXlsiwEkWdS2fUw2vD/DYAvuvgwYMHDx4sfkt8fDw3WvU2/AtCmkV3jWy0En6H8k/xxaWd99L1dWhGAtFCD/+UALwQ72y6HRoa2qhRo+K3cD3YvBD/gvC+kdVLMCn/Cgf/4BaXvlkL33SXee0YLC4FwPXS0tKuXbsmEAhiY2NDQ0NzcnK4ptjFicVibudfWlqaUCgMDAxECN24cePJ62VVq1ZlGCY3N7f4LcW35z969EgsFpdYQRobGxsbG1vxmg0Gw+LFi2fNmlXxT3EV/uXJTQOCrfQ+QET8/+LSzy8xsy8y42vjb9TCxdC5FIAXM3fu3IULF7Zu3drpdJ4+ffqbb77Jy8v75ptvEEK5ubkGg4HLv5iYmF9//dVsNteqVSs8PPzatWsIoaFDh1IUhRBKSkqqVq2aTCZDCH355ZdJSUmffPJJ0cHxv/76a9HbBoOhZs2atWvXPnfu3IuUXVhYuGjRIgjCCoElo74EFpcC4CpWqzU9PV2pVM6cOfP69eu1atVCCDmdTqPRqNVquV5oy5Yt27Jly+HDh4s+a9OmTc2aNbt58+aZM2eaNGly4cIF7vbAwMAVK1a0b9+eezcpKal79+7r169/8nk3bNjQvn378+fPX758+annLuXm5jocjpCQEISQ3W7Pzs72hvnSZ1ty4nA4uD8WSjCZTC6q5+luGSAIfVDLYGxXZ/LPbuTdQhS12TnpJP3IXBnrmQHwAYmJiW+++WbdunW7d+/OTWAWrXMWCARFvWBKtXr16pEjRw4dOnT16tXlP4vFYrlw4cKjR4+efITRo0cPHjx4zZo1ZX2uwWDAMGzq1KmJiYnvvfceRVHjx4+Pjo7u169fzZo1i9LXUyo6Irx3797w4cOTk5MZhmnWrNn69es1Gg1C6MiRIyNGjDCbzWq1esOGDU2bNnVntQhBu22fxi0ufWDCl99gYHEp4AU77TiZdrYynzFaExmmCClxY0pKyvXr10UiEcuyvXr1io+P7969e8uWLXv16hUaGlrWQ6WkpFy+fLl3796NGzdu0qTJokWLyjkL6dy5c2+//fbNmzcTExM3bdrENSO9du3a3bt3e/ToUaNGjY4dO86dO5c7lbdUGo3m1q1bCKGlS5fevn07JSVFJBLt3Llz9OjRly9ffuYXwnUqGoTvvvtuTEzM0aNHKYrq2bPn559/vmDBAqfTOWTIkAULFgwZMmTFihXDhg27desWjrs3pZINbE24RujTqsmxuY2J9+KJb5KY9r9THcPweY3xMF5tPQT+w2AvPPzw78rZkI0QwjBMSAieDMJRo0ZxCYRh2LZt2/bv3//nn3+uWLFi6tSpv/32W1H/lxJWr149cOBAuVxes2bNmJiY7du3DxkypNR7vvrqqxMnTkQIFRYWdurU6YsvvpgzZw5C6Mcffxw8eLBIJKpXr15YWNjvv//er1+/soovOqpi27Zt3D56hBDLsjdu3MjPz3/m18J1KhqEWVlZXbt2xXFcKBQ2atQoNTUVIbR//36SJF955RWE0JgxY2bOnHnixImWLVu6r1yu3Tb8TvQHWhH6X338nbr43Mt0vW3UlDrEO3VhKQ3wOkHSgFkt3/N0FYibouNgGNapUyfuDIpJkyZ9+umnpQYhTdPr1q1TqVTcPTMyMlatWlVWEHKrZhBCSqVy2LBh27dvRwg5HI4NGzbo9XruEfLy8latWlVOEBZN0ubl5RUWFt69e5d7d/bs2eWclVgJKhqEM2bMeP/99wUCgdVq3bp1608//YQQunPnTq1atbgvgCCI6OjoO3fulBOExTNfKpWWM4Iuy60Ctia02/YnMhLNbki8EYvPPMvU3Ex91ggfEQ0T4wCUiaZpgvj/Pxijo6OPHz9e6j3/+OMPHMeXLl3KvetwOPr373///v3q1auX/xS3bt3S6/UIod27d6tUqsWLF3O3W63WgQMHPn78uJzJWE5cXFyVKlWmT59e/Eaj0Vj+Z7lPRYMwNDRUIpFs27bNarUGBwfrdDqEkMFgKD6hrFAoyhresixrs9m4NbucV1999eOPPy7r6bgWa08uzLmcSdSQ4UajZ04x5p3KXMTkViqEvmmIjmbhMy4KVt1k59ZzxqkrbymNz7yMHuQzLdZsNpunS3iKW7duvfLKK4MHD46IiEhNTZ0/f35ZGxJWr149ZMiQ4k3LunTpsnbt2lJ/M48bNy4iIkKv158+fXrt2rVHjhxBCK1atWro0KHFH6FNmzbr168vkXBPmjVrVps2bXAcb9KkSUZGxqlTp1asWFH+p7Asa7FYniMsxWJxOScScyoahK+88sr777/PHcwxa9ascePG7d27NzAw0GAwFN0nPz8/KKjkcZEcDMPEYnHFZ4HL6jX6wE7XCcAUCnEFHwf4UpPMlxSoSwRalcz0O4b3DMfnNCICKusbwZdeRo/wmSDkekx7m3HjxnGbJRBC0dHRc+fOPXr06P79+3U63ZYtW9q0aVN0z+bNmxdNoiYmJvbp06f44/zvf/9LSkri3v7www+LN87u0aPHwYMHb926FR4efvXq1cjISIZh2rVrV6Jl2qefflo04VmcWCyeO3du0QqSmJiYixcvrlmzZtu2bQEBAdx8rFqtLmcTIYZhUqnUXT+JbAUwDCMSic6fP8+9+9tvv8XExLAse+zYMZ1OR1EUy7JWq1Uul1+6dKnUR7BYLBKJpCLPxXE6nRaL5cnb+/1FbbpDV/xx/FxhYaGnS3CLXBs78QQV9JNjyVWaYtz+dL76MlYmu91ut9s9XYULLF++/M033/R0Ff6oY8eOf/31l5sevEIjQgzD2rRp88UXXyxevNhmsy1ZsqRdu3YIoRYtWoSEhHz88cdvvfXWokWL6tSp89TdlC8I2m0DhJBWhJY2J0bXxCedpNemMEubE62C4cIxAF4hIyNj/vz5xW8hCGLBggWeqqciKjo1unbt2v/9739dunQRCAQdOnTgBrDcOt0pU6a0atWqXr16mzdvdmOl0G4b/Fc9HXakB7nrITPiCF1Hg5YlEtXk8L0BgIfJ5fISi1TdvafuxVU0CIODg0u9mBkVFbVr1y6XllQmaLcNntQzHO8Qis+/QjfcRo2Pw2ckELDFAgAPksvlnTt39nQVz8bbg7o4aLcNSiUl0awGxKV+5N1CVHMztS6F8XRFAAA+4VMQQrttUI4qMmxdW2JtW2LRVab979TVPGhVCgCoEF4FIbTbBk/TLgQ734ccXAPvso964zid4+2bvgAAnsenC27QbhtUBImjsbH4gAj8kwt07S3OD+oR4+Pg1F/gMklJST/88IOnq/A7aWlp7ntwPgXhLQML1whBBcEWC+AOzZo1u3jx4vnz5z3y7BRFkSSffmm7UOvWrYuaBrgcxlZK03Sr1arT6SwWSwXv/2RnmQIHqvqzs3CkAH6ZVZzRaISWKAihXQ+ZiSeZ595iAS/ji/OZzjKeBd+KbsKbmUZupQykIHgOPcPx6/3JhgFYw23UrAu0jfZ0QQAAb8KfIDSwsXBAK3hesMUCAFAWPgUhLBkFL4jbYrEOtlgAAIrhTRDeLIDd9MA12sIWCwBAMbwJwlsGNgamRoGLcFssrvUXiAlUe4tz6TWGhsEhAP6KN0F438hGKiAIgStxWyz+7EZuvc803k4dz4AwBMAf8SYIbTSCdtvAHbhTLD5piA8/Qvf8k3pggjgEwL/wIwhZhBgWQXMQ4D7cFouWerzxdthiAYB/4UcQOmgk4EelgMekJJqegF/oC1ssAPAv/IgXJ4OEcMgcqBSwxQIAf8OPy25OhgcjQtbpoAuyaWM+nZ9NG/PpghzEMsJqsaIadQmVztPVgWfDbbFYfoPpuJcaUJVckAiXqAHwWfz44faWIGRo2lhA52fTxjy6IIcuzKMNuYwhlzbk0oZclnIQqgBCqSXUAbhSS6oDWISsl44VbF2Oi2WiGnVFUXWFNeqSWr2nvwxQISSOJsThQ6LwCcfoetuotW2I5kFwmRoAH8SPIHQwrBCv3N9BDGO5fIzKekQbcmlDHm3IYQrzGHMhLlcR6kBCqSHUgbhSI46Kx5VaQh1AKHW4VF7K47Tth1jWmZlqv3PFduOsYdcqRJCiqLqiGnVFNeqSgWGV+kWBZ6cToe+bOo8ViAYeoPtXx+Y3IUQwSw+Ab+FHEFbyiNCZfj//58UYKRBFJwjDahBxTXCljlDpCIUGPUffbwwTBIcLgsNRix4IISo7zX7nmv3OlcI/NyKKEtWoK6xRVxRVV6APf54HB5Wie1XsYl/yzeN0w+3UujZEgwD4lwLAd0AQ/gdLU6bDW42Htqq6j5A17+aOZCIDw8jAMFmzLggh2pDruJdkv3vNcvoPKjtNEB4jjqkvjIgTVovBCH780/iPQDH6rSPx6z2m2x/U6Jr47IaEV0zXAwBeGD9+21ZOEDoeJuf//CWpC9ZP+6ZylrcQKp2kXitJvVYIIdqQa/6/E28AACAASURBVL9z1XHnWsHZ/bQhTxgZJ6pRVxLXlNRXrYRKQAUNjMBb6vHXj1GtdlNr2xDQCB4AH8CPIHQwSOjOIGSdjsK96yznDqr6jJU2aOvGZyobodJJG7Tlnp0xGex3r9lvX8n+5j1BSIS8dS9x7aYwceolQqRoVxdyxU2m1S5qal1iWjxeydevAQCuxY8T6s9ksxNO0Kd7uyW27Xev5/+yWBAaoRkwHpd71wkXLE3Zrp40HtnOFObJErvJmnXFZcqKfzqcZ+0SZb2M94zsqCO0AEerWhPhz37wvV+BE+pdAn6i3YQfI0I3TY2yDnvhHxss5w6o+78tiU90/RO8MIwgublTR2qK+eTejM/GiOs0VbQfKAip7unSAIpQYIdeIr++zjTaTn3WiBgbC9cMAeAlfgShww1BaLt5Pn/TUnGtRvr3f8DFMhc/uqsJq0YLq0Yruw43n/g9Z/lMMjhc3qq3JK4pwuGXryfhGJpUB+8Yho04Qm9/wKxsRYRKYWgIAM/w49eo09XXCKmsR3nr52mHvKMZNNH7U7AIodQouw4L/nidrFlX44FNGXPGGA9uYSwmT9fl7+I02KleZJtgvPF2eut96FAKAM/wY0To4qlRls3/ZbGy2whRdD3XPWjlwQiSW1bjfHTbdGJPxuxR4jpNFe0GCEIjPF2a/xLgaHoC3jYEG3mE/uUO+11LQivydE0AgIrhy4iQFbhuZZ7x0G8IJ+QtXnLVA3qKoEqUZtDE4A9WCoKr5az4OGvpVOulY4iBA4Q8pmkQdrEvGalEdX+jdj+Ebt0A8AM/gtBBu+z0CSrrkfHAZs3gyT6zGwGXqxUdBoX8b42y40DTiT3pn4wo3PcTYy70dF1+SkKiuY2JdW2J8SfoN47TJqenCwIAPA0/gtBlU6P/TIoOJwNCXfFw3gTHxXHNAt/6IuD1T+nCvIw5Y/I2LKDT73u6LD/VIRS70p9ECMVvpY5mwNAQAK/mX0FoPLIN4bi8RQ8XPJa3ElSpoRk0MXjmj4KgqpYNc3N//NSZ+dDTRfkjpQB935L4JpEYcoiedJK2w4w1AN6KH0Hoku0TVHaa8a9fNIOn+MykaDlwuUrRabB8yjJxrUY530zP27CANuR6uih/1K0qdrEv+ciMGm6nLuTA0BAAb8SPIHTB9gmWzf95sbLrMB+cFC0bRpCyxO7BH6widSGZ88cZdq1ibGZPF+V3uG7dHzfAu/9BzThLO2F7BQBehjdB+IIjQuORbSzLylv2dFFFfIKJJMquw/TvLWespozPXjUe2MxSsISjsg2MwC/1E1zPZ1vuom4WwNAQAC/iF0FI5aYb92/SDp3qD5OiZSFUOs2giYET5jtSUzI+f9V8Yg+qlDazoEiwBO3sTL4ag7feTc27zDDw8gPgHfgRhC90jZBl8zcsUnYe4leTomUR6MN1oz7QjXjfcu5A1peT7CmXPF2Rf8EQGhuLn+5N7kllWu+m7hRCGALgefwIQifDPncQmo5sY1lG3qqXSyviN2H1WoETFym7Dsnf/HX2t+870+54uiL/wnXrHhiBN99J/XATrhkC4GF8CUIkfK7OMlRueuH+zdohfj0pWhZxXDP9jO+l9Vrl/PBR7po5VG6GpyvyI1y37kMvkd/fZLr/QT22wNAQAI/hTRA+z4iQZfM3fqns/AoZGOb6mnzCv8tKfxRWjc76cmL+5q8Yk8HTRfmROA12uhfZJhivv43acBuGhgB4hi8HoenodpahYVL0qTChWNFhUPDMlbhEnjn/TeOBzazT4emi/AWJo+kJ+K7O5GeXmEEH6Fy7pwsCwP/wIwgdz76PkMrNKPxrE0yKVhwuU6p6jgmcsNCRmpIx51XziT2IgTFKJWkSiF3oQ0YqUfxv1K6H8LIDUKn4EYTPMSI07PpR0X4ATIo+KzIwTDfqA92omZZzBzIXjbfdPO/pivwF1617cwdiyilmxGHo1g1A5fHRIGQYe/JFaaMO7irI1/27rHR4wW/f5q39gjEVeLoif9FCj13oS0pIFL+VOpIOK2gAqAz8CEIH/WxB6HiUQqgCCKXGbRX5BUnd5vr3lpP6qhlz3zSf2OPpcvxFUbfuYYehWzcAlYEfQehknu08QnvyJVHN+m4rx49gAqGy67DAt74wn9yX8/2HdF6WpyvyF92qYlf7k7l21GAbdR66dQPgTrwJQvJZlrzYki+Ka9ZzWzl+RxAaETRlsTi2UeaXE4wHNkNvtsqhFqKf2hIz6+Hd/6AWXIGObAC4C1+CkBUSFU1C1ulwPLglrFHXrSX5HZyQt+kTNHmJ7eb5rK+mwhmHlWZoFH6uD7n9AdPrTwo2VwDgDvwIwmfaPuG4lyQIjcDFUndW5KfIgJDAt+bKGnfM/npa4b6fWJrydEV+oaoMO/IS2TAAa7CN+jsTRoYAuBg/gvCZVo3CvKh7YZgssbt+2rfOtLtZiyY4HiZ7uiC/QOJoVgPim0Si335q1gUa5kkBcCEfDEJ78iVRNKyUcS9CpdO9+pGy85DclR8XbF3OOmyersgv9AjHzvUh96exff6i82CaFAAX8bUgZKwmZ+ZDYfVYN1cEEEJIUq8Vd95v5vxxcKJT5agqww6/RDYIQA22USdgmhQAV+BHEFb8GqH99hVRRG2MFLi5IvAPXK7WDp2m7jcub+OX+Zu/YmwWT1fk+7hp0q8T8b77qXmXYZYUgBfFjyCs+IgQdhB6hLh2E/305RgpyPzideuVE54uxy/0DMfP9SF3PGD6/EXnwzQpAC/A14LQlnwJVsp4BC6WqfuN045837B7Ve6aOXCcUyWoKsOO9iBrqVH9bdSpLBgZAvCc+BGEjooFIW3IZUwFgrAa7q8IlE4UWUc/7VtSF5Ixdyx0ZasEJI7mNia+ao73/gumSQF4TvwIQmfFrhHaky+JohLg3CXPwgRCVc8xgeM+N53Yk/PD/+iCbE9X5Pt6VcPP9ia3P2D6wjQpAM+OL0HICvCnx5stBeZFvYUgrIb+naWiGnUzF04wHdkOXdncLVyOHetBxsI0KQDPji9BWKGpUXvyJREEoffACUWHQUGTv7Re/Tv763epzFRPF+TjuGnSpc3xvn9RS6/B6b4AVBQ/gtBBP/30CSrrEUIITuL1NmRAaODb86WNOmR9/S50ZasEvavhp3qTP99l+v5FFzg8XQ0AfPAMQXj69OnOnTtHREQ0adLkyJEj3I3Jycldu3aNiIjo2bPngwcP3FNkhUaEtuRL4hjYOOGVMEyW2D3ona8c929kfTnRkZri6YJ8XDU5duQlMlyO6m+jzmTDNCkAT1HRILx69Wr37t0HDhx46NChZcuW6XQ6hBDLsn379m3atOmpU6diY2MHDRrkpiop9unHMNmTL8K8qDcjtfqAN+coO72Su+Ijw65VLOX0dEW+TESgpc2JL5vivf6EaVIAngJjK7aK4eWXXw4PD1+wYEHxG48dO9avX7+MjAyCIGw2W2Bg4NGjR+vXL2VYZrVadTqdxVLRtiMURTmdTolEghCiGCRZ43SOKbdZDMM8/vBl/fTvCJWugk/hD4xGo0Kh8HQVJdHGfMPOlY57NzSDJ4ui4j1dztN558tYQfeN7MsH6Soy7MfWhFrosTIcDgdCSCj0XAU+gdffit6soiPCc+fOVa1addCgQa1bt54zZ47T6UQIXb9+vV69egRBIITEYnHt2rWTkpJcXmJF5kUdabcJpRZSkBcIhUY7dJq679i8nxbkb/6KtVs9XZEvq67AjvYgq8hQ0x3UpVyYJgWgFGQF75eWlrZkyZLVq1crFIrXXnvNbDZ//vnnubm5xf88UavV2dmlbxpjWdZms2k0mqJbRo0aNXv27LKejhsR0jSNEDJSGIkEJpOpnPLs187gEXHl38cPmc1mzGt3VVarIx+/0PrHT+lfjJX0fkMQ7b3T2l79MlbMnLqoqQbvspecWot+K4au/AJgROgSPvCtWPnEYjFJPiXpKhqEarV67Nixbdq0QQh99NFHU6dO/fzzzzUaTfHsKSws1Gq1pX46hmFisfju3btFt8hksnJ+KopPjdpsSEg45XJ5OeXZ7ifJW/WSlHsfP8SybPmvm4fJ5YqhU+13rub/soQJi9QMnIDLlJ6uqRTe/jJWzJBY1CiEHXQAP1Mg/LEVoarcSIIgdAnf+Fb0QhWdGo2KipLJZNzbcrncZrMhhCIjI2/evMldZaRpOiUlpUaN8tqbaYqp+I+Ek3nK3gmWphwPbopq1KngAwKvIqpRV//eclIXkjnvTcvZ/Z4ux5fVVGGnepNhUtR0B3UlD6ZJAfhHRYNwzJgx69atKywsdDqd3333XadOnRBCHTt2ZBhmw4YNCKEVK1ZotdrExESXl/jUtjKOu9fJ4HBcAn8o8RXXlU035iPjgV9zVnxMF+R4uiKfJSbQ0ubEnEZ4572wmhSAf1Q0CEeOHNm8efPq1auHhYUhhBYtWoQQIknyl19++fjjj9Vq9dKlSzdu3OiO+eunLpaxJ18Sw5H0/CesHhs07RtRZFzmwvHmE3ugK5v79I/Aj/QgVyUzww/TJtjGAvxeRbdPcBiGYRjmyQuPDoej/KnOF9k+caOA7b+fThpQ5uXMrCWTVS+NFkUnVPDB/QdPF1s70+/n/7IEI0jNK1O8oVUQT1/Gp7LRaPoZ+o9H7OYORLzWvUsw4BqhS/jqt6LHPVuLNRzHS11+49bv7/JHhKzd6kx/IKxey30FgEomCKkeNHmxtFH7rKXvGA9sRgzM4LkFN036YX28wx6YJgV+jQe9RssPQvvty8JqsZgA/tL0LVxXtslLbDfPZy+bRudlebognzUsCj/ek/wxmRlxmDZDI1jgl3gQhOWfymu7dRFajPoqMiAk8K25koRWmV9OsJw94OlyfFaMCjvTm5SQqNF26iqsJgX+hwdBWP6pvHD0ko/DMHmbPoHj5xsPbcldM4exQM8EtxAT6PuWxPsJePs91MpbME0K/As/grCsESFtzKcLc4VhUZVbEahsguBqQVOWEkpt1sK37XevebocnzUiGj/Wg/zqOkyTAv/C7yC0J18URSUgnAdfBXhBmECo7jdOPWhi3rq5hl2r4FxDN4lVY6d7kSICNd5OXcuHaVLgF3gQIQ6aFRKlr+22J18S1YQLhH5EHNtQ/+4yZ/r97KVTqZzHni7HN0lItKIVMSMB77iH2nAbpkmB7+NBEJY3Iky5LIYLhH4Gl6sDXv9E1qxL9tKp5hN7PF2OzxoRjR96iZx7mRlxmLbA8Bv4NB4HIZWdxtIUGVSl0isCnoZhssTuAW/PM/29O3f1HMZi9HRBvqmWGjvZi2QQarydug7TpMB38SEIWSQobWbUnnxJFNOg0ssB3kIQHB40eQkZEJK1cLz9zlVPl+Ob5AL0U1tiegLeYQ/18x2YJgW+iQdB6KBLP33ClnxJ7MWH2IFKwHXrVr88KW/d3IKty2EFjZuMiMYPvkTOuQTTpMA38SAIS58aZVn7nSvQXxQghMQxDfTTvqVyM7OXvkNlp3m6HN9UW42d7EVSLGqxi0oxwDQp8Cl8DUJn2h1cpiLUAZ6oCHgdXK4KeH2WvFWvrCVTTEe2e7oc36QQoI3tiCl18Fa7qV9gmhT4kIqeUO9BpbZYs6VcFsNwEPyXtHFHQdXovHVz7feTNAMn4lI4otL1RkTjdbXYoAP0nlT2+5aEhAe/QgB4Cn6MCJ9ssUZlpgrCIj1RDvBqguBqQe98RWqDMxe+Zb99xdPl+Kb6OuxCX9LJoha7qNuFME0KeI8fQfjkiJA25BIqmBcFpcBIgarnGO3gKXk/zYcVNG6iEKCf2xGT6+Atd1Gb7sI0KeA3XgQhK8BL7p+gC7LhAiEoh6hmff20b+n8rOwlsILGXUZE43u6kB+eY944TttpT1cDwPPiRRCWOiLMIVQ6T5QDeAOXKXWvfixvDSto3KhBAHa+L2lwoMRd1B2YJgX8xIMgdDxxjZB12FnKicuUHqoI8Im0ccegiYvMZ/fnrvqUMRd6uhwfpBSgX9oTb8TizXdSm2GaFPAQD4LwyREhbcghlDAcBBVF6qsGTVkiCI3MnD/OduOsp8vxTWNj8X1dyZnnmDeO0w5IQ8Ar/AzCghy4QAieCUaQyq7DtMOm5W/6qmDrcpZyeroiH9QgALvQl8y3o8Sd1F0jTJMC3uBBEDpoCELgGqLoevrpyxlTQdaiCc70+54uxwcpBWhzB2JsLN5sB7XlHgwMAT/wIAidTMleo7QhB/ZOgOeDS+TaEe8r2g/IXvae6ch2xMLAxfXGxuJ7upIzzjKTTsI0KeABPgThE6dP0IZcWDIKXoS0ccegd76yXDyS8/2HdGG+p8vxQY0CsAt9yccW1GIndQ+mSYF340MQln6NMNBD5QAfQeqCgyYuFFavlbXwbVsSrKBxPW6adFgU3nQHte2Bp6sBoGw8aBQI1wiBu+CEsuswUY26eRsXSeo2V/V6DSMFnq7Jp2AITaqDtwjGBu2njmWyXyaW0i4RAI/jwXelk2GFxH/mRmE3PXAhUXSCfvpyxlyYtXC88/E9T5fjgxoFYKd7sI+tWMtdME0KvBEvgvC/I0KGZixGQqHxWEHA5+BimXb4dEWHgdnfzoAVNO6gEqKf27BDa+AtdlF7UuHlBd6Ff0FIG3JxuRrhPKgc8Ms/K2guHcv5/kO6MM/T5fgabpp0Rydywgl60knaCatJgdfgQZyUaLEGFwiB+5BafdCE+cLqtbIWjrddP+3pcnxQ40DsTB/ydiHbchd1H6ZJgXfgQRA+OSKETYTAjXBC2XWYbsxHBdu+y9/8Feuwe7ogX6MTod1dyCE18MRd1F6YJgVegIdBCAcwAfcTVo8Nevcb1mHL+nICnQ4raFyMmybd3ol8G6ZJgRfgQRA6SgRhYR4sGQWVABdLtcPeU3QYZF47x3R0B6ygcbkmgdiZ3mRKIdtxD/XYAi8v8BgeBKGzxDXCfBgRgsojbdxR8eZc66Vj2ctn0oZcT5fjawLE6PcuZL/qeKPt1L5HkIXAM/gRhP+9RgiLZUClwtQBgePni2PqZ3050Xr1pKfL8TXcNOkv7cnXj8E0KfAMXgQhK8D/f0M9XQAdt0Glw3FFh0G6Vz8y7FyZt2EBa7d6uiBf0zoYu9iXTDawnfbCNCmobDwIQgdd7PQJlqWN+XCNEHiEMDxG/+4yTCDKnP+W416Sp8vxNQFitKcr2bca3mg79QdMk4JKxIMgLD41SpsKMJEUGkICT8FEEs2giarer+Wuml247yfEwESeK3HTpD+3I189Rs84S9OQhqBS8CEIix3DxBhyCTUMB4GHSeJbBE371vHgVtZX71A56Z4ux9e0CcEu9iUv5bId91DpFk9XA/wAH4Kw2IiQKsghVHAAE/A8QqkJGPuptH7brCWTLWcPeLocXxMoRnu7kl2r4A22Of9Mg4EhcC9vD0IWIZpBZNHUaEEOCUtGgZfAMHmbPoHj5xsPbcldM4exmDxdkE/BEJqegP/cnhxzFKZJgXt5exA+uXcCh6lR4E0EwdWCpiwldSFZC9+2377q6XJ8TdsQ7Gxv8kwW220flQFrdYF78C4IodEo8DqYQKjqOUb98qS89XMLti5nacrTFfmUECn6qzuZqMfqbXVuvgurk4DreXsQ/mfvBEJ0QTZMjQLvJI5poJ/2LZ2Xmb3kHSrrkafL8SkEhmY1IHZ1Jj++wAw6QOdCI3TgUt4ehE903Ibd9MB74XKV7rVZ8ta9spa+Yzqy3dPl+JrGgdjFvmSkEsX/Ru18AEND4DLeH4T/bStjyIX+asDLSRt3DJq4yHLuQM73H9KF+Z4ux6eICTS3MfFrB2LqaWbEYdro9HRBwCd4exAWP5WXsZowDMdEEo9WBMDTkfqqQVOWCKvFZi1825Z0xtPl+JpEPXa+LykhUfxW6nA6LCcFL8rbg/A/bWVgNz3gkX8O+P1fwdblcMCvyykF6PuWxPIWxPDD9BvHaQusTwIvgPR0AU/xnyAsyCHU3rub3uQwm50Wi9NiKvZ/CSkmcRIhpBDKEUIETkhIMUJIRIiEhAAhJBNIcQzHMVwmkHq2fuAOwuq1gt5dVrDlm6wvJ2iHTxeE1fB0RT6laxXsan9ywgm68XZqbVuiUQD29M8B4AneHoSO/4wIcyqz3Xa2JafAXmh2WswOi9lpMTvNZqfV7DAbHSaL02J2Ws1Os9nJfchicphlAuk//wmlUoFULpDJBFILZaUZGiFkdJgQQjRDWykbQshG2520EyFkdloYlmFYxuy0IIQEOCkmxQghMSkW4KSIFCmEcoVQphDKuf/kQrlSJC96l7tFgHv7v6M/w8Uy7bD3rJeO5Xz/obxNX0X7gQiD39cuoxai9W2JX+8xPf6gRtXEP21ICL19ngt4HW//BVr8VF66wI0nEdIs/dDwKCX/bkre3eT8u7fz74oIoUaskQuk0n/jTSaQKkTyELleJpBKBVL5P4EnlQqk3IDvxTlpp422I4RslM3JUDbKbnSY/v8/uym1MK3oXZPDZHSYCu0mASFQFkWjSC4XyhVCuVIoFzBkiDpYK9HoJBqtWCMmRS4pEjwHSb1WwvCYvA0LbLcuaIdOgxNUXGtgBN46GH/9GN14O7W+LRGvhT81wDPgQRAWv0YorBLlqkd20I47BfdT8u5y4XfP8DBQqovWREZrI4eHDozWRqpESlc9V8UJCIGAEKB/p1IryOK0Fo9Go8Nc6DAaHabHpozrBbfybQW51vw8Wz6GsECpTi1Wc7moFat1Uq1GrAqQ6DRitUaswjH4W9qNCG1Q4Pj5pqM7Mhe8pe43Ttqgracr8il6CdrZmViXwnTcQ02tS0yLx3FIQ1AxvArCgmyiTrMXf8zb+Xd3pOw79OB4sCwoWhsZranROaJdlCaCu3rHR1KBRCqQ6GUlL6AajUaFQlH0rpWyZVtyC/7NxTxr/uXMa/k2Q441N89WUGgvVImUGrE6QKLTiFU6iUYr0egk2iBpgF4WpJWoMQS/V14Yhsnb9BFG1M5bP89246xmwHhYBe1aI6Lx1sHY6KP0zofM2jZElBK+acHTeXsQOkqsGn2BCSUH7Tj04PiOlH051rxeUV3W9fhGK9G4pkqekJDicGVYuDKs1I8yLJNvM+Tb8nOsefk2Q64lL82YcSUrKdOcnWnJNtpNellg0X/BMn2wLDBIFhgoDYArlM9KGF5T/95yw64fM+eP0w6bJoyI83RFPqW6Ajv4ErniJtNqF/VJQ+L1WBgZgqfw9l9hT1wjfJ5Vo2nG9N23/9x7d3+kuvqgWr1bVW1GYMTTP83P4Biuk2h0Ek2UJvLJjzoZKtuSk2vNy7XmPzZlXMu+8ee9Q7nWvHRTppAQhsqDQ+R6nUQbINVyb1dRhMI62HJgAqG63zhrVELuqs9kLV5Sdh6CcJiXdhkMobGxeMtgbORheut95sdWRJgM0hCUyfuD8J/OMqzTwTpsuFTx1E8pwrDsxcwrW27uvJl3u2tE++VdFobI9W6r1McJcDJUHhwqD37yQ9zFyHRT5mNTRo4l71bu7XRT5iPjYydDBUi0IXJ9qDw4VBEcKg/WSTQ6iTZYpsdh2SRCCCFJfKIwolb+z19mLX1HO/w9MiDU0xX5lNpq7GQvctFVpsF2amFTYngU/KkBSuf9QfjP1Og/eycq/At0R8reddc2V5GH9KrZ7dPW78P0nfsohPIYbVSMtuQ6JoO9MMOclWXOzjBnZ5gzk3KSM81ZmeZsG2ULluv10v+faA2R64OkgQFSrR+O1AmFJuD1T80n92YvnarsNlyW2N3TFfkUEkfTE/BuVbHhh+ldD9jlLQkdLJ0GT/D2eCg6fYIueIYDmH69uXNHyt6F7T+JUIW7sThQLpVIqRIpnwxIG2XPMGdxoZhpzjrz+GKGOSvDnFVgM2glGr0sMFgWFCwL0ssCgv59W0gIPfIlVBIMkyV2F0bG5a2fb7t1UfPyxGea+QBPFa/FTvcmZ12g43+jlrfAe1WDoSH4D28PwmIjwuwKbiLcmbLvt1u7vur0RZAU2nN7IzEpqq6qWl1VtcTtFEPnWHMzzFmZpqwMc3ZSTvKhh39nmrMzTFkasSpUEVJFERIqDw5ThIQpQkLlwT52DVIQXC1o8uLCfT9lLnxbO+RdUVS8pyvyKVy37l7h7Mgj9JZ77DctCIXA0zUBr8GfIKzYbvo/7x1af23z0k6fQwryDokT3PgPBZX8EHcZ8n7Bw/uG1EMPjnPXIHEM5xbmFF2DDJHrQ+R6/m7z4A74FcfUz/tpviS+har36xjh7T+h/JKoxy70Jd89TcdvpVa3JtqG8PVbBbjWM/+YffzxxwKB4MMPP+TezcvLmzNnzvXr1+vVqzdz5kyl0sWb0B3/PyLMfepSgqOpJ7+7uObLDp+VuqYD8FeplyGLL9K5lXv70IPj9w0PjQ5zgERbXVW1uiq8KB35tTxHVLO+ftq3+b8syVo4XjtihiCkuqcr8ikKAfq+JfHHI3b4Ybp7VWxxM0IKf2z4PYxln+EQk7Vr186cOTMoKOjixYvcLe3atQsNDX3ttdeWLVvmdDp37txZ6idarVadTmexWCr4RBRFOZ1OiUSy5BrzwMQubkbkrpotbdhOktCyrE858/jC5yeXLGg/K7q01f/+qcSGen9QPB0fGzO4t3OseSXWr1ZTVQ1XhlVwbY6nXkbL2f0F239Qdh4ib92b7+1JHQ4HQkgo9KJrvQUONOEEfSGH5VG3bj/8ia4cz/C3UHp6+vz58z/44IMVK1Zwt1y6dOncuXNZWVkSiaRBgwZ6vf727dtRUS7rgob+c42wvMUy5zMuf35y8dy2H0EK+rlSx45O2pltzeXS8V7Bw3Ppl0pNxxC5vpqyithrGgxJG3cUVq+Vt36eLfmS5uVJhNK/+j+4G3TrBkWeIQjffvvt2bNn22y2olvOnz/foEEDiUSCEFKpVHFxcRcuXHBXEBaUuVjmctb1tAeNBwAAIABJREFU2X8v/Kz1B7G6aBc+NfAZAkLAbYJsiBKKbnQyVLopM82Y/tiU/siYfjnzepopPcOUpZWoiy/JqaIIUSKFAnnmz3AyMCxw0pfGPzZkzn9T2XWYPPEl2HfvWly37rHH6cbbqXVtiQTo1u2XKhqEP//8s9Pp7Nev38aNG4tuzMzM1Gq1Re/qdLr09PRSP51lWbvdXr9+/aJbBgwYMHHixLKejpsapSjKZCMRQkYjxZgLzSyOGY0l7nkr/87sMwunNRhfTRxmfOKjfs5kMnm6BK+mwZQapbKOMqboFoZlsq256ebMdHNmemHmlYykDHPmY3OGhJSEyoJDZfoQmT5UFhwqDw6R6aVkJbUJxVr1k8Y2Me1aaTyxT9L7dSKUf4caeuHUaBEpQj81Qz/fJzr+Tk6IoSfGUoS3piH8RD8HsVgsEDxliXCFgjAnJ+ejjz46cuRIidtlMlnxAaLFYilr/hrDMIFAsHLlyqJbQkNDy5nsLrpGiJG0jMQUMqGBYZQabYm73Sm4P+fc4unNJyWGNa7IF+KH4IrCs1IpVVHoPxPsRqMRiTDucuO9goensy+k3ytl2Wp1VdUIVTW5UOaWshS10KRF1kvHCjbME8U0UPd5A5d54HSU5+bNQcgZWxd1rs6OPor/kSnw5m7d8BPtDhUKwqtXrz548KBOnToIIYfDYbPZtFrto0ePwsPD7927x92HZdn79++Hh5e5gR3H8YYNGz5rfdzUKOOwY4KSP0IPC9PeO/TJxIavQwoCdyu69Ng2vEXRjcUX5lzJStqV8kepmzqqq8J1LmrvLqnXShRTv3Dv+sz541Q9RksbdeD7IhqvUtStu/lO6sN6xMQ60K3bXzzbqlGE0MaNGxcsWMCtGjWbzVWqVNm5c2erVq327t07ZsyYBw8elPpH33OvGp1yiq4mxyaEG7IWvBUy++eiO6QZ0yft/+DN+iM7Vm/zTPX7FVhj5hLP9DLmWfPTTOmPCh8/MqWnGf/5j8CJKoqQKorQsH//H6YIUQqf/5/G+ehO/q9fY6RAPXC8ILjacz9OpfH+EWFxSQXsyMO0Toy8rVs3/ES7yQvtoJHJZMuWLevTp0/NmjVTUlJWrFjh8m90bkTIOm2Y8D8tAn+4tK5/TA9IQeBttBKNVqKpG1i7+I35NsMj4+M0Y/oj4+O/H53h3iAwgkvEKorQKsqQMHloFUWIUlShX3OCKjWCJi+2nDuQ880MSf3WqpdGwbmGLlTUrbv+NmpRM+jW7fueeUTITY0W3zhfWFh4586d6OhoubzMQ9Wfe0Q49jjdOAAbpU7NW/O5fsb3/zyjw/jKjrGbeq901/UYXwF/P7qEm17GJ7c83jM8dNCOMHlIxdvl0IX5hl0r7SmX1X3fkCS0cnmRrsKvEWGRK3ns8MN0jArzkm7d8BPtJs88IhQKhSW+m5VKZfHloK71z4jQYS8+Ijz04HjT0IaQgoDXSt3yaLAXphnTHxnTHxkfX8q89vvtv9KM6QxiQmT6YLk+mDsSWR7E9aKTC2WEUqMdOs1++2rBlq/NZ/ar+79FauGsMZeBbt1+wtubC5UahPvuHhhdd4gHqwLATbgjO2oHxBS/sdBuTDdnZpiyMs1Zaab08xmXM8yZGeYsHMO5RAyWB+l7ddc8fCD7dnz1ht1COg3HSGgp7RrQrdsfeHsQOmgkwBFrt2GCf/p9pBamZZlzGoXU82xhAFQapUihFCmePNDKTjtyrXlFk6vXZbbHCcGPc3fZf9kRIFaHasJLHInM647knlW8W/eq1kQ76NbtW7w9CJ0MEhL/GRHuu3ugU0RbHIM5CuDvRITwyY45CCHDteP3dn2XH+YwNgxNtxceyj2ebsrMseZyHcmLLkDqxJoAqY7vR3ZUmqJu3SOgW7fP8fZ/SSfDCnCc/XcfIcOyf947PK/dx56uCwDvparTMqFmY+PBX00/b1J2GixvNZxrzOagHTnFRpC3cm/nWvNyrfkl2q4WBSS/Tu2oHF2qYFf7kxNP0I22U+v4060blM/7g5DbPmHHhGKE0MXMKxqxOlLNg41TAHgQJhQpuw6TNmibv2WZ+ex+zaAJwvAYYRkjyIoEpE6iDZBq+XislcuphWgddOv2Ld4ehA4GCXHEOmyYQIQQ2nf3YJfI9p4uCgB+IIOqBI77wnLuQO7KWeI6zVW9XsPF0ifvVlZAFj+147Ex40pWEncqclFA6iTaAIm22LmPQf5zzQK6dfsSbw/Cf0eEDkwoslG2U2nn3m44xtNFAcAfGCZt3FEc16xw3/rML15X9Rgtbdyxgp9a6qkd6L8BmWPJ405FLj8g3fCFeZ5egnZ0ItalMJ32UFPrEu/G417brRuUjy9BaMcEokMP/44Pqq0WqTxdFAA8g0vl6n7jpI07Fvz6tfn0n+qB4wX6MtsCP1WZAclQ2ZacdFNmjiU315ZfPCB1Yo1OogmQ6oo3CtDLAit4NrI3GxGNtw7GRh+ldz5kvLlbNyiHtweh4999hLhUse/uwQExPTxdEQB8JawaHTR5senYruyv3pU2bKfqMZq79O4qApzkArLE7U6GemxIzzBnFTgMJQJSIZQFSHQlOunwLiChWzffeXsQOhkkxBHjtOdhjvsFD5vBQRMAvAickLfpI6nf2rDrx8z549T93hLXdvvPlAAnQ2T6EJm+RFMqJ0NlmbMzzFkZ5qwMU9aVrKQ/7h7KMGcV2AwBUq1eFhQiCwqWBwXL9MGywGC5XiNWiwgvbdKGITQ2Fm8ZjI08TO99xHhbt25QPh4EITciPGi7265aSwHu7QUD4P0IpVY7dJo95XL+lmXk37s1A94mNB64jCfASa7teInbuSnWDFMWl5EXMq9kmDIzzFl51gISJzRitVqsUouUKrFKI1JpxCqVSKkWq9RilUasVomUHgzLom7dDbZTC5tCt27e8PZcKQrCQ8aUDxt94OlyAPAdougE/XvLTYe3Zi54W966t6LTYIzwil8IZU2xIoQsTmueLb/AVmiwFxrshXnW/GxLbkr+3QJ7YZ61wGAvLLAZPBuWJI6mJ+DdqmIjDtO7HrDftiACXDn9DNzCK77vy8FdI7Q5LVnOwhhdyRZTAIAXgRGkosMgSUKrgt++yVo0QTNwvDAiztNFlUcqkEgFkiqK0HLuw4WlwV7I5WWeraAoLPOtBQX/hiUXilxYqkVKrVj9T1iKVBrJi4ZlvBY7xXXr3ur8rgUB3bq9nLcHoZNhhTj2mDYHCtXQBQoAdyADQgLe+Mx2/VTu+nmiGnXVvcfich6vza5IWNpph9FhMjpMedb8HEuu0WE2OkwPC9NyrXnc7dmWXCdDKYQyhVCuk2h1Eq1CKFcI5QFSrU6iUQjlcqE8UKIr5wycf7t149Ct2/t5fxAiAY5yGWugWO3pWgDwZeK4ZsHR9Qr/2Jgxd6yy0yvy1r2R77aPERFCkUQbINFGqMrcRlJqWHJtdyoelvEa3eV+svfPQrdur8aPIMxmrUGS59/2BACoCEwoVvUcwzVms5w7qB40Xli1pqeL8hjXhmUDoXz6QY1epm0XqtCIn2FkCSqBtwehg0ZCAuVg9iBpgKdrAcAvCMIigyYuspw7kLviY0m91sruI0ttzAbQM4blQ2P+V1eyN9wz96tmzrU+28gSwtKtvD0IKRaRGMrFqHryQE/XAoDf+KcxW9PCfT9lzh2remlUxRuzgRKKh2WbKujXe8yEE/Somvinrf+/W7fFac23FRTYDdwan3xbQYHNcCUro8BuKLAZ8m0Gg70Qx3CFQK4SK0WEUEyK5EKZiBCKCO4NkYgUKoRyIS4QkSKFUC4khCJCKBfKxIRIRIpkAvhTpjxeHYQUgzCEcAzl4nRwaWupAQDug0sV6n7jJPVa/197dxoYRZWvDfxUV1dv6YQkndBpAkYSEpQQSIdNAyEQFhlQhMjijICgGVZ1xlH0ul3FdRQHR3ldGBTceGcGr7IYkAQYRUbgwkAQBUNYBKQhO9m607XfDwVtZDOkO6nq7uf3qVNUOn8qJ3lyTp1zqu6Tpe7dW2ImLdDbu6ldVNC77G7dygSfS5dUtuThm8/UnpUNhBU4VmQbuSZO5JTuJidyjWzTmcZyVuQ4kWvkmliRZQXOzXuaBS8ncm7eY9KbDDRjZSJMeqORNipf0UAbLHqzhTEbaaOZMVn0ZiNtMDPnj5j0RisTYaQNRr3RaogI4emKmg5C5am8hJBqRuocdbU5YADQTozJ6faFbzZt/7xy6ULrkFsjR0xRHg4Kbda23botjNluiY+MjGzz1/WlJityTVxTI9fEiTwrsI2cm5M4VuCqPbW+f/VlKifyrMg2sE28xBtoQ6QhwkgbDTTTot9pNdIGA22INFiNeoNBZ4g0RBj1P59z/p9oQ6TB2ubi25W2g1AmDEWIJNYwcufI0NzAHiAIKBuzZebUF66oeHlu9KQFphv6qV1T0FN2677na3H9Ken9oXRqp3bvbymDtP68w0VRqmRkI9vESTwrcI1cUyPb5DuHu3Bz9PxpnJsTOT+jNIKJaI9nYWo7CCXC6Ei9u5YiFMa4AdRFd7LF3rXQe2h33SdL9QlJMZPvo6Nx594v10dSW8fpl5dK2Z8Hx27dfkapIInNQrOH97Ai1yx43ZyHFTlWZJs4t+9FU3OTkqkX4tbNiqzywiuwgiQ4rPZV45cFNg41HYTKlNHKhrNxIvZlANAEU6+B9tTMxq2rK16ZHzlicuTwSUSHH8+2U3brzkmgZoTBbt16Ha1MiG3zO0iyzIpswDuFmm7BvCQzOqq88WycoOnABggrFGOIGjOt859eZ498W/GX+7kTP6hdUdC7MZraOV4/3KFzrhH+ViqpXY526SjKrA/85q0aD0LC6EhVU1W8hL2JALRFH9clbu4LncZOr/ngpdpViyV3g9oVBTdlt+4tY/VvHZKmbBWrvWoXFE6CIAgrPVVxBLPUALTIlH6T/dG3dWZr+Uu/d+/YSGRZ7YqCm7Jbd3IU6fMZv+4kuoYdRNNByEnEoCOVzbVxFB5kAqBROlNEdP68+HkvuncXVy1dyJ89oXZFwU3Zrft/Rugf/l9pxldiI692QWFA00F4vkfoPRenw5RRAE1jElM6/+G1iJvHVL/9WN1nb8tss9oVBbdsO/Vtvj7GSPp8Jnx5Fv3s9hUEQVjF1cfrNboMEwB+RlGWASPtjy6TBb78pd979mxRu6DgZtGT12+m3xlM3/2VOOffokdQu6DQpekg5CTC6ORa3m1DEAIECV1EVMyUB2KnP9q49ZPq5U8LtRVqVxTcbulKHbhDTwjpv1YoOafp39jBS9OXVZAIQ+osOsZgMKtdCwBcA2NKRueFb5rSnJVLHmjY9LEs4E5X20UbyLIh9NNZuslfM/d+LR5vxEhpgGk6CDmJ6El1ZzqCMmCyDECQoWi9NXeCfeFb/NmTtUvu547sV7ui4DY1Wbd/HJfWiRq4VpiyVTxSjzgMGE0HIS/JOrnGRpl0jFHtWgCgLehONtusJyJvu6fxf/5fzfsviI3n1K4oiFn18qN9dT/eyfSLo25eL0zZKpYhDgNB00HIiYQiVfHERAxYRwgQxAw3Dox96E0mIani5XlN29YSCSvk2i6SIb44zEYcBoKmg5CXCCXW2GQDeoQAwY4yGKPGTOv8h780H/zfiiX3cydL1a4ouCEOA0jbQSgTWaqOl/QUghAgJOjjE+PnvxQ18s6a9xZhYzb/IQ4DQttBKBFRqrIJOsqAIAQIHebMHPtjy3Vma8Ur8zx7tmBjNj8hDv2k6SDkRCII1TaOIAgBQozObI3Onxf3+0VN2z+vevMRvvyU2hUFPcRhm2k6CFlRFMWmTrxMMVg+ARCCmK49Oj/4V0vW8KqlD2NjtoBAHLaBpoOwkath9DEUx6JHCBCyKCoie6z90WVSc1PFK/O9h/aoXVAoQBxeE00HYYO3ysjYZAQhQKijo2Ji71oY89sH69Ytr17+tFhbqXZFoUCJwxMt4vAw4vByNB2Ebq7azMTLPItZowDhwNijj/2Rt0xpzool9zds+lgWsc90AFhbxOGQzxGHl6HxIKyyMHHoEQKED2Vjts4P/pU7dbji5blsGTZmCwwlDn+ciji8DE0HYTNfYzHGybwXPUKAsKK3OeJmPxd9e0HtP5bUrlosNdWpXVGIQBxelqaD0MtXWQ3xMsehRwgQhkzpNyU8uoyO6FTx8jz3rk1YbhgoiMOLaDoIWb46ymjDPUKAsEUZzZ0mzI6b96Jnz9by52c1fvmp1NykdlEhAnHoo+kg5ISaaCaGEELRerVrAQDVMF26x9+/2DbrSaHKVf7szNpVi/mzJ9QuKkQgDomWg9ArsDLhonQGdAcBgBDCdO0RM+WBhCfeZRKSqpc9Vfn6Q837txNJVLuuUOCLw8F2avgGYcpWsbQujOJQu0FY1VxN6+KMEqaMAsDPdNboyBFTEp5aGZk7oWnHxvIXChq3rpY8jWrXFQqsDPlDb13ZZGawncrbGEZxqOEg9NRQtM0k4gYhAFyMovXmzJz4+S/FznxcqCkvf24WxksDJQzjULv33iqbq3V0PHqEAHAVhm6phm6pUb+Z7tm9uXrZk3SMPTJ3grnPYKLT7l/5QUGJw3t76t47LOVtFIbYdc/2090QTaldV7vQblupbq4lVJxBxCJCAPgVdGRM5IgpCU+9H5k7oXHb2vIX7sV4aUCESe9Qu0FY1VwtEZtB4tEjBIDWUMZLO//hL+fHS1+499zqNzBe6r+Qj0MtB2GNRMUxIoZGAeDaGLqlxkx5wP5fy/S2hOplT16YXyqpXVdwC+E41HAQemoEEseIXjyMEADa4JfjpWswXhoQl8bhD8Efh9oNwurmGp7YGAE9QgBouwvjpUt+MV5aflLtuoKbEodHpjCD7dSI4I9DjQZhA9dIUzQnmxGEABAQvxgvfftxjJf6L0L/izi8rVjYXxOUcajRIKzyVMebbbxEaKwjBIDAOT9e+t8fnB8vfVEZL8X+pW3ni8ORXXRji4IyDjUahJWe6niLjZcILbA6BCEABNTP46V3P86Xnyx/fhbGS/3UMg7HFYnBFYetXVB/9OjRDRs2/PTTT926dZs+fXpsbKxynOf5lStXHj58OCMjY/r06TRNB6QsRsf0st2wpobQAkfMCEIAaBeGbqmxdy0UG8559myufucJfedu1uyxWI/fZkocFtyge7dUGlckZsWR5/rRmTatL8Nv7Tf71ltvPXToUNeuXbdv356RkVFZWakcnzlz5ocfftijR4+33nrrvvvuC1RZ/RMy7+w5AT1CAOgAdNT5+aXW7LGNX32G8VI/KXFYNkUfLL1DSm7dsy69Xq/JZCKEyLLsdDrnz58/e/bs48ePp6ennz592mazuVyulJSUEydOJCQkXPrpzc3NNpvN4/G0sixBEHieT1mj30Pejr6+R0T22Nb/l8CnsbExMjJS7SqCHi6j/ziOI4QYDAa1C2kV7qcjTV+v9R7cbc7MsQ6dwCRcp3ZF5wVjU3QL5N1S6ZUDUlYcebYf7dRk77C1PUIlBQkhFEWxLGu1Wgkh27dvdzqdNpuNEJKYmJiamrpjx44AFsdLRCdgsgwAdChlvNT+X3/T2xKq33m86q3HML+0zVr2Dm8tEm8rFkq01zu85k2333zzTUmSJk6cSAgpLy+Pj4/3/ZPdbj9z5sxlP0uWZZ7nCwoKfEdycnKmTp16pa+i9Ag50SpzHoHSeb3ea60TCCFer5dhGLWrCHq4jP5TeoRScGWJwcwMHh9901j+0O76f31S9/kK08BRpoG3UOYItSoK3qZIEzKnB5l2PVl5hBq3SXba5P/uK/eN6YhEZBjmVyevXFsQrl27dtGiRcXFxWazmRCi1+tF8eenYvI8f6VvEkVROp2uf//+viM9e/a8yneUoihCCC8TWuD1ZkuQfu9VxzAMLp3/cBn9p9yCCcrLyDCGrNyIrFz+9FH39nXnFs819R1iGTJelfHSYG+K0Qx5sA+Z3Yu8d1jO/5I4bfKiLF1mbPt+UV0r5j1dQxBu3Lhxzpw5X3zxRWZmpnKkS5cuLpfLd4LL5UpMTLzSp9M0PXfu3FZ+LVmWJUniJUKJLG20BGoyarihaRqXzn+4jP5TLmBQX0Y6qacp6RGx4Zx7x4batx9junSPzL3d1GsQoTrupldoNMUo+nwcvlsqjd+siXuHrb1HuHnz5lmzZq1fv75fv36+g6NGjTp8+HBZWRkhpKSkpKqqKjc3N1CVyYSIEiEcSwXJDXYACHl0VEzUmGmOZz6yZo9t2PyP8/uXNmN+6TXT1L3D1vYIp06dyjDMggULlA9nzZq1YMGCuLi4xx9/fOTIkaNHj/7iiy8WLVoUwBlNvEQYHZE5TJYBAG2h9Iw5M8ecmaPMLy1/bpY5M8eaO4Gxa2V+abBoue7w1iIxK44syqKz4jq6d9ja5RMlJSUtb3QnJCT4RkH37dtXWlqakZGRkZFxpU9vw/KJumb++k/1R3/6fdz8l/Q2Rys/EVoKxsnWGoTL6L/gWj5xrcSGWveOjU3/Lmzv8dLQbootF1p0cBy2Ngj91IYgrHILvdbRh47O6PzwW3RUTLuWF6pC+8emw+Ay+i+0g1AhC3xzybbGrz6TOTbiplsissfqzNbAfolwaIoegSzv8DjU7jZCvEQMNJE5PH0CAIIApWcsA0baF74VO+O/+PKT5c/NOrf6Db7ilNp1BRnL+T1L9SO76MZvFm8rFvZVt3tvTctBKDM6SuJYHYIQAILHhfX4y+io2Kqlj1S99Zj34C7SIWNvIUOJw7LJHRSHWg5CYqYEiqKILuinCwNAuKGjYpX5pREDRtRv/LD8xYKmbWtltlntuoJJh8WhpoMwUvZiyigABK+W46Xc6SNnF804t/oNoeIntesKJr44HO7QjSsS7vyXKAU6Da95i7UOw8uUVcYNQgAIBeef91Rf4975ReXSh5kuyR2/Hj+oWfTkTxm6uTfqik5LukBfMw0HoUSshEMQAkDIoDvZosZMixw5tblkW/2GD+rWLrcOuS3i5jGUwaR2acHBoicTrw/8QKZ2g5CTSAThMDQKACFGGS+1DBjJHj/Y9PXahqJV5sycyNyJens3tUsLU9oNQl4iEbIXfygBQKgyJqcbk9MxXqo6bQehhHuEABDiMF6qOk0HoUXGRqMAEBYuHS+19M+z5k7U2xLULi30aTgIZcqCWaMAEGZ+MV762h8wXtoBNByEkmzG0CgAhKXz46UjpjTv/7q+8H1lvJRkDCEkxPcaVYWWg5CYJI5iQnmXXgCAq6AYg2XASEv/EeyRb5u2r/MWreJvyDKmZRnTMvWxdrWrCx3aDkKZxe1iAAh3FGVMyzSmZdafPqF3HWbL9jdseJ8ymo09naa0TGOPvrqIKLVLDG7aDUJOJCYJO24DAJyn62SL6HpLxKBbCCFCzVn2cImn5Otz/3xDb0swdO9lTO5turE/ZTSrXWbw0W4QCjIVKXEUg28qAMDF9DaHPtsRkT2WSBLnOsaWlTTt2Fj7///C2LsZ05ymnk5DSgZFa/c3vKZo9zLxEjFJLGWIVrsQAAAN0+kM3VIN3VIjR0yROS934gfv4ZL6z1fwlacNSTeYejqNaU5Dt1S1q9Q07QYhJxGDhHWEAACtRRlMxjSnMc1JCJGa6tij33nLSppWPkcE3pDc25TmNPUaSEfHqV2m5mg3CHmJGEVssQYA0BY6a7Q5M8ecmUMu3FD0lpXUF67Uma3Gnk5TmtN4Q5bOFKF2mZqg6SA0CCyWTwAA+Ak3FK9Ou/9zXhkaxaxRAIBAaXlDkee4Hw9e5oZi1x7htouNdoNQkAkjcBSDoVEAgMCjGMPPNxTdDdyJQ+zxQ3Wr3xDOVRl7ZJjSnKYb+tOxndUusyNoNwg5iTAieoQAAO1OFxFlSr/JlH4TIURsqOWOH/SWlTQUraL0hvM3FNMydZaQ3d1Nu0HIS0SPIAQA6Fh0VOxFs2w8+78+t/p1vc1x/oZi9/QQm72h+SDE8gkAAJW0mGUjcq7jbFlJw5bV/Kky5rq0ULqhqO0g5NEjBADQAB196bL9utVvCFUu5rqewb5sX8NBKFM6wavDOkIAAC1puWxfbDzHHfveW1bStOI5IgqG5HRTmtOUPojuZFO7zGug4SCUCC3gMUwAANpFR8Zcbtn+Cp058sKy/X46k0XtMn+FdoNQFnhCUURHq10IAAD8uuBdtq/FmhQyz0pYRAgAEHR+sQ84y504pPFl+9oNQh3PSpgyCgAQzCiDseWyffbIt96yEvf7L8gca0jpbUpzmm7sT8eovGxfw0EocrIeQQgAECJ0EVGXuaG44X2dKeLCsn2nzmLt+MK0G4Q0z8roEQIAhKKfbyi2XLb/z9f1cReW7Sf3pvRMBxXTMV+mDSieJVhECAAQ6q6wbP/5Dlu2r90g1IssQY8QACB8tFy2zzZzJ0vPL9s/V2ns0ceU5jT2zNLbEgL+ZbUbhDqBo3CPEAAgLFFG88/L9uuqvIdL2LL9DZs+0lki7Y8uC2wHUbtBSAvYaBQAAAgdHR8xaHTEoNFElsXGuoAPk+oC+3YBxIgsZcQ6QgAAuICi6KiYgL+rdoNQL7I6TJYBAIB2pt0gZBCEAADQ/jQchAJLGxGEAADQvjQchCJLY7IMAAC0M+0GoUHmaAyNAgBAO9NuEBpFVm/AwwgBAKB9aTsIsXwCAADamUaDUJSJWfbqMFkGAADamUaDkBOJWeaxswwAALQ3jQYhLxOLxFKYLAMAAO1Mq0EoIQgBAKAjaDcIzbKXYjBZBgAA2pdGg5ATZTN6hAAA0P40GoSCTBlFjmKwjhAAANqXRoOQl4hJZnUGDI0CAED70m4QGiU8mBcAANqdRoOQ43mRoolOo+UBAEDI0GjSiCzL0egOAgBAu9NqEHIsq8NMGQAAaHcaDUKBZXkdeoQAANDuNBqEEoehUQAA6Ah6/99ix44d33//vdPxmsPNAAANbklEQVTpHDBggP/vppA4L69HEAIAQLvzt0f45JNPTps27cCBA3fcccfLL78ckJoIIRLP8ugRAgBA+/OrR1hZWblkyZLvvvsuJSVlzpw52dnZ8+bNi4qK8r8siWVl9AgBAKD9+dUj3LJly4033piSkkIIycjIcDgc27ZtC0hZMs8J6BECAED786tHeObMmcTERN+HiYmJLpfrSicLgvDiiy/6PszOzh48ePCVTha5Zok28DzvT3nA8zyuof9wGf2nXECKotQuJLihKbYBTdO6X9ubxa8gFEWxZcvW6XSiKF72TFmWCSG1tbW+I1VVVVc6mRByXY8eVVGdrnICtIYoiriG/sNl9J9yAXEZ/YSm2Aa/moLEzyB0OByVlZW+DysqKrp06XLZMymK0uv1r776aivfudeNyXyPbiYTNt32C8/zuIb+w2X0n/LLyGDALhl+QVNsJ37dI8zNzd2/f39VVRUh5NSpU8eOHbvKaCcAAIAG+dUjTEpKuvPOO8eNGzd16tSPPvpo9uzZnTt3DlRlAAAAHcDfdYTvvvvuH//4x9ra2ieeeOKvf/1rQGoihLhcrm+//TZQ7xa2tm/f3tTUpHYVwU0QhC1btqhdRdA7evRoWVmZ2lUEva1bt2KyTHuglGks7a25udlms3k8nlae/7e//W3Xrl0rVqxo16pC3uDBg//85z/n5OSoXUgQ+/HHH4cPH37ixAm1CwluzzzzjCAIzz//vNqFBLeUlJSioqIePXqoXUio0eheox0TzwDQMfATDVqm0SAEAADoGAhCAAAIax13jzAqKmrYsGGtPN/lctXX1/fq1as9iwp9//nPf3r06BEdHa12IUHM6/Xu3bsX64L8dPz4cVmWle0Yoc127NjhdDrNZrPahQSTiRMnzp8//+rnBOAxTK1hNpvfe++9Ky23v5Tb7Xa73ViM4aeffvrJ4XDo9R30XQ5JsiyfPHny+uuvV7uQ4FZXVyfLckxMjNqFBLcTJ04kJSVhp7pr0r179189p4N6hAAAANqEe4QAABDWEIQAABDWEIQAABDWEIQAABDWtDifsLm5edOmTR6PZ/To0fHx8WqXEzRKS0tPnz7t+3DkyJG+1zt37iwtLc3Kyurbt68apQWB06dPl5WVpaen2+1238ErNcXS0tKdO3d27949NzcXU/h8BEE4ePBgTU1NXl6e7+C+fft8DyK1WCzZ2dnKa0mStmzZcubMmaFDhyYnJ6tQriaJorh79+4jR4507tw5Ly+v5YOrvv/++z179qSmpg4ZMsR30Ov1FhUVNTY2jho1qmXThWsja0xDQ0NGRsbIkSPvuuuu+Pj4H374Qe2Kgsb8+fNTUlJGXiBJknL8kUceSU5Onjt3rsPhWLp0qbpFalNGRkZkZKTRaPzHP/7hO9jY2JiRkTFixIhp06a1bIp///vf4+LiZs+enZ6efvfdd6tTsfbs2LHDbDbHxcUZDIaWx/Py8vr06aO0yZkzZ/qOjx8/3ul0FhQU2Gy2jRs3dni9GpWXl9e3b9+ZM2cOGDAgLS2tsrJSOb58+XK73T537tzU1NT7779fOeh2u51O57Bhw2bMmGGz2b777jv1Cg9umgvCpUuX5uTkKL/EH3zwwRkzZqhdUdCYP3/+c889d9FBl8tlMplOnjwpy/KuXbtiY2M9Ho8a1Wna8ePHRVHs27dvyyB88803Bw8eLIqiLMsPPfTQ9OnTZVkWRTE5OXnNmjWyLNfW1nbq1OnAgQNqla0p9fX1FRUV+/fvvzQIP/nkk4tO3rZtm8PhaGxslGV55cqVTqez4wrVtqNHjyovJEnKzs5WfqJZlrXb7f/6179kWT579qzFYvnxxx9lWV6+fPnAgQOVJvrYY49NnTpVtbqDnObuERYWFubn5yvDTZMmTSosLFS7omDicrm++OKL0tJS35FNmzY5nc7rrruOEDJo0CCLxfLNN9+oV6BGde/eXXmEekuFhYV33HGHctzXFA8dOnT27Nlx48YRQmJiYvLy8tBEFVFRUVfaAaOsrKyoqOjkyZO+I4WFhWPGjLFarYSQO+64Y//+/S1H9cOZb/MdiqISEhI4jiOE7NmzRxAEZWeuhISEm2++ecOGDeTCb8uLmii0geaC0OVyJSYmKq8TExNra2ubm5vVLSlY6PX6AwcOvP3220OGDMnPzxcEgRDicrm6du3qOycxMdHlcqlXYzC5qCmeO3fO7Xa7XC673c4wjO84rufVmc3mrVu3vvbaaxkZGQ8//LBysGWzjIyMjIyMxGW8yL59+7Zu3XrXXXeRC03Rdzc6MTHxzJkz5JIm6na76+rq1Co4qGlusowoir6/zWmaJoQov9DhVy1ZskS5YnV1df3793/vvffmzJkjimLL2Rx6vR7Xs5UubYqiKF50PWma9nq96tQXJNatW6dcvWPHjmVlZY0fP37o0KFolld38uTJ/Pz8JUuW9OzZkxByaatTLhd+WwaK5nqEDoejsrJSeV1RUaH8tahuScFC+UkghERHR//mN78pKSkhv7yehJCKiorW7/ga5i5qilarNSoqyuFwVFdXS5LkO+5wONSrMQj4mmVKSsqAAQP27dtHfnltWZatq6tDs/RxuVwjRoz405/+dM899yhHLv0pVlrdRU3UaDTabLaOLzgEaC4Ihw0bVlxcrLwuLi5u/QMrwEeW5ZKSEuW+YG5u7u7du+vr6wkhx44dc7lcN910k9oFBodhw4YVFRUpr31NsVevXiaTaefOnYQQjuO++uqr4cOHq1hkEPF4PKWlpUqzHDZs2NatW0VRJIRs3rw5KSkpKSlJ7QI1oaKiYtSoUQUFBQ888IDvYL9+/Zqamr777jtCiMfj+eabb5RWd1ETxWKeNtPcptvl5eWZmZmTJ0+22+2vvvrqhg0b8AScVho2bFhOTo7Vat28efPRo0f37t2r/Hk4ZcqUs2fPTpo0acWKFaNGjXr11VfVrlRzXn/99UOHDn366aeZmZkpKSkPPfSQMnO9T58+kyZNcjgcixcvLiwsVNZvLV68+J133rnvvvs2b97sdru3bdumdvma0NDQsHDhwpqamjVr1hQUFMTGxr700ktnzpyZPn360KFDGYb55JNPLBbLV199xTCMIAgDBgzo3r37kCFDlixZsmjRonvvvVft/4EmDB069MiRI+PHj1c+HDRokNIvfOqppz799NPZs2evX7/earWuX7+eEFJTU9OnT5/bbrstKSnplVde+eyzz/BnWdtoLggJIS6X66OPPmpubs7Pz8cC8NbbuHHj3r17WZZNTk6eOnVqRESEcpzn+Y8//vjw4cP9+vWbNGkS/ma81MaNG1vOWrz11luVkTpfU5w4cWJmZqbvhA0bNnzzzTfdunWbOXMmHg6n8Hg8H3/8se9Dq9X6u9/9juf5NWvWHDx4kBCSnp6en5/veyhYY2PjypUrKyoq8vLyRowYoU7R2rN69eqWE15SU1N92bZmzZrdu3enpKTMmDHDt9D+7NmzH374odvtnjBhQlZWlgoVhwQtBiEAAECH0dw9QgAAgI6EIAQAgLCGIAQAgLCGIAQAgLCGIAQAgLCGIAQAgLCGIARQ3w8//PDBBx/4dm4DgI6EIAToaP/85z8v2o9m06ZNM2fOxI7JAKqgn3nmGbVrAAgv+fn5dXV1t99+u+8ITdNJSUnYKxJAFZp7DBNAGBo4cODAgQMvPc6ybH19fcsH3no8HlEUL/tIFp7na2pqoqKiLBZLO9YKEHIwNArQoVJSUo4dO7Zq1arY2NjY2Nh58+YRQpYvX+5wOHieV85xOBwvvvjiQw89FB0dbbfbU1NT9+3bV11dffvtt1ut1ujo6FtuuaWmpsb3nm63e8GCBbGxsQ6HIyoqaty4ccqDWwGgNRCEAB1qxYoVDodj+PDhq1evXr169X333UcIaWpqKi8v9238W1dX98Ybb5w6daq4uHjTpk2SJN15551Tp07t3bv3rl27Vq5c+e9///upp55STpYkacKECatXr37ttdcOHjy4bt26o0ePjhkzxherAPArZADoWKmpqffcc0/LI0uWLCGEsCyrfGgymdLT0wVBUD58//33CSH33nuv7/yCggKHw6G8LiwsJIQUFhb6/lV5JvO6deva978BECrQIwTQopEjR/qe7Z6amkoIGT16tO9f09LSysvLlT5fcXExwzB6vX7LBdXV1RaL5fvvv1elcoCgg8kyAFoUHR3te608fC4mJqblEVmWeZ5nGKaiokIUxd/+9rctP91oNLZ8rB0AXAWCECC4derUyWQyVVVV+XqQAHBNMDQK0NGsVmtzc3Og3i03N9fj8WzYsCFQbwgQbhCEAB0tPT39yy+/XL9+/d69e0+cOOHnu02ePDkrK6ugoGDVqlXV1dX19fV79+599NFH9+/fH4hiAUIfghCgoz377LO9e/eePn16//79n376aT/fjWGY4uLivLy8u+++Oz4+Pjo6esCAAV9++WVERERAqgUIeZR8YekSAKhFlmVJkvy8yVdfX3/48GGTydS1a9fY2NhA1QYQ8hCEAAAQ1jA0CgAAYQ1BCAAAYQ1BCAAAYQ1BCAAAYQ1BCAAAYQ1BCAAAYe3/AEV0XjY7RzgHAAAAAElFTkSuQmCC", "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ], "text/html": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "sol = solve_prob(p_best)\n", "plot(saveat, sol, label=[\"pSTAT5A_rel\" \"pSTAT5B_rel\" \"rSTAT5A_rel\"], xlabel=:time)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Loss function includes parameters transformation to log10 scale:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "function loss_func(p_init)\n", " p = exp10.(p_init)\n", "\n", " sim = solve_prob(p)\n", " σ = p[7:9]\n", " # loss\n", " return loss(sim,data,σ)\n", "end\n", "\n", "function loss(sim,data,σ)\n", " loss = 0.0\n", " obs = names(data)[2:end] \n", "\n", " for i in 1:length(obs)\n", " loss_i = loss_component(sim[i],data[!,i+1],σ[i])\n", " loss += loss_i\n", " end\n", " return loss\n", "end\n", "\n", "function loss_component(sim,data,σ)\n", " loss_i = 0.0\n", " \n", " for i in eachindex(sim)\n", " loss_i += ((sim[i]-data[i])/σ)^2 + 2*log(sqrt(2π)*σ)\n", " end\n", " return loss_i\n", "end;" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we are ready to perform identifiability analysis with *CICOBase*" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " 5.692108 seconds (12.60 M allocations: 585.771 MiB, 2.88% gc time, 93.04% compilation time)\n", " 0.447164 seconds (2.62 M allocations: 91.556 MiB, 3.74% gc time, 4.74% compilation time)\n", " 0.243965 seconds (1.43 M allocations: 49.670 MiB)\n", " 0.145805 seconds (941.54 k allocations: 32.667 MiB)\n", " 0.697688 seconds (4.27 M allocations: 144.614 MiB, 5.28% gc time)\n", " 0.141389 seconds (846.00 k allocations: 29.317 MiB, 7.34% gc time)\n", " 0.181994 seconds (1.12 M allocations: 38.955 MiB)\n", " 0.184535 seconds (1.11 M allocations: 38.560 MiB, 6.97% gc time)\n", " 0.431702 seconds (2.89 M allocations: 100.730 MiB, 2.10% gc time)\n" ] } ], "source": [ "α = loss_func(log10.(p_best)) + 3.84 # chisq with 1 df\n", "\n", "# search CI with CICOBase\n", "num_params = length(p_best)\n", "\n", "intervals = Vector{ParamInterval}(undef,num_params)\n", "p_log = log10.(p_best)\n", "\n", "tbounds = fill((-7.,7.), num_params)\n", "sbounds = (-5.,5.)\n", "for i in 1:num_params\n", " @time intervals[i] = get_interval(\n", " p_log,\n", " i,\n", " loss_func,\n", " :CICO_ONE_PASS,\n", " loss_crit = α,\n", " theta_bounds = tbounds,\n", " scan_bounds = sbounds,\n", " scan_tol = 1e-2,\n", " local_alg = :LN_NELDERMEAD,\n", " silent = true\n", " )\n", "end;" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
9×8 DataFrame
RowparamsLSatusUStatusLBoundUBoundLCountUCountInitValues
SymbolSymbolSymbolUnion…Union…Int64Int64Float64
1Epo_degradation_BaF3BORDER_FOUND_BY_SCAN_TOLBORDER_FOUND_BY_SCAN_TOL-1.70803-1.41987523494-1.56892
2k_exp_heteroSCAN_BOUND_REACHEDBORDER_FOUND_BY_SCAN_TOL-3.1462941036-4.9997
3k_exp_homoBORDER_FOUND_BY_SCAN_TOLBORDER_FOUND_BY_SCAN_TOL-2.48322-1.98024237289-2.2097
4k_imp_heteroBORDER_FOUND_BY_SCAN_TOLBORDER_FOUND_BY_SCAN_TOL-1.85648-1.68961171179-1.78601
5k_imp_homoBORDER_FOUND_BY_SCAN_TOLSCAN_BOUND_REACHED0.190413128774.99011
6k_phosBORDER_FOUND_BY_SCAN_TOLBORDER_FOUND_BY_SCAN_TOL4.159944.272431431684.19774
7sd_pSTAT5A_relBORDER_FOUND_BY_SCAN_TOLBORDER_FOUND_BY_SCAN_TOL0.437580.7684131722430.585755
8sd_pSTAT5B_relBORDER_FOUND_BY_SCAN_TOLBORDER_FOUND_BY_SCAN_TOL0.7207680.9901632311860.818983
9sd_rSTAT5A_relBORDER_FOUND_BY_SCAN_TOLBORDER_FOUND_BY_SCAN_TOL0.4001020.6661052049290.498684
" ], "text/latex": [ "\\begin{tabular}{r|cccccccc}\n", "\t& params & LSatus & UStatus & LBound & UBound & LCount & UCount & InitValues\\\\\n", "\t\\hline\n", "\t& Symbol & Symbol & Symbol & Union… & Union… & Int64 & Int64 & Float64\\\\\n", "\t\\hline\n", "\t1 & Epo\\_degradation\\_BaF3 & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & -1.70803 & -1.41987 & 523 & 494 & -1.56892 \\\\\n", "\t2 & k\\_exp\\_hetero & SCAN\\_BOUND\\_REACHED & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & & -3.14629 & 4 & 1036 & -4.9997 \\\\\n", "\t3 & k\\_exp\\_homo & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & -2.48322 & -1.98024 & 237 & 289 & -2.2097 \\\\\n", "\t4 & k\\_imp\\_hetero & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & -1.85648 & -1.68961 & 171 & 179 & -1.78601 \\\\\n", "\t5 & k\\_imp\\_homo & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & SCAN\\_BOUND\\_REACHED & 0.190413 & & 1287 & 7 & 4.99011 \\\\\n", "\t6 & k\\_phos & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & 4.15994 & 4.27243 & 143 & 168 & 4.19774 \\\\\n", "\t7 & sd\\_pSTAT5A\\_rel & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & 0.43758 & 0.768413 & 172 & 243 & 0.585755 \\\\\n", "\t8 & sd\\_pSTAT5B\\_rel & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & 0.720768 & 0.990163 & 231 & 186 & 0.818983 \\\\\n", "\t9 & sd\\_rSTAT5A\\_rel & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & BORDER\\_FOUND\\_BY\\_SCAN\\_TOL & 0.400102 & 0.666105 & 204 & 929 & 0.498684 \\\\\n", "\\end{tabular}\n" ], "text/plain": [ "\u001b[1m9×8 DataFrame\u001b[0m\n", "\u001b[1m Row \u001b[0m│\u001b[1m params \u001b[0m\u001b[1m LSatus \u001b[0m\u001b[1m UStatus \u001b[0m\u001b[1m LBound \u001b[0m\u001b[1m UBound \u001b[0m\u001b[1m LCount \u001b[0m\u001b[1m UCount \u001b[0m\u001b[1m I\u001b[0m ⋯\n", " │\u001b[90m Symbol \u001b[0m\u001b[90m Symbol \u001b[0m\u001b[90m Symbol \u001b[0m\u001b[90m Union… \u001b[0m\u001b[90m Union… \u001b[0m\u001b[90m Int64 \u001b[0m\u001b[90m Int64 \u001b[0m\u001b[90m F\u001b[0m ⋯\n", "─────┼──────────────────────────────────────────────────────────────────────────────────────────────────────────────────\n", " 1 │ Epo_degradation_BaF3 BORDER_FOUND_BY_SCAN_TOL BORDER_FOUND_BY_SCAN_TOL -1.70803 -1.41987 523 494 ⋯\n", " 2 │ k_exp_hetero SCAN_BOUND_REACHED BORDER_FOUND_BY_SCAN_TOL \u001b[90m \u001b[0m -3.14629 4 1036\n", " 3 │ k_exp_homo BORDER_FOUND_BY_SCAN_TOL BORDER_FOUND_BY_SCAN_TOL -2.48322 -1.98024 237 289\n", " 4 │ k_imp_hetero BORDER_FOUND_BY_SCAN_TOL BORDER_FOUND_BY_SCAN_TOL -1.85648 -1.68961 171 179\n", " 5 │ k_imp_homo BORDER_FOUND_BY_SCAN_TOL SCAN_BOUND_REACHED 0.190413 \u001b[90m \u001b[0m 1287 7 ⋯\n", " 6 │ k_phos BORDER_FOUND_BY_SCAN_TOL BORDER_FOUND_BY_SCAN_TOL 4.15994 4.27243 143 168\n", " 7 │ sd_pSTAT5A_rel BORDER_FOUND_BY_SCAN_TOL BORDER_FOUND_BY_SCAN_TOL 0.43758 0.768413 172 243\n", " 8 │ sd_pSTAT5B_rel BORDER_FOUND_BY_SCAN_TOL BORDER_FOUND_BY_SCAN_TOL 0.720768 0.990163 231 186\n", " 9 │ sd_rSTAT5A_rel BORDER_FOUND_BY_SCAN_TOL BORDER_FOUND_BY_SCAN_TOL 0.400102 0.666105 204 929 ⋯\n", "\u001b[36m 1 column omitted\u001b[0m" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "ENV[\"COLUMNS\"]=120\n", "res = DataFrame(\n", " params = [:Epo_degradation_BaF3,:k_exp_hetero,:k_exp_homo,:k_imp_hetero,:k_imp_homo,:k_phos,:sd_pSTAT5A_rel,:sd_pSTAT5B_rel,:sd_rSTAT5A_rel],\n", " LSatus = [k.result[1].status for k in intervals],\n", " UStatus = [k.result[2].status for k in intervals],\n", " LBound = [k.result[1].value for k in intervals],\n", " UBound = [k.result[2].value for k in intervals],\n", " LCount = [k.result[1].counter for k in intervals],\n", " UCount = [k.result[2].counter for k in intervals],\n", " InitValues = p_log\n", ")" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Julia 1.11.3", "language": "julia", "name": "julia-1.11" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.11.3" } }, "nbformat": 4, "nbformat_minor": 4 }