PLS.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "350b0787",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import scipy.linalg as sp\n",
    "import time\n",
    "from sklearn.cross_decomposition import PLSRegression"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3a8eead7-a6e8-458e-82c0-4ac3e62ccf75",
   "metadata": {},
   "source": [
    "Here, we create a $n\\times m$ matrix $X$ and a $n\\times p$ matrix $Y$, where\n",
    "* $n$ is the number of samples;\n",
    "* $m$ is the number of human genetic variants;\n",
    "* $p$ is the number of pathogen genetic variants."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 84,
   "id": "8e011e10",
   "metadata": {},
   "outputs": [],
   "source": [
    "S = (v.X).T @ v.Y\n",
    "n = S.T @ S\n",
    "eigenvalues, eigenvectors = np.linalg.eig(n)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 88,
   "id": "03dc05a2",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "50"
      ]
     },
     "execution_count": 88,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "S.size"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 56,
   "id": "635f2794",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 55.68946576, 186.73702662,  83.76003807,  36.47353508,\n",
       "        36.47353508])"
      ]
     },
     "execution_count": 56,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "n @ q"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 176,
   "id": "0acd9367",
   "metadata": {},
   "outputs": [],
   "source": [
    "def transfer(mat, model = 'additive'):\n",
    "    if model == 'additive':\n",
    "        return mat/2\n",
    "    elif model == 'force':\n",
    "        tmp = np.where(mat >= 0.5, 1, mat)\n",
    "        return np.where(tmp < 0.5, 0, tmp)\n",
    "    elif model == 'dominant':\n",
    "        return np.where(mat == 2, 1, mat)\n",
    "    else:\n",
    "        tmp = np.where(mat == 2, 1, mat)\n",
    "        return np.where(tmp == 1, 0, tmp)\n",
    "    \n",
    "def affect(mat, weight = np.array([1]), model = 'additive', encode = True):\n",
    "    mat = transfer(mat, model)\n",
    "    dosage = mat*weight\n",
    "    if len(weight) != 1:\n",
    "        dosage = np.sum(dosage, axis = 1)\n",
    "    return transfer(dosage, model = 'force')\n",
    "class Human:\n",
    "    def __init__(self, n_pop, n_var, ploidy = 2):\n",
    "        self.n_pop = n_pop\n",
    "        self.n_var = n_var\n",
    "        self.ploidy = ploidy\n",
    "        self.X = np.random.randint(low = 0, high = ploidy + 1, size =(n_pop, n_var)) # nxm matrix X\n",
    "\n",
    "class Virus:\n",
    "    def __init__(self, human, n_var, ploidy = 1):\n",
    "        self.X = human.X\n",
    "        self.n_pop = self.X.shape[0]\n",
    "        self.n_var = n_var\n",
    "        self.ploidy = ploidy\n",
    "        self.Y = np.random.randint(low = 0, high = ploidy + 1, size =(self.n_pop, self.n_var)) # nxp matrix Y\n",
    "# Returns the means of the matrix for centralization\n",
    "def getMean(mat):\n",
    "    return np.tile(np.mean(mat, axis = 0),(mat.shape[0],1))\n",
    "\n",
    "# This function calculates the variance in the original matrix\n",
    "def getVariance(mat):\n",
    "    totalVar = 0\n",
    "    for i in range(mat.shape[1]):\n",
    "        totalVar += np.var(mat[:,i])\n",
    "    return totalVar*(mat.shape[0])/(mat.shape[0] - 1)\n",
    "\n",
    "# Plot variance graph \n",
    "def plotVar(varXAry, varYAry, varXTotal, varYTotal):\n",
    "    numFactor = len(varXAry)\n",
    "    numFactorAry = np.arange(1, numFactor + 1)\n",
    "    plt.plot(numFactorAry, np.cumsum(varXAry/varXTotal), color = 'c', label = 'varX')\n",
    "    plt.plot(numFactorAry, np.cumsum(varYAry)/varYTotal, color = 'm', label = 'varY')\n",
    "    plt.axhline(y = 1, color = 'g', linestyle = '--')\n",
    "    plt.xlabel('number of factors')\n",
    "    plt.ylabel('percent variance explained')\n",
    "    plt.legend()\n",
    "    plt.title('variance explained in $X$ and $Y$')\n",
    "    plt.show()\n",
    "    \n",
    "# The SIMPLS algorithm\n",
    "def SIMPLS(X, Y, numFactor):\n",
    "    X = np.matrix(X)\n",
    "    Y = np.matrix(Y)\n",
    "    \n",
    "    varXAry = []\n",
    "    varYAry = []\n",
    "    varSAry = []\n",
    "    \n",
    "    numSample = X.shape[0]\n",
    "    numHuman = X.shape[1]\n",
    "    numPathogen = Y.shape[1]\n",
    "    \n",
    "    Y0 = Y - getMean(Y) # Centralize Y\n",
    "    S = (X.T @ Y0) # Calculate S\n",
    "    for a in range(1, numFactor+1):\n",
    "        eigenvalues, eigenvectors = np.linalg.eig(S.T @ S)\n",
    "        dominant_index = np.argmax(np.abs(eigenvalues))\n",
    "        q = np.matrix(eigenvectors[:, dominant_index])\n",
    "        r = S @ q\n",
    "        t = X @ r\n",
    "        t = t - np.mean(t)\n",
    "        normt = np.sqrt(t.T @ t)\n",
    "        t = t/normt\n",
    "        r = r/np.sqrt(r.T @ r)\n",
    "        p = X.T @ t\n",
    "        q = Y0.T @ t\n",
    "        u = Y0 @ q\n",
    "        v = p\n",
    "        if a > 1:\n",
    "            v = v - V @ (V.T @ p)\n",
    "            u = u - T @ (T.T @ u)\n",
    "        v = v/np.sqrt(v.T @ v)\n",
    "        S = S - v @ (v.T @ S)\n",
    "        R, T, P, Q, U, V = r, t, p, q, u, v\n",
    "        if a == 1:\n",
    "            T_mat = T\n",
    "            R_mat = R\n",
    "        else:\n",
    "            T_mat = np.hstack([T_mat, T])\n",
    "            R_mat = np.hstack([R_mat, R])\n",
    "        varXAry.append(np.diag((P.T @ P)/(numSample - 1))[0])\n",
    "        varYAry.append(np.diag((Q.T @ Q)/(numSample - 1))[0])\n",
    "        # If the information is already extracted, stop the algorithm immediately\n",
    "#         if abs(np.sum(varYAry) - getVariance(Y)) < 1e-10:\n",
    "#             B = R @ Q.T\n",
    "#             h = np.diag(T @ T.T) + 1/numSample\n",
    "#             return B, tAry, varXAry, varYAry, varSAry\n",
    "    B = R_mat @ Q.T\n",
    "    h = np.diag(T_mat @ T_mat.T) + 1/numSample\n",
    "    return B, tAry, varXAry, varYAry, varSAry, rAry"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 313,
   "id": "97399d18",
   "metadata": {},
   "outputs": [],
   "source": [
    "h = Human(20, 100)\n",
    "v = Virus(h, 10)\n",
    "# v.Y[:,2] = v.X[:, 5]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 314,
   "id": "f5e82f5f",
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1, 0, 1, ..., 1, 0, 2],\n",
       "       [0, 2, 2, ..., 1, 2, 2],\n",
       "       [1, 2, 0, ..., 0, 0, 1],\n",
       "       ...,\n",
       "       [2, 1, 0, ..., 0, 1, 0],\n",
       "       [2, 2, 1, ..., 2, 2, 1],\n",
       "       [2, 1, 1, ..., 1, 0, 1]])"
      ]
     },
     "execution_count": 314,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "v.X"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 315,
   "id": "05adf66b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0, 0, 1, 1, 0, 0, 0, 0, 0, 1],\n",
       "       [1, 0, 1, 1, 1, 1, 1, 1, 0, 0],\n",
       "       [1, 0, 0, 0, 0, 0, 1, 0, 0, 1],\n",
       "       [1, 0, 0, 0, 0, 1, 1, 0, 1, 0],\n",
       "       [1, 0, 0, 0, 1, 0, 0, 1, 1, 0],\n",
       "       [1, 0, 0, 1, 1, 0, 0, 1, 0, 1],\n",
       "       [1, 0, 0, 0, 1, 1, 1, 1, 0, 1],\n",
       "       [0, 1, 1, 1, 1, 0, 0, 1, 0, 0],\n",
       "       [1, 0, 1, 0, 0, 0, 1, 0, 1, 0],\n",
       "       [1, 0, 0, 1, 0, 1, 1, 0, 0, 0],\n",
       "       [1, 0, 1, 1, 1, 1, 0, 1, 1, 0],\n",
       "       [1, 0, 1, 0, 0, 0, 0, 1, 1, 0],\n",
       "       [0, 0, 0, 1, 0, 1, 1, 1, 1, 0],\n",
       "       [1, 1, 1, 0, 0, 0, 1, 1, 1, 0],\n",
       "       [1, 1, 0, 1, 1, 1, 0, 0, 0, 0],\n",
       "       [0, 0, 0, 1, 1, 1, 1, 1, 0, 1],\n",
       "       [1, 1, 1, 1, 0, 1, 0, 1, 0, 1],\n",
       "       [1, 1, 1, 1, 1, 1, 1, 1, 0, 0],\n",
       "       [0, 1, 0, 1, 1, 0, 0, 0, 0, 1],\n",
       "       [0, 0, 0, 0, 1, 0, 0, 0, 1, 0]])"
      ]
     },
     "execution_count": 315,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "v.Y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 158,
   "id": "c1f6d23e",
   "metadata": {},
   "outputs": [],
   "source": [
    "varX1Total = getVariance(v.X)\n",
    "varY1Total = getVariance(v.Y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 159,
   "id": "68016d9e",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "matrix([[0.5       +0.j, 0.4794412 +0.j, 0.50509717+0.j, 0.24950829+0.j,\n",
       "         0.49950829+0.j],\n",
       "        [0.5       +0.j, 0.56330078+0.j, 0.48430577+0.j, 0.25151397+0.j,\n",
       "         0.50151397+0.j],\n",
       "        [0.5       +0.j, 0.54724313+0.j, 0.48828696+0.j, 0.25112992+0.j,\n",
       "         0.50112992+0.j],\n",
       "        [0.5       +0.j, 0.41001489+0.j, 0.52231011+0.j, 0.24784782+0.j,\n",
       "         0.49784782+0.j]])"
      ]
     },
     "execution_count": 159,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(v.X - getMean(v.X)) @ B + getMean(v.Y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 160,
   "id": "fb21d847",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[]"
      ]
     },
     "execution_count": 160,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "varS1Ary"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 161,
   "id": "58fbee4f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAHHCAYAAABDUnkqAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAA9hAAAPYQGoP6dpAAB+yElEQVR4nO3deVxU1fvA8c8wwLCvIoKKoCwC7lq5ZFaWpn0t27Ss1MxftpppWVZmlmWr2qqtmmVWZpaVpW2aaVYqprKjuCGK7PsyM+f3hzqJgjI4wzDD8369eL2cO+fe+9x7cebh3ueco1FKKYQQQgghHISTrQMQQgghhLAkSW6EEEII4VAkuRFCCCGEQ5HkRgghhBAORZIbIYQQQjgUSW6EEEII4VAkuRFCCCGEQ5HkRgghhBAORZIbIYQQQjgUSW6EEEII4VAkuRFCCCGEQ5HkRggLWrJkCRqNhn379tk6FLt1PufQ1uf/6aefRqPRnLWNtWJ8/vnn0Wg0JCQknPHe2LFjcXZ25vvvv7foPhujIecIYPjw4Xh4eJCbm3vGewcOHCA0NJSLLrqI8vJya4Qp7JwkN0JY0NVXX82ff/5JSEiIrUNpkezh/Fsrxvvvvx8/Pz+ef/75WsufeuopPv74Y9544w2uvvpqi+7TmmbOnElFRQVvvPFGreVFRUWmxOe7777Dw8PDRhGK5kySGyEspLy8nKCgIPr27YtOp7N1OC2SPZx/a8Xo4+PDAw88wFdffUVKSgoAixcv5tlnn2X69Oncc889Ft2ftfXr14/Bgwfz1ltvme7OVFdXc91113HkyBF++OEHgoKCbBylaK4kuREO7+uvv0aj0fDLL7+c8d7ChQvRaDTs3LmTjIwM7rjjDqKiovDw8KBt27aMGDGCXbt2nbHeyVvr27dv58Ybb8Tf359OnTrV+cihods9uc3ExERuueUWfH19CQ4OZsKECRQVFZ0RQ0pKCrfccgvBwcHodDrCwsIYO3YsVVVVpjbp6emMGTOG1q1bo9PpiI2N5a233mrwuTvX+pWVlfTs2ZPIyMhaMR45coQ2bdpw6aWXYjAYah1fQkIC119/PT4+Pvj6+nLbbbdx7Nixs8bR0HNY1/k357yac76+//57evTogU6nIyIigldeeaVB5/R8YzybKVOm4OHhwdy5c/n555+ZNGkSo0aN4oUXXjjnutb4PW3sOTpp5syZ5OXl8cEHH6CUYsKECfz555+sXr2aqKgos7YlWhZnWwcghLX973//o3Xr1ixevJjBgwfXem/JkiX06tWLbt268fvvvxMYGMgLL7xAUFAQ+fn5fPTRR1x00UUkJCQQExNzxravv/56br75Zu6++27KysooKCg4o83hw4fN2u4NN9zA6NGjufPOO9m1axczZswA4MMPPzS1+ffff7n44otp1aoVzzzzDFFRUWRnZ7N69Wqqq6vR6XQkJSXRv39/wsLCePXVV2nTpg1r165l8uTJ5ObmMmvWrLOet4as7+bmxhdffEHv3r2ZMGECK1euxGg0cuutt6KUYvny5Wi12lrbve666xg1ahR33303iYmJzJw5k6SkJP766y9cXFzqjMXcc1iXc51Xc87XL7/8wrXXXku/fv347LPPMBgMvPTSSxw9evSccZxPjOcSEBDAPffcw/z58/n666+56KKLWLp0aYNqXCz9e2qJczRo0CAGDhzIvHnzOHToEMuXL2fFihX079+/wdsQLZQSogWYOnWqcnd3V4WFhaZlSUlJClBvvPFGnevo9XpVXV2toqKi1EMPPVTrvVmzZilAPfXUU7WWL168WAEqMzOz3ljq2+7Jbb700ku12t97773Kzc1NGY1G07LLL79c+fn5qZycnHr3M3ToUNWuXTtVVFRUa/n999+v3NzcVH5+fr3rmrv+559/rgC1YMEC9dRTTyknJye1bt26WuudPL7Tz+WyZcsUoD755BOl1Pmdw7rWbeh5Ned4L7roIhUaGqoqKipMy4qLi1VAQIA618fq+cTYENu3b1eACg0NVXl5eQ1e73Tn+3t6PufoVGvXrlWAAtT8+fMbfTyiZZHHUqJFmDBhAhUVFXz++eemZYsXL0an0zFmzBgA9Ho9zz//PHFxcbi6uuLs7Iyrqyvp6ekkJyfXud0bbrjhnPs2d7vXXHNNrdfdunWjsrKSnJwc4Hhtz4YNGxg1alS9NQeVlZX88ssvXHfddXh4eKDX600/w4cPp7Kyki1bttQbs7nrjxo1invuuYdHHnmEOXPm8Pjjj3PllVfWue1bb7211utRo0bh7OzMb7/9Vm88jbk2pzvbeTXneMvKyvjnn3+4/vrrcXNzM23P29ubESNGNCiWxsTYEMXFxdxxxx0A5ObmmtWTyJK/p5Y8R0OGDMHDw4MrrriCKVOmmLWuaLkkuREtQnx8PBdccAGLFy8GwGAw8Mknn3DttdcSEBAAwNSpU5k5cyYjR47k22+/5a+//uKff/6he/fuVFRU1LndhvR4MXe7gYGBtV6fLDw92bagoACDwUC7du3q3WdeXh56vZ433ngDFxeXWj/Dhw8HqLOL7fmsP2HCBGpqanB2dmby5Mn1brtNmza1Xjs7OxMYGEheXl696zTm2pzubOfVnOMtKCjAaDSecRx1HZu5znXtz6ampoYbbriBvXv38vXXX2M0Gnn55ZcbvG9L/p5a8hzt3buX8vJy+vXrZ9Z6omWTmhvRYtxxxx3ce++9JCcns3fvXrKzs01/5QJ88sknjB079oyutLm5ufj5+dW5zYbUMjRmu2cTEBCAVqvl0KFD9bbx9/dHq9Vy++23c99999XZJiIiwmLrl5WVcfvttxMdHc3Ro0eZOHEi33zzTZ3rHTlyhLZt25pe6/V68vLyzviyPJWlz+HpzDlef39/NBoNR44cOaNNXcuaysSJE1m/fj3fffcdQ4cO5bbbbuO9997jiSeeoHXr1udc35Ln2JLnaPv27QD06tXLrPVEyyZ3bkSLccstt+Dm5saSJUtYsmQJbdu2ZciQIab3NRrNGd1zv//+e7Kyss5rv5berru7O4MGDWLFihX13n3x8PDgsssuIyEhgW7dutGnT58zfs6WTJi7/t13382BAwf46quv+OCDD1i9ejXz58+vc9vLli2r9fqLL75Ar9dz6aWX1huPta7NSeYcr6enJxdeeCFfffUVlZWVpm2UlJTw7bffWiQecz355JMsXbqURYsWMXToUAAef/xxqqur670Op7PkObbkOTqZ3PTu3dvsOETLJXduRIvh5+fHddddx5IlSygsLOThhx/Gyem//P5///sfS5YsoXPnznTr1o1t27bx8ssvn/XxT0NYY7vz5s3j4osv5qKLLuKxxx4jMjKSo0ePsnr1at555x28vb157bXXuPjiixk4cCD33HMP4eHhlJSUkJGRwbfffsuvv/561n00dP3333+fTz75hMWLFxMfH098fDz3338/jz76KAMGDODCCy+std2vvvoKZ2dnrrzySlNvqe7duzNq1KgmPYeNPV6AZ599lquuuoorr7ySadOmYTAYePHFF/H09CQ/P99iMTXEu+++y3PPPceTTz7JnXfeaVoeFRXFTTfdxNtvv82jjz56zrsvlj7HljpH27Zto1WrVrRv375RcYgWytYVzUI0pXXr1pl6XqSlpdV6r6CgQN15552qdevWysPDQ1188cVq48aNatCgQWrQoEG12p7sMXLs2LFay+vqCdPQ7ZqzTaWO9/a66aabVGBgoHJ1dVVhYWFq/PjxqrKy0tQmMzNTTZgwQbVt21a5uLiooKAg1b9/fzVnzpwGna9zrb9z507l7u6uxo0bV2u9yspK1bt3bxUeHq4KCgpqHd+2bdvUiBEjlJeXl/L29la33HKLOnr0qEXO4dl6IjXkvJpzvlavXq26detmOvcvvPCCaV9nc74xnur7779XWq1W3XbbbXW+v2vXLqXRaNTs2bPPGpNS1vk9bew5OlVQUJAaMmRIg9sLoZRSGqWUauqESgjR8jz99NPMnj2bY8eO0apVK1uHI4RwYFJzI4QQQgiHIsmNEEIIIRyKPJYSQgghhEOROzdCCCGEcCiS3AghhBDCoUhyI4QQQgiH0uIG8TMajRw+fBhvb+8GDZ0vhBBCCNtTSlFSUkJoaGitAVjr0uKSm8OHD8tIl0IIIYSdOnjw4DlHzm5xyY23tzdw/OT4+PjYOBohhBBCNERxcTHt27c3fY+fTYtLbk4+ivLx8ZHkRgghhLAzDSkpkYJiIYQQQjgUSW6EEEII4VAkuRFCCCGEQ2lxNTcNZTAYqKmpsXUYdsnV1fWc3fSEEEIIa5Hk5jRKKY4cOUJhYaGtQ7FbTk5ORERE4OrqautQhBBCtECS3JzmZGLTunVrPDw8ZKA/M50cJDE7O5uwsDA5f0IIIZqcJDenMBgMpsQmMDDQ1uHYraCgIA4fPoxer8fFxcXW4QghhGhhpDDiFCdrbDw8PGwciX07+TjKYDDYOBIhhBAtkSQ3dZBHKedHzp8QQghbkuRGCCGEEA5FkhshhBBCOBRJbkS9/v33X3Q6HatXr661fOXKlbi5ubF7924bRSaEEELUz6bJze+//86IESMIDQ1Fo9Hw9ddfN3jdTZs24ezsTI8ePawWX0tWXV1N9+7dmTlzJnfddRd5eXkA5OTkcPfddzN79my6dOli4yiFEEI0N4llZWSUl9s0BpsmN2VlZXTv3p0333zTrPWKiooYO3YsgwcPtlJk9uWdd96hbdu2GI3GWsuvueYaxo0bx549e7j22msJDg7Gy8uLCy64gJ9//rlW2/DwcObMmcP48ePx9fXl//7v/wCYMWMGYWFh3HfffQBMmjSJqKgoHn744aY5OCGEEHbj06NHuXDbNm5MTKTChj1mbTrOzbBhwxg2bJjZ602aNIkxY8ag1WrNutvTGEopyk9LGpqKh5NTg3oe3XTTTUyePJnffvvNlPAVFBSwdu1avv32W0pLSxk+fDhz5szBzc2Njz76iBEjRpCamkpYWJhpOy+//DIzZ87kySefNC3TarV89NFH9OrVizFjxrB27Vp27NiBVqu1/AELIYSwS1VGIw9lZLDw8GEAglxdqTAacbfRd4XdDeK3ePFi9uzZwyeffMKcOXPO2b6qqoqqqirT6+LiYrP2V2404rVxo9lxWkLpwIF4NuAXIyAggKuuuopPP/3UlNysWLGCgIAABg8ejFarpXv37qb2c+bMYdWqVaxevZr777/ftPzyyy+v845MbGwsU6ZM4YUXXuDFF18kOjraAkcnhBDCEeyrqOCmpCS2lpSgAWZ26MBT4eFobTgsiF0VFKenp/PYY4+xbNkynJ0blpfNnTsXX19f00/79u2tHKVt3HrrraxcudKUyC1btoybb74ZrVZLWVkZ06dPJy4uDj8/P7y8vEhJSeHAgQO1ttGnT586t11aWsrnn3+Oh4cHG22U6AkhhGh+vs/Lo9e2bWwtKSHA2Zk1XbsyOyLCpokN2NGdG4PBwJgxY5g9e7ZZdw5mzJjB1KlTTa+Li4vNSnA8nJwoHTjQrFgtxcOMmbVHjBiB0Wjk+++/54ILLmDjxo3MmzcPgEceeYS1a9fyyiuvEBkZibu7OzfeeCPV1dW1tuHp6Vnnth955BFcXV3ZvHkz/fr1Y+nSpYwdO7bxByaEEMKuGZRiVmYmz534I/lCb29WxMcT5uZm48iOs5vkpqSkhK1bt5KQkGB6lGI0GlFK4ezszLp167j88svPWE+n06HT6Rq9X41G06BHQ7bm7u7O9ddfz7Jly8jIyCA6OprevXsDsHHjRsaPH891110HHL8Ts2/fvgZt96effuL9999n48aNdO/eneeff54pU6Zw5ZVXEhISYq3DEUII0UzlVFdzS1ISvxYWAnB/27a80qkTulP+IDdUGtC62e67026SGx8fH3bt2lVr2dtvv82vv/7Kl19+SUREhI0iaz5uvfVWRowYQWJiIrfddptpeWRkJF999RUjRoxAo9Ewc+bMM3pW1aW4uJg777yThx9+mL59+wIwefJkVq5cyV133cW3335rtWMRQgjR/PxRWMjopCQOV1fj6eTE+zEx3BwcjLHGSOHGQvLW5JG/Jh+XVi70+K2HzeK0aXJTWlpKRkaG6XVmZiY7duwgICCAsLAwZsyYQVZWFkuXLsXJyemMcVVat26Nm5ubjLdywuWXX05AQACpqamMGTPGtHz+/PlMmDCB/v3706pVKx599NEGFVZPmTIFX19fZs+ebVrm5OTE4sWL6d69uzyeEkKIFkIpxfxDh5i+Zw8GINbDgy9aRRP4QwWJaxLJX5ePoei/rt8anQZDmQGtZwvsLbV161Yuu+wy0+uTtTHjxo1jyZIlZGdnn1H0Kuqn1Wo5fKIb3qnCw8P59ddfay07OW7NSXU9pvrwww/r3E9kZCRlZWWND1QIIYTdKNLruSMlhW+O5hKTChN3eXDx307kbt9B7intnAOdCbgqgMCrAwkYEmCzxAZAo5RSNtu7DRQXF+Pr60tRURE+Pj613qusrCQzM5OIiAjcmklRlD2S8yiEEI4hYX8hryxNJPz3Gi78G3xPu+nv1duLwOGBBAwPwOcCHzRa6/WSOtv39+nspuZGCCGEENallKJ0Ryn5a/JJ/jobl+2V/N8pJZpaHy0BQwMIGB5AwFUB6No0vsOONUlyI4QQQrRg+iI9BT8XHC8G/iGf6uzjw4ScvO9+LFJL7LVtaDciCJ/+Pji5NP8h8iS5EUIIIVoQpRTlyeXkfX+8Z1PRH0Uo/X8VKlVusLUX/H0RXHxje6b164iTjQflM5ckN0IIIYSDM5QZKPitgPw1+eStyaNqf1Wt992j3Sm4zJ250YX82cWIn6cLy+PiGOzvb6OIz48kN0IIIYQDKs8oNyUzhesLUVX/3Z3R6DT4XepH4PBAfK/y5xntEV45eBCAAT4+fB4fT9vzGADX1iS5EUIIIRyAscpI4e+FpoSmIq2i1vu6MN3xbtrDA/C/zB+tp5bsqiquTkpiY1ERAFPbteOFjh1xMWP6n+ZIkhshhBDCTlUeqCT/h+PJTMHPBRjL/+vapHHW4DvQl4DhAQQOD8Qj1gPNKbUz6wsKuDkpiaM1NXhrtSzu3JkbgoJscRgWJ8mNEEIIYSeMNUaKNxebpjko2117QFXXNq6mZMb/Sn+cfc78mjcqxUsHDvBEZiZGoKunJyvj44ny8Giio7A+SW6EEEKIZqzqSBX5P+aT/33+8WkOiv+b5gAn8OnrYxpIz6uHV627M6crqKlhbEoK3+XlATAuOJi3o6PxsIMJos0hyY2oU1paGj169OD999+vNU+V0Wjk4osvJjg4mFWrVtkwQiGEcEzKoCj+p9hUO1O6rbTW+86BzgQOO57MBAwJwCXQpUHb3VZSwo2JieyrrESn0fBWdDQT2rQ5azJkryS5EXUKDw/nhRde4IEHHuCyyy4jJCQEgFdffZWMjAy+/vpr2wYohBAOpCavhvy1x5OZ/B/z0efpa71/cpqDwKsD8e7jbdY0B0op3s3OZnJ6OtVK0dHNjS/j4+np7W3pw2g2JLlxAO+88w7PPPMMBw8exOmUCvdrrrkGf39/nnrqKaZOncqWLVsoKysjNjaWuXPncsUVV5jahoeHM3HiRDIyMli1ahUjR45kyZIlfPPNN/zf//0f3333HSkpKTz11FMsX76c1q1b2+JQhRDCIZw6zUHe93kU/1UMp05z4KslYMj5T3NQZjBwT1oaHx89CsC1gYEs6dwZP5eG3e2xVzJx5inqmvBRKVWr+rwpOXk4Neh2YX5+PiEhIaxZs4bBgwcDUFBQQJs2bfj2228JDg5my5Yt9O/fHzc3Nz766CNeffVVUlNTCQsLA44nNwUFBcycOZORI0cCx2f/PnDgAF27duWVV17hgw8+oHPnzixZsuSs8cjEmUIIcab6pjk4ybOrp6kY2Kff+U9zkFpezo2JiewuK0MLPN+xI4+0b2+3j6Fk4kwLMpYb2ei10Sb7Hlg6sEFTxgcEBHDVVVfx6aefmpKbFStWEBAQwODBg9FqtXTv3t3Ufs6cOaxatYrVq1dz//33m5ZffvnlPPzww7W2HRYWxoIFC5g4cSJt27Zl7dq1Fjo6IYRwbEopypPKTT2bTp/mwMnDCf8r/E3FwG7tLffH4IqcHCakplJqMNDG1ZXP4uIY5Odnse03d5LcOIhbb72Vu+66i7fffhudTseyZcu4+eab0Wq1lJWVMXv2bL777jsOHz6MXq+noqKCAwcO1NpGnz596tz2HXfcwcyZM5k8eTK+vr5NcThCCGGXGjLNwclkxu8SP5x0lh0sr9poZPqePbyWlQXAIF9fPouLo40djzbcGJLcnIOThxMDSwfabN8NNWLECIxGI99//z0XXHABGzduZN68eQA88sgjrF27lldeeYXIyEjc3d258cYbqa4+7Zaop2e923d2dsbZWX5dhBDidOea5sD/Mv/jtTPDAvCItN5YMgcrKxmVlMSW4mIAHgsL49nwcJztfLThxpBvq3PQaDQNejRka+7u7lx//fUsW7aMjIwMoqOj6d27NwAbN25k/PjxXHfddQCUlpayb98+G0YrhBD2y1BpoOj3ItPjpor006Y56KAz9Wzyu8wPrYf1v0N+ys9nTHIyuTU1+Gq1LI2N5ZpWray+3+ZKkhsHcuuttzJixAgSExO57bbbTMsjIyP56quvGDFiBBqNhpkzZ2I02qZIWggh7JFpmoPv8yj4xbxpDqzJqBRz9u/n6X37UEAvLy9WxMfT0d29SfbfXEly40Auv/xyAgICSE1NrTXw3vz585kwYQL9+/enVatWPProoxSfuG0phBDiTOec5iDklGkOrqh7mgNry62u5rbkZNYWFABwV0gIr0VG4uZgow03hnQFP4V0YbYMOY9CCHtUdaSK/B/yyV9zlmkOTsyq7dX97NMcWNtfxcXclJjIwaoq3J2cWBQdzdg2bWwWT1OQruBCCCHEOZxrmgOXVi4EXBVg9jQH1qSU4s2sLKbt2UONUkS5u7MyPp6uXl62Dq1ZkeRGCCFEi3GuaQ68+3ibHjeZO82BtZXo9dyVlsZnOTkA3BgUxAcxMfhIT9YzyBkRQgjhsJRSlCaUmmpn6pvmIPDqQAKuCsA12NV2wZ5FYlkZNyYmklJejrNGw8sdO/Jgu3Z2O9qwtUlyI4QQwqE09TQH1rbs6FHuSk2l3GikrasrX8TH018GVD0rSW7q0MJqrC1Ozp8Qoimdc5oDz1OmORhm2WkOrKnKaOShjAwWHj4MwBX+/iyLjaW1a/O8u9ScSHJzCpcTs6SWl5fj3sLHCDgfJ0c+1kp3RCGElRjKDBT8eso0Bwfqmebg6gD8Blp+mgNr21dRwU1JSWwtKUEDzOzQgafCw9HKY6gGkeTmFFqtFj8/P3JOFGt5eDTdQEyOwmg0cuzYMTw8PGS6BiGERTWXaQ6s7fu8PG5PTqZAryfA2ZllsbFcFRho67Dsinz7nKbNiXECTiY4wnxOTk6EhYVJYiiEOC8Nmubg6kAChzfdNAfWZFCKWZmZPHdiUuMLvb1ZER9PmIwXZjZJbk6j0WgICQmhdevW1NTU2Docu+Tq6opTC5yoTQhx/mrya8j5Iof8Nflnn+bg6kA8OjvO3fWj1dWMSUri18JCAO4LDeXVyEh08lnaKJLc1EOr1UrNiBBCNCF9kZ5tfbZRmVlpWtYcpjmwtj8KCxmdlMTh6mo8nZx4LyaGW4KDbR2WXXO83xIhhBB2KeOhDCozK3ENdaXtfW2bxTQH1qSUYt6hQzy6Zw8GINbDg5Xx8cR6eto6NLsnyY0QQgiby/02lyOLj4AG4j6Lw2+gn61DsqoivZ47UlJYlZsLwC2tW/NudDRe0hHDIuQsCiGEsKnq3GpS/y8VgPbT2jt8YvNvaSk3JiaSUVGBq0bDgshI7g4Nddg7VLYgyY0QQgibUUqRfk86NUdr8IjzIPzZcFuHZFVLsrO5Jz2dSqORMJ2OL+PjueAcM1wL80lyI4QQwmZyPsvh2JfH0DhriF0ai9bNMTtyVBgMPJCezgdHjgAwLCCAj2NjCXSx/UzjjkiSGyGEEDZRdbiK9PvSAejwZAe8e3vbOCLr2FNRwY2JiewoLUUDPBMezuMdOuAkj6GsRpIbIYQQTU4pRerEVPQFerx6exH2eJitQ7KKr48dY3xKCkUGA0EuLnwaG8sVAQG2DsvhSXIjhBCiyWW/n03+D/lodMcfRzX3mbnNpTcamZGZySsHDwLQ38eHL+LjaavT2TiylkGSGyGEEE2qIrOCPVP3ANDxuY54xjnWuC7ZVVWMTkpiY1ERAFPbteOFjh1xkdGGm4wkN0IIIZqMMipSxqdgKDXgO9CXdlPa2Toki/qtoIBbkpI4WlODt1bL4s6duSEoyNZhtTiS3AghhGgyh147RNHvRTh5OtF5SWc0WscoqjUqxYsHDvBkZiZGoKunJyvj44nysN/Zye2ZJDdCCCGaRFlyGXtn7AUg8tVI3Du62zgiyyioqWFsSgrf5eUBMC44mLejo/GQ+QltRpIbIYQQVmfUG0kZm4KqUvgP9SfkrhBbh2QR20pKuDExkX2Vleg0Gt6MiuLOkBAZbdjGJLkRQghhdQfmHqBkawnOfs50/qCz3X/5K6V4NzubyenpVCtFRzc3voyPp6e3Y47VY28kuRFCCGFVJdtL2P/MfgCi3oxC19a+u0OXGQzck5bGx0ePAnBNYCAfde6Mn4w23GxIciOEEMJqjFVGkscmo/SKVte3ovWY1rYO6byklpdzY2Iiu8vK0ALPd+zII+3b2/2dKEcjyY0QQgiryXwqk/LEclxauxC9KNquk4AvcnK4MzWVUoOBNq6ufBYXxyA/P1uHJeogyY0QQgirKNpUxMGXj4/QG/1ONK5BrjaOqHGqjUYe2bOH17OyABjk68tncXG0kdGGmy1JboQQQlicocxA8rhkUBA8NpigkfY5kN3BykpGJSWxpbgYgMfCwng2PBxnGW24WZPkRgghhMXtmb6Hyj2V6NrpiHwt0tbhNMpP+fmMSU4mt6YGX62WpbGxXNOqla3DEg1g09Tz999/Z8SIEYSGhqLRaPj666/P2v6rr77iyiuvJCgoCB8fH/r168fatWubJlghhBANkv9TPoffPgxAzOIYXPzsqxeRUSme2bePoTt3kltTQy8vL7b36SOJjR2xaXJTVlZG9+7defPNNxvU/vfff+fKK69kzZo1bNu2jcsuu4wRI0aQkJBg5UiFEEI0RE1hDakTUgEIvS+UgCsCbByReXKrqxm+cyez9u1DAf8XEsKmnj3p6O4Yoym3FBqllLJ1EAAajYZVq1YxcuRIs9aLj49n9OjRPPXUUw1qX1xcjK+vL0VFRfj4+DQiUiGEEPVJHpfM0aVHcY90p8+OPmg97WcKgr+Ki7kpMZGDVVW4OzmxMDqacW3a2DoscYI53992XXNjNBopKSkhIKD+vwyqqqqoqqoyvS4+URQmhBDCso59fYyjS4+CE3T+qLPdJDZKKd7MymLanj3UKEWUuzsr4+Pp6uVl69BEI9l1uferr75KWVkZo0aNqrfN3Llz8fX1Nf20b9++CSMUQoiWofpYNWl3pQHQ/pH2+Pb3tXFEDVOi13NLUhKTMzKoUYobg4LY2ru3JDZ2zm6Tm+XLl/P000/z+eef07p1/SNezpgxg6KiItPPwYMHmzBKIYRwfEop0u5Oo+ZYDZ5dPImYHWHrkBoksayMC7dv5/Njx3DWaJjfqRNfxMXh42zXDzUEdvpY6vPPP+fOO+9kxYoVXHHFFWdtq9Pp0MlAS0IIYTVHlx0l96tcNM4aOi/tjJOu+f/dvOzoUe5KTaXcaKStqytfxMfT39c+7jaJc7O75Gb58uVMmDCB5cuXc/XVV9s6HCGEaNEqD1WSfn86AB1mdcC7Z/OeFbvKaGRKRgaLDh/vqn6Fvz/LYmNp7WqfoyeLutk0uSktLSUjI8P0OjMzkx07dhAQEEBYWBgzZswgKyuLpUuXAscTm7Fjx/Laa6/Rt29fjhw5AoC7uzu+knELIUSTUkqROjEVQ5EB7wu8CXsszNYhndW+igpuSkpia0kJADM7dGBWeDhaO57vStTNpvcOt27dSs+ePenZsycAU6dOpWfPnqZu3dnZ2Rw4cMDU/p133kGv13PfffcREhJi+nnwwQdtEr8QQrRkh985TMHaApzcnI4/jnJuvo+jvs/Lo9e2bWwtKSHA2Zk1XbvyTESEJDYOqtmMc9NUZJwbIYQ4fxV7Kvin+z8Yy4x0mt+J9lOaZ09Ug1LMyszkuRN/KF/g7c2K+Hg6uLnZODJhrhYzzo0QQoimpwyKlPEpGMuM+A7ypd3kdrYOqU5Hq6sZk5TEr4WFANwXGsqrkZHoZNJLhyfJjRBCCLMcnH+Qoj+K0Hpp6by4Mxqn5vdo54/CQkYnJXG4uhpPJyfei4nhluBgW4clmogkN0IIIRqsLLGMzCcyAeg0vxPuEc1rziWlFPMOHeLRPXswALEeHqyMjyfW09PWoYkmJMmNEEKIBjHWGEkem4yqVgQMCyDkzhBbh1RLkV7PHSkprMrNBeCW1q15NzoaLxmUr8WRKy6EEKJB9j+3n9LtpTj7OxPzfgyaZtTT6N/SUm5MTCSjogIXjYYFkZHcExrarGIUTadByc3UqVMbvMF58+Y1OhghhBDNU8m2EvbP2Q9A1NtR6EKbz8jvi7OzuTc9nUqjkTCdjhXx8VwovWFbtAYlNwkJCbVeb9u2DYPBQExMDABpaWlotVp69+5t+QiFEELYlKHSQPLYZDBA0E1BtB5d/3x+TanCYOCB9HQ+ODGg67CAAD6OjSXQxcXGkQlba1By89tvv5n+PW/ePLy9vfnoo4/w9/cHoKCggDvuuIOBAwdaJ0ohhBA2k/lkJuVJ5bgEuxD1dlSzeNSTUV7OjYmJ/FtWhgZ4Jjycxzt0wKkZxCZsz+xB/Nq2bcu6deuIj4+vtXz37t0MGTKEwyfm62iuZBA/IYRouMKNhewYtAMUdFndhVYjWtk6JL4+doxxKSkUGwwEubjwaWwsVwQE2DosYWXmfH+bPZJRcXExR48ePWN5Tk4OJSfm6xBCCGH/9KV6UsangII2d7SxeWJTYzTyyJ49XJeYSLHBQH8fHxL69JHERpzB7N5S1113HXfccQevvvoqffv2BWDLli088sgjXH/99RYPUAghhG3seXgPlXsr0YXpiFwQadNYDldVcXNSEhuLigCY2q4dL3TsiIuMNizqYHZys2jRIh5++GFuu+02ampqjm/E2Zk777yTl19+2eIBCiGEaHr5a/PJficbgM6LO+PsY7uRQ34rKOCWpCSO1tTgrdWyuHNnbggKslk8ovlr9MSZZWVl7NmzB6UUkZGReNrJ6I9ScyOEEGdXU1DDP13+ofpwNW0faEvU61E2icOoFC8eOMCTmZkYga6ennwZH0+0h4dN4hG21SQTZ2ZnZ5Odnc0ll1yCu7s7SqlmUUEvhBDi/KQ/kE714Wrco93p+EJHm8RQUFPD2JQUvsvLA2BccDBvR0fjodXaJB5hX8x+WJmXl8fgwYOJjo5m+PDhZGcfv205ceJEpk2bZvEAhRBCNJ1jK4+RsywHnKDzR53RejR9MrG1uJhe27bxXV4eOo2G96KjWdy5syQ2osHMTm4eeughXFxcOHDgAB6n3BocPXo0P/74o0WDE0II0XSqj1aTdncaAGGPhuHb17fJY3jv8GEGJCSwr7KSjm5u/NmrFxNlGgVhJrMfS61bt461a9fSrl27WsujoqLYv3+/xQITQgjRdJRSpE5KpSa3Bs9unoTPCm/S/dcYjTyUkcFbJ8ZKuyYwkI86d8ZPRhsWjWB2clNWVlbrjs1Jubm56HTNZ64RIYQQDXf046PkfZOHxkVD7NJYnHRN18U6v6aGUYmJ/FJYCMCciAhmhIXJaMOi0cz+7b3kkktYunSp6bVGo8FoNPLyyy9z2WWXWTQ4IYQQ1ld5sJL0B9IBCH86HK/uXk227+SyMi7avp1fCgvxdHJiVXw8T8g0CuI8mX3n5uWXX+bSSy9l69atVFdXM336dBITE8nPz2fTpk3WiFEIIYSVKKMidUIqhmIDPn19aD+9fZPt+4e8PG5OSqLYYKCDTsfqrl3p5tV0iZVwXGbfuYmLi2Pnzp1ceOGFXHnllZSVlXH99deTkJBAp06drBGjEEIIKzm86DAFPxfg5O5E54864+Rs/cdRSilePXiQ/+3aRbHBwEBfX/7p3VsSG2ExjR7Ez17JIH5CCHFceUY5W7tvxVhuJPK1SNpNbnfulc5TldHI3WlpLDlyBICJISG8FRWFq0yjIM7B6oP4FRYW8vfff5OTk4PRaKz13tixYxuzSSGEEE1IGRQp41Iwlhvxu8yPtve3tfo+j1ZXc/3u3WwuLsYJmB8ZyQNt20o3b2FxZic33377LbfeeitlZWV4e3vX+qXUaDSS3AghhB04+OpBijcXo/XW0nlxZzRO1k0wEkpKuHb3bg5WVeGr1fJFfDxDZDZvYSVm3wecNm0aEyZMoKSkhMLCQgoKCkw/+fn51ohRCCGEBZXuKiVzZiYAkQsicevgZtX9fZmTw8UJCRysqiLa3Z2/eveWxEZYldl3brKyspg8eXKdY90IIYRo3ozVRlLGpqCqFYH/C6TNHW2sty+leHb/fp7etw+AIf7+fBYXh78MzCeszOw7N0OHDmXr1q3WiEUIIYSV7Z+zn9IdpTgHOBP9XrTV6l3KDAZGJyWZEpsp7drxfdeuktiIJmH2nZurr76aRx55hKSkJLp27YrLab+o11xzjcWCE0IIYTnF/xSz//nj0+REL4xG18Y6o8ofrKzk2t27SSgtxUWjYVF0NBNCQqyyLyHqYnZXcKezdNfTaDQYDIbzDsqapCu4EKIlMlQY2NZrG+Up5bS+uTVxy+Ossp8/i4q4bvdujtbUEOTiwlfx8Vzs52eVfYmWxapdwU/v+i2EEKL5y3wik/KUclzbuBL1VpRV9vHRkSPclZpKtVJ08/RkddeudHCzbrGyEHVp1Dg3Qggh7EfhhkIOLTgEQMwHMbgEWLbuxaAUj+3dyysHDwIwslUrPu7cGS9n+YoRttGg37zXX3+du+66Czc3N15//fWztp08ebJFAhNCCHH+9CV6UsangIKQiSEEDg+06PaL9XpuSUpizYmhQJ7s0IHZ4eEy8aWwqQbV3ERERLB161YCAwOJiIiof2MaDXv37rVogJYmNTdCiJYk9a5Ust/Lxi3cjT47++Dsbbm7KRnl5VyzezfJ5eW4OTmxOCaGm4ODLbZ9IU5l8ZqbzMzMOv8thBCi+cpbk0f2e9kAxCyOsWhi82tBATclJpKv1xPq6so3XbrQR/5gFM2EPBAVQggHVJNfQ+rEVADaTWmH/6X+Ftv221lZTE5PxwBc6O3Nqi5dCNVZp1u5EI3RqOTm0KFDrF69mgMHDlBdXV3rvXnz5lkkMCGEEI2Xfn861dnVuMe4E/F8/eUE5qgxGnkwI4OFhw8DcGvr1rwXE4O7VmuR7QthKWYnN7/88gvXXHMNERERpKam0qVLF/bt24dSil69elkjRiGEEGbIWZFDzvIc0ELs0li07ueffOTV1HBTYiK/FRaiAZ6PiODRsDCZ0Vs0S2ZPvzBjxgymTZvG7t27cXNzY+XKlRw8eJBBgwZx0003WSNGIYQQDVR1pIq0e9IA6DCjAz4Xnn8dTFJZGRdu28ZvhYV4abV806ULj3XoIImNaLbMTm6Sk5MZN24cAM7OzlRUVODl5cUzzzzDiy++aPEAhRBCNIxSirS70tDn6fHq4UWHmR3Oe5vf5+XRd/t29lZWEu7mxp89ezKiVSsLRCuE9Zid3Hh6elJVVQVAaGgoe/bsMb2Xm5truciEEEKY5ciSI+R9m4fGVUPnpZ1xcjX7I95EKcXLBw4wYtcuSgwGBvn68k+vXnTx8rJgxEJYh9k1N3379mXTpk3ExcVx9dVXM23aNHbt2sVXX31F3759rRGjEEKIc6jcX0nGgxkARDwTgVfXxichlQYDd6Wl8fHRowDcFRLCG1FRuJ5lbkEhmhOzk5t58+ZRWloKwNNPP01paSmff/45kZGRzJ8/3+IBCiGEODtlVKRMSMFQYsCnnw/tH27f6G0dqariusREthQXowUWREZyX9u2Ul8j7IrZyU3Hjh1N//bw8ODtt9+2aEBCCCHMk/VWFoW/FuLk4UTnjzqj0TYuEdleUsK1u3dzqKoKf2dnVsTHM9jfcuPjCNFUZBA/IYSwY+Vp5ex99Pi0N51e6oRHlEejtvNFTg7jU1KoMBrp7OHB6i5diPJo3LaEsLUGJTf+/v4NviWZf2LyNCGEENZl1BtJGZeCscKI/xX+hN4Tav42lOLpfft4dv9+AK4KCOCzuDh8ZUZvYcca9Nu7YMECK4chhBDCXAdfPkjxlmK0PlpiPoxB42Te46gyg4Gxycl8daKn67R27XixUye0Ul8j7FyDkpuT49oIIYRoHkp3lrJv1j4Aol6Pwq29m1nr76+s5Npdu/i3rAxXjYZ3oqMZHxJihUiFaHqNuu9oMBhYtWoVycnJaDQaYmNjufbaa3GW25hCCGF1xmojyWOTUTWKwGsDCR4bbNb6m4qKuH73bnJqamjt4sKqLl3o7+trpWiFaHpmZyO7d+/m2muv5ciRI8TExACQlpZGUFAQq1evpmvXrhYPUgghxH/2zd5H2b9luLRyIeadGLO6aS/OzmZSWho1StHd05PVXbsS5mbeXR8hmjuzR2SaOHEi8fHxHDp0iO3bt7N9+3YOHjxIt27duOuuu6wRoxBCiBOK/yrmwAsHAIheFI1rsGuD1jMoxbSMDCakplKjFDe0asWmXr0ksREOyezk5t9//2Xu3Ln4nzL2gb+/P8899xw7duwwa1u///47I0aMIDQ0FI1Gw9dff33OdTZs2EDv3r1xc3OjY8eOLFq0yMwjEEII+2QoN5A8NhmM0PrW1gTdENSg9Yr0ev63axfzDh0C4KkOHfgiPh5P7fnPFi5Ec2R2chMTE8PRE0NynyonJ4fIyEiztlVWVkb37t158803G9Q+MzOT4cOHM3DgQBISEnj88ceZPHkyK1euNGu/Qghhj/bO2EtFWgWuoa5EvRHVoHXSy8vpu307P+bn4+7kxBdxccyOiMBJekQJB2Z2zc3zzz/P5MmTefrpp01zSW3ZssU0K3hxcbGprY+Pz1m3NWzYMIYNG9bgfS9atIiwsDBT1/TY2Fi2bt3KK6+8wg033GDWcZRVl6GtPvOvFq2TFjdnt1rt6uOkccLdxb1RbctrylFK1dlWo9Hg4eLRqLYVNRUYlbHeODxdPRvVtlJficFosEhbDxcPU41Alb4KvVFvkbbuLu44aY7n69WGamoMNRZp6+bshtZJa3bbGkMN1YbqetvqnHU4Ozmb3VZv1FOlr6q3ravWFReti9ltDUYDlfrKetu6aF1w1bqa3daojFTUVFikrbOTMzpnHXB8YsfymnKLtDXn/72tPiMKfi8gY2EGuECndzvh4u9Sb9uTfi3IZ2xKKoXKmXY6Hd906UKsm/NZ45DPCPPbymfEcU3xGdFQZic3//vf/wAYNWqU6Rfq5H+qESNGmF5rNBoMhvp/eRvjzz//ZMiQIbWWDR06lA8++ICamhpcXFzOWKeqqso0izlgSr5CXw2FOh41D48azvdjvje9bv1K63o/FAd1GMT68etNr8NfCye3vO6Z0fuE9uGf//vH9DrurTj2F+2vs21cUByJ9yaaXl/w3gUkHUuqs20H3w7sm7LP9PqSJZew9fDWOtu28mjFsUeOmV4PWzaMDfs31NnWw8WDssf/+wC84YsbWJO+ps62AGrWfx+st6+6nS+Tvqy3bemMUtMH3aTvJvHRvx/V2zbn4RyCPI/fep+6dipvb61/uo/MBzMJ9wsH4IlfnuCVP1+pt+3ue3YT3zoegOc3Ps/sDbPrbfv3xL+5oO0FALy25TWm/zy93ra/jfuNS8MvBeDdbe9y/w/319v2u1u+4+roqwFYtmsZd3xzR71tv7jxC26KvwmAVcmrGPXlqHrbLr52MeN7jAdgbcZa/rf8f/W2fXPYm9x34X0AbDywkcs+uqzeti9d8RKPDHgEgO3Z27nw/QvrbTtr0CyevvRpAJKPJdNlYZd62z7c72FeHvIyAAeKDhDxWkS9be/tcy9vXf0WALnlubR+pXW9bcd1H8eSkUuA41/+XnPrn0jyxrgbWXHTCtPrs7W16WfEEyfe2x9HIg37jEAXTN8ha1gVH08bnY4L3rtAPiOQzwh7/YxoKLOTm99++83cVSzmyJEjBAfX7vIYHByMXq8nNzeXkDrGaJg7dy6zZ9f/SymEEI7MU6vlt+7dcZP6GtGCaFR9zzuamEajYdWqVYwcObLeNtHR0dxxxx3MmDHDtGzTpk1cfPHFZGdn06ZNmzPWqevOTfv27Tl87HCdj83ksVTdbeWWs9xylsdS5re1xGdE7g+5JN2UBBrotrYbfv396v2MyK2u5rbkZP4oLgJgTngEU9q3b/T/e/mMkM+I5vQZUVxcjK+vL0VFRecsezH7zs3MmTN5+umn0Z72V0BRURF33303y5cvN3eTDdamTRuOHDlSa1lOTg7Ozs4EBgbWuY5Op0On052x3NPVs9Z/tvo0pE1j2p76wWTJtqd+OFqy7akf5pZsq3PWoePM63O+bV21rg1+Rmutti5aF9OHgiXbOjs54+zasP+65rTVOmkb/DtsTlsnjZNV2mo0Gqu0Bev9vzf3M6Imr4asu7Nwr3Gn3bR2tL20bb1td5eWcs3uVDIrq/F29eLT2Fj+16rVGW3lM+I4+Ywwv21z+IxoKLN7Sy1dupQBAwawZ88e07L169fTtWtX9u3bZ8nYztCvXz9++umnWsvWrVtHnz596qy3EUIIe5Z2bxrVR6rxiPUgYk79tUirc3Ppl5BAZmUlHd3c+LNXrzoTGyFaCrOTm507dxIeHk6PHj147733eOSRRxgyZAjjx4/njz/+MGtbpaWl7NixwzQ+TmZmJjt27ODAgeMDVM2YMYOxY8ea2t99993s37+fqVOnkpyczIcffsgHH3zAww8/bO5hCCFEs5bzeQ7HvjgGWui8tDNatzNrZpRSvLB/PyN376bUYOAyPz/+7t2beE/L/hUshL0x+7GUr68vn332GU888QSTJk3C2dmZH374gcGDB5u9861bt3LZZf9VXk+dOhU4PlHnkiVLyM7ONiU6ABEREaxZs4aHHnqIt956i9DQUF5//XWzu4ELIURzVpVdRdq9aQB0eLIDPn3OrC+oNBiYmJrKspwcAO4NDWVBZCQuTmb/zSqEw2lUQfEbb7zBo48+ynXXXce2bdvQarV8+umndO/e3RoxWpQ5BUlCCNHUlFLs+t8u8tfk49XLi15beuHkUjthya6qYuTu3fxdUoIWeCMqinva1l2PI4SjMOf72+wUf9iwYcyePZulS5eybNkyEhISuOSSS+jbty8vvfRSo4MWQggBRz48Qv6afDQ6DbFLY89IbLYWF9Nn2zb+LikhwNmZn7p3l8RGiNOYndzo9Xp27tzJjTfeCIC7uzsLFy7kyy+/ZP78+RYPUAghWoqKfRVkTMkAIGJOBJ7xtWtnPjt6lIE7dnC4uppYDw/+7t2by06Z508IcZxFx7nJzc2lVTOv0JfHUkKI5kgZFTsu30HRhiJ8L/alx/oeaLTHx28xKsVTmZk8d6IG8eqAAD6Ni8PH2eyySSHsllUfSwFs3LiR2267jX79+pGVlQXAxx9/TEpKSmM2J4QQLV7WG1kUbSjCydOJzks6mxKbUr2eGxITTYnNI+3b803XrpLYCHEWZic3K1euZOjQobi7u5OQkGAa/bekpITnn3/e4gEKIYSjK0spY+9jewHo9Eon3DsdHzxvX0UFAxIS+Do3F1eNho86d+alTp3QyozeQpyV2cnNnDlzWLRoEe+9916tgfP69+/P9u3bLRqcEEI4OqPeSMq4FIyVRvyH+BM6KRSAjYWFXLB9OzvLygh2cWF9jx6MrWOKGSHEmcy+r5mamsoll1xyxnIfHx8KCwstEZMQQrQYB188SMnfJWh9tcR8EINGo+GD7GzuSUujRil6ennxTZcutHdr+JQFQrR0Zt+5CQkJISMj44zlf/zxBx07drRIUEII0RKU7Chh3+x9AES9GYVzqCtT0tOZmJpKjVLcFBTExp49JbERwkxmJzeTJk3iwQcf5K+//kKj0XD48GGWLVvGww8/zL333muNGIUQwuEYq4yk3J6CqlG0ur4VulEBXL1rF6+d6KQxOzycz+Pi8NSeOe2CEOLszH4sNX36dIqKirjsssuorKzkkksuQafT8fDDD3P//fdbI0YhhHA4+57eR9nuMlyCXNDOa0/fhATSKirwcHJiaWwsNwQF2TpEIexWo8e5KS8vJykpCaPRSFxcHF5eXpaOzSpknBshhK0VbS4iYWACGKFmaQdGRWRRqNfTXqdjdZcu9PD2tnWIQjQ75nx/N3qgBA8PD/r06dPY1YUQokUylBlIGZcCRsi/0Zub2u/HqId+Pj6s6tKFYFdXW4cohN2TUaCEEKIJ7X1sLxUZFZS1cWLsuBKMwPg2bVgUHY1OZvQWwiIkuRFCiCZS8EsBWW8eLxieNc1IhRe82qkTD7Vrh0YG5hPCYiS5EUKIJqAv0rNrfDIAX18L6Rdp+S4ujmGBgTaOTAjHI8mNEEI0gV/v3o3roWqyQuGnKW5s6dWVWE/Pc68ohDBbox7wfvzxxwwYMIDQ0FD2798PwIIFC/jmm28sGpwQQtg7pRQL30/C9bNCjBr46TkvNg7oLYmNEFZkdnKzcOFCpk6dyvDhwyksLMRgMADg5+fHggULLB2fEELYrQqDgTs3JRI8PQeA9IlevDe2FwGnzMsnhLA8s5ObN954g/fee48nnngC7SkjZ/bp04ddu3ZZNDghhLBXWVVVXJKQQIcncgkogKoYVya+3hMX6RElhNWZ/b8sMzOTnj17nrFcp9NRVlZmkaCEEMKe/V1czAXbtuH7dSmDfgecof/yrmjdZCoFIZqC2clNREQEO3bsOGP5Dz/8QFxcnCViEkIIu/Xp0aNckpBATVY1U18/viz8qXC8e8qow0I0FbN7Sz3yyCPcd999VFZWopTi77//Zvny5cydO5f333/fGjEKIUSzZ1SKJzIzeeHAAVCw8DUXPEpq8O7jTdhjYbYOT4gWxezk5o477kCv1zN9+nTKy8sZM2YMbdu25bXXXuPmm2+2RoxCCNGslej13JaczOq8PADe3BxAxKZ8NDoNnZd2xslF6myEaEqNnjgTIDc3F6PRSOvWrS0Zk1XJxJlCCEvKrKjgmt272V1Whk6jYbFbBG0H78NYZqTTvE60f6i9rUMUwiFYdeLMzMxM9Ho9UVFRtGrVyrQ8PT0dFxcXwsPDzQ5YCCHs0YbCQm7YvZs8vZ42rq6siotHd81eisqM+F7iS7sH29k6RCFaJLPvlY4fP57Nmzefsfyvv/5i/PjxlohJCCGavXcPH+aKf/8lT6+nt5cX//TqRdsPiynaWITWS0vnJZ3ROMl8UULYgtnJTUJCAgMGDDhjed++fevsRSWEEI5EbzQyOT2dSWlp6JVidFAQv/fsif9eA3sf3wtAp3mdcI9wt3GkQrRcZj+W0mg0lJSUnLG8qKjINFqxEEI4ooKaGkYlJfFzQQEAz4aH80SHDii9ImFcCqpKETAsgJCJITaOVIiWzew7NwMHDmTu3Lm1EhmDwcDcuXO5+OKLLRqcEEI0FyllZVy0fTs/FxTg6eTEV/HxPBkejkaj4cDcA5RsLcHZ35mY92PQaORxlBC2ZPadm5deeolLLrmEmJgYBg4cCMDGjRspLi7m119/tXiAQghhaz/m5XFzUhJFBgMddDpWd+1KNy8vAEq2lbD/2eMTCEe9FYUuVGfLUIUQNOLOTVxcHDt37mTUqFHk5ORQUlLC2LFjSUlJoUuXLtaIUQghbEIpxfyDB7l61y6KDAYu9vXl7969TYmNodJA8rhklF4RdGMQrW+2n2ExhHBkZt+5AQgNDeX555+3dCxCCNFsVBmN3JOWxuIjRwC4s00b3o6OxvWUiS/3PbWP8sRyXIJdiFoYJY+jhGgmGpXcFBYW8vfff5OTk4PRaKz13tixYy0SmBBC2EpOdTXX797NpuJinIB5kZFMbtu2VvJS+EchB185CEDMuzG4tnK1UbRCiNOZndx8++233HrrrZSVleHt7V3rP7tGo5HkRghh1/4tLeWaXbs4UFWFr1bL5/HxDA0IqNVGX6onZVwKKGgzvg2trmlVz9aEELZgds3NtGnTmDBhAiUlJRQWFlJQUGD6yc/Pt0aMQgjRJFYdO0b/7ds5UFVFlLs7f/XufUZiA7B3+l4q91aiC9MRuSDSBpEKIc7G7Ds3WVlZTJ48GQ8PD2vEI4QQTU4pxZz9+3lq3z4ArvT35/O4OPxdXM5om78un8MLDwPQ+cPOOPs26um+EMKKzL5zM3ToULZu3WqNWIQQosmVGwzckpRkSmwebNuWNV271pnY1BTWkDIhBYC297fFf7B/U4YqhGggs//kuPrqq3nkkUdISkqia9euuJz2AXDNNddYLDghhLCmQ5WVjNy9m22lpbhoNLwdFcXE0NB622dMzqA6qxr3KHc6vtixCSMVQphDo5RS5qzg5FT/zR6NRtPsp2AwZ8p0IYTj+qu4mJG7d3OkuppWLi6sjI/nEj+/etsfW3WMxOsTwQl6/tET336+TResEMKs72+z79yc3vVbCCHszSdHjjAxNZUqpejq6cnqLl0Id69/osvqnGrSJqUBEDY9TBIbIZo5qYQTQrQYBqV4fO9eXjp4fHyaawMD+Tg2Fm/n+j8KlVKk3Z1GzbEaPLt5Ev50eBNFK4RorEYlN2VlZWzYsIEDBw5QXV1d673JkydbJDAhhLCkYr2eMUlJfH9iyIrHw8J4NiICp3OMKnz0k6PkrspF46IhdmksTjqz+2EIIZqY2clNQkICw4cPp7y8nLKyMgICAsjNzcXDw4PWrVtLciOEaHb2VlQwYtcuksrLcXNy4oOYGMYEB59zvcpDlaQ/kA5A+KxwvLp7WTtUIYQFmP0nyEMPPcSIESPIz8/H3d2dLVu2sH//fnr37s0rr7xijRiFEKLRfiso4IJt20gqLyfU1ZXfe/RoUGKjlCL1zlQMRQa8L/Km/aPtmyBaIYQlmJ3c7Nixg2nTpqHVatFqtVRVVdG+fXteeuklHn/8cWvEKIQQjbIoK4shO3eSr9dzgbc3//TuzQUN7CV5eNFhCtYV4OTmROxHsTg5y+MoIeyF2f9bXVxcTPNJBQcHc+DAAQB8fX1N/xZCCFuqMRq5Ly2Ne9LT0SvFmNat2dCjB6E6XYPWr9hTwZ6H9wDQ8YWOeMTIiOxC2BOza2569uzJ1q1biY6O5rLLLuOpp54iNzeXjz/+mK5du1ojRiGEaLD8mhpuSkzk18JCNMDzERE8GhZWa5Lfs1EGRfK4ZIzlRvwu86PtA22tG7AQwuLMvnPz/PPPExISAsCzzz5LYGAg99xzDzk5Obz77rsWD1AIIRrq39JSLty2jV8LC/HSavm6Sxce69ChwYkNwMF5ByneVIzWW0vMhzFonBq+rhCieTB7hGJ7JyMUC+F4ygwGntm3j1cPHsQAhLu5sbpLF7p6mde7qSyxjK29tqKqFTHvxxByZ4h1AhZCmM2c72+bV8i9/fbbRERE4ObmRu/evdm4ceNZ2y9btozu3bvj4eFBSEgId9xxB3l5eU0UrRCiuVmTl0eXf/7hpROJzfWtWvF3r15mJzbGGiPJY5NR1YqAqwNoM6GNdQIWQlhdg2puevXqxS+//IK/vz89e/Y86y3e7du3N3jnn3/+OVOmTOHtt99mwIABvPPOOwwbNoykpCTCwsLOaP/HH38wduxY5s+fz4gRI8jKyuLuu+9m4sSJrFq1qsH7FULYv8NVVUzJyGDFsWMAhOl0vBkVxYhWrRq1vf1z9lO6vRTnAGdi3osx61GWEKJ5aVByc+2116I70ctg5MiRFtv5vHnzuPPOO5k4cSIACxYsYO3atSxcuJC5c+ee0X7Lli2Eh4ebBgqMiIhg0qRJvPTSSxaLSQjRvBmU4p3Dh5mxdy/FBgNa4KH27ZnVoQNeZ5lG4WyKtxaz/7n9AES/HY0upGG9qoQQzVODPglmzZoFgMFg4NJLL6Vbt274+/uf146rq6vZtm0bjz32WK3lQ4YMYfPmzXWu079/f5544gnWrFnDsGHDyMnJ4csvv+Tqq68+r1iEEPbh39JSJqWm8ldJCQAXenvzTnQ0Pby9G71NQ4WBlLEpYICg0UG0Ht3aUuEKIWzErJobrVbL0KFDKSwsPO8d5+bmYjAYCD5tpNDg4GCOHDlS5zr9+/dn2bJljB49GldXV9q0aYOfnx9vvPFGvfupqqqiuLi41o8Qwr6UGQxM37OH3lu38ldJCd5aLW9GRbG5V6/zSmwAMp/MpDy5HNc2rkS/FW2hiIUQtmR2QXHXrl3Zu3evxQI4/bm2UqreZ91JSUlMnjyZp556im3btvHjjz+SmZnJ3XffXe/2586di6+vr+mnfXsZQl0Ie7ImL4/4v//m5RMFwzcGBZFy4YXc17Yt2vOsiyn8vZBD8w8BEPN+DC6BLhaIWAhha2Z3BV+3bh2PPvoozz77LL1798bT07PW+w3tXl1dXY2HhwcrVqzguuuuMy1/8MEH2bFjBxs2bDhjndtvv53KykpWrFhhWvbHH38wcOBADh8+bBp/51RVVVVUVVWZXhcXF9O+fXvpCi5EM3d6wXAHnY63oqO5OjDQItvXl+jZ2n0rlZmVtLmzDZ3f72yR7QohrMOcruBmV99dddVVAFxzzTW17rCcvONiMBgatB1XV1d69+7NTz/9VCu5+emnn7j22mvrXKe8vBzn0woGtVqtaf910el0pmJoIUTzZ1CKRYcP8/gpBcNT27dnVng4nif+v1vCnof3UJlZia6Djsh5kRbbrhDC9sxObn777TeL7Xzq1Kncfvvt9OnTh379+vHuu+9y4MAB02OmGTNmkJWVxdKlSwEYMWIE//d//8fChQsZOnQo2dnZTJkyhQsvvJDQ0FCLxSWEsI0dJSVMSkvj7xMFwxd5e/NOTAzdzRyz5lzyfswj+91sADov7oyzT+N6WQkhmiez/0cPGjTIYjsfPXo0eXl5PPPMM2RnZ9OlSxfWrFlDhw4dAMjOzq41Gef48eMpKSnhzTffZNq0afj5+XH55Zfz4osvWiwmIUTTKzMYeHrfPuafqKvx0WqZ27Ejk0JDz7uu5nQ1BTWk3pkKQNsH2+J/2fn1/BRCND+Nnn6hvLycAwcOUF1dXWt5t27dLBKYtcj0C0I0L9/l5nJfejoHTtTGjQoKYn5kZINn8DZX0q1J5Hyag3uMO30S+qB1t9yjLiGE9Vi15ubYsWPccccd/PDDD3W+39CaGyFEy5ZVVcWD6emszM0FjhcMvx0dzXALFQzXJefLHHI+zQEniP0oVhIbIRyU2V3Bp0yZQkFBAVu2bMHd3Z0ff/yRjz76iKioKFavXm2NGIUQDsSgFG8eOkTs33+zMjcXLfBI+/YkXnihVROb6qPVpN2dBkDYjDB8LpI7t0I4KrPv3Pz666988803XHDBBTg5OdGhQweuvPJKfHx8mDt3rowWLISo146SEu5KS+MfKxcMn04pRepdqejz9Hh29yT8qXCr7k8IYVtmJzdlZWW0bn18ePKAgACOHTtGdHQ0Xbt2NWvSTCFEy1Gq1/P0vn0sOHTIVDD8QseO3GWFguG6HF16lLzVeWhcNMQujcXJ1eyb1kIIO2J2chMTE0Nqairh4eH06NGDd955h/DwcBYtWlTnIHpCiJbt2xMFwwdPKRheEBlJSBONP1V5oJL0yekAhD8Tjlc3694lEkLYntnJzZQpU8jOPj4+xKxZsxg6dCjLli3D1dWVJUuWWDo+IYSdyqqqYnJ6Ol+dKBgOd3Pj7agohlmxruZ0yqhImZCCodiAT18fwh4Ja7J9CyFsp9FdwU8qLy8nJSWFsLAwWrVqZam4rEa6ggthXQaleDsriycyMyk5McLww+3b81R4OB4WHGG4IbLeyiL9/nSc3J3o828fPKI8mnT/QgjLsWpX8A0bNtQayM/Dw4NevXqZH6UQwuEknCgY3nqiYLivjw/vREfTzcoFw3UpTy9nzyN7AOj4UkdJbIRoQcyuqrvyyisJCwvjscceY/fu3daISQhhZ0r1eqZmZNBn2za2lpTgq9WyMCqKTT172iSxUQZFyrgUjBVG/Ab70fbetk0egxDCdsxObg4fPsz06dPZuHEj3bp1o1u3brz00kscOnTIGvEJIZq51bm5xP3zD/MPHcIIjA4KIvnCC7m7bVucmqAnVF0OvnKQ4j+L0fpo6fxhZzROtolDCGEb51Vzk5mZyaeffsry5ctJSUnhkksu4ddff7VkfBYnNTdCWMahykomZ2Sw6kTBcMSJguGrmrBguC6lu0rZ1mcbqloRsziGkPHSi1MIR2DVmptTRURE8Nhjj9G9e3dmzpzJhg0bzmdzQgg7YFCKt04UDJcaDDhrNDzcvj0zO3Ro8oLh0xmrjaSMTUFVKwKvCaTNuDY2jUcIYRuNTm42bdrEsmXL+PLLL6msrOSaa67h+eeft2RsQohmZntJCZNOKRjud6JguKsN6mrgeFfv8pRyijYVHf/ZWETl3kqcA52JficajY0eiwkhbMvs5Obxxx9n+fLlHD58mCuuuIIFCxYwcuRIPDykJ4IQjqpUr+epfft47URdja9Wy4udOvF/ISFNWldjKDdQ8k+JKZkp/rMYfYG+VhsnNyc6f9gZXZumGSRQCNH8mJ3crF+/nocffpjRo0fbxbg2Qojz801uLg+cMsLwza1bM79TJ9o0wQjDVUeqKN5UbEpmSreXovS1ywSd3J3wvtAb3wG++A7wxaefDy7+LlaPTQjRfJmd3GzevNkacQghmpmDJwqGv26igmFlVJQllh2/I7P5eEJTubfyjHauIa7Hk5gBPvgO8MWrhxdOLjJXlBDiP+dVUCyEcDwGpXgzK4snTykYfqR9e560cMGwocxA8d/Hk5jiTcUU/VmEochQu5EGPLt64tv/v2TGLdxNammEEGclyY0QwmRbSQmTUlPZVloKQP8TBcNdLFAwXJVVdfzx0ubjyUxJQgmclss4eTrhc5HPf4+Y+vrg7CsfU0II88inhhCCEr2emZmZvJGVhRHwc3bmxY4dmdjIgmFlUJTtLvuvF9OmIqr2V53RzrWtqymR8R3gi2d3T5yc5RGTEOL8SHIjRAv39bFjPJCRwaETBcO3tG7NPDMLhvWlekr+OqUX05ZiDMWn35YBr25epsdLvv190YXp5BGTEMLizE5uOnbsyD///EPgaUWFhYWF9OrVi71791osOCGE9RysrOSB9HS+ycsDoKObGwujoxkSEHDOdSsPVv5XK7OpiNJ/S8FYu43WS4tPXx9TMuNzkQ/OPvL3lBDC+sz+pNm3bx8Gg+GM5VVVVWRlZVkkKCGE9eiNRlPBcJnRiLNGw/QTBcPudRQMG/VGynaV1Upmqg6e+YhJF6b7r1ZmgA+eXeQRkxDCNhqc3Kxevdr077Vr1+Lr62t6bTAY+OWXXwgPD7docEIIy9paXMyktDS2nygYHuDjwzsxMcR7epra6Iv1FG/5b2yZkr9KMJSe9geNFry6e9VKZtzauTXloQghRL0anNyMHDkSAI1Gw7hx42q95+LiQnh4OK+++qpFgxNCWEbxiYLhN08pGH6pY0cmtGlDzcFqjn5z1JTMlO0qO/MRk48Wn37/9WLyvtAbZy95xCSEaJ4a/OlkNB7/tIuIiOCff/6R0YmFsANKKb4+McJwVnU1TgZ4MM+fSQd9Mb5ewF+b91GdVX3Gem7hbv8V/g7wxTPeE41WCn+FEPbB7D+9MjMzrRGHEMLCDlRWMm17Kvv/KGDIbrggyYnOyQpNeQFHKTC10zhr8Orp9d+ov/190YXKvExCCPvVqPvKv/zyC7/88gs5OTmmOzonffjhhxYJTAhhHqUUlZmVFGwsZPNP2ZT8Wcw9meBkmorp+P9VZz9n0yMmnwE++Fzgg9bTciMPCyGErZmd3MyePZtnnnmGPn36EBISImNUCGEjxhojpQml/40ts6mY6iPHHzG1OvED4NRRR9DFfqYpDDzjPNE4yf9bIYTjMju5WbRoEUuWLOH222+3RjxCiHrUFNSYJpQs2lREyT8lGCtq3zmtcYb0KEjvpuHCK9pw3dXhuIXIIyYhRMtidnJTXV1N//79rRGLEOIEpRQVeypM48oUbSqiPKn8jHbOAc6U93Hjq44V/NnZQEpnuKl9a16NjCTY1dUGkQshhO2ZndxMnDiRTz/9lJkzZ1ojHiFaJGOVkZLtJbXuzNTk1JzRzj3K3VQrU97Hnamag3xbkA9AJzc3vouO5soGjDAshBCOzOzkprKyknfffZeff/6Zbt264eLiUuv9efPmWSw4IRxVTV4NRZv/q5Up/qcYVaVqtdG4avDu4308mel/vBeTa2tX9EYjr2dl8VTmLsqMRlw0Gh4NC+PxsLA6RxgWQoiWxuzkZufOnfTo0QOA3bt313pPiouFOJNSior0ilqFv+UpZz5icmnlcjyJOTG2jFdvL7RutZOVf4qLuSstjR0nRhi+2NeXd6KjiTtlhGEhhGjpzE5ufvvtN2vEIYTDMFQaKN12Si+mzcXU5J75iMmjs0etZMY92r3ePxCK9XqeyMzkrawsFODv7MzLnTpxR5s2OMkfFUIIUUujx0/PyMhgz549XHLJJbi7u6OUkjs3okWqPlZduxfT1hJU9WmPmHQafC44ZYbsfj64tjp3wa9Siq9yc5mcns7h6uPdvG8LDubVTp1oLQXDQghRJ7OTm7y8PEaNGsVvv/2GRqMhPT2djh07MnHiRPz8/GR+KeHQlFKUp5TXmiG7Ir3ijHYuQS7/jfg7wBfvXt446cybIXt/ZSX3p6fzXV4eAJHu7iyMiuIKKRgWQoizMju5eeihh3BxceHAgQPExsaalo8ePZqHHnpIkhvhUAwVBkq2lvyXzGwuQp+vP6OdR5xHrRmy3TvV/4jpXPRGI69lZfFUZiblJwqGHztRMOwmBcNCCHFOZic369atY+3atbRr167W8qioKPbv32+xwISwheqj1abHS0WbiijdXoqqqf2IycndCe8LvfHt72t6xOQS4FLPFs3zV3Exk1JT+besDICBJwqGY6VgWAghGszs5KasrAwPD48zlufm5qLTyUiowr4Yq40UbS4i/8d88n/Mp+zfsjPauLZxrTVDtlcPL5xczXvEdC5Fej1P7N3L24cPo4CAEwXD46VgWAghzGZ2cnPJJZewdOlSnn32WeB492+j0cjLL7/MZZddZvEAhbC0ir0V5K89nswU/lqIodRQ633PLp61khm3CDerFcsrpVh57BiTMzLIPlEwPDY4mFc6dSJICoaFEKJRzE5uXn75ZS699FK2bt1KdXU106dPJzExkfz8fDZt2mSNGIU4L4ZyA4XrC013Z04vAHZp7ULA0AACrgrA/0p/XIOaJqnYV1HB/enpfJ9/fIThKHd3FkVHc7m/f5PsXwghHJXZyU1cXBw7d+5k4cKFaLVaysrKuP7667nvvvsICQmxRoxCmEUpRXlSuSmZKdxYWGv0X42zBp/+PgRcFUDA0AC8eng16SzZNUYjrx06xKx9+0wFwzPCwpghBcNCCGERGqWUOnczx1FcXIyvry9FRUX4+PjYOhxhITWFNRT8XEDB2gLyf8yn6lBVrfd1YToChh1PZvwv98fZt9FDPJ2X0wuGL/H1ZZEUDAshxDmZ8/1t9if84sWL8fLy4qabbqq1fMWKFZSXlzNu3DhzNymE2ZRRUbKtxFQ7U7ylGE4pnXFyc8J3kO/xuzNXBeAR42HTQSaL9Hoe37uXhacUDL9yomBYBr8UQgjLMju5eeGFF1i0aNEZy1u3bs1dd90lyY2wmuqj1eSvO57MFKwrOGNKA4/OHqZkxvcSX7Tutn/Eo5Tiy2PHeFAKhoUQosmYndzs37+fiIiIM5Z36NCBAwcOWCQoIQCMNUaK/yw21c6UJpTWel/rrcX/Cn9T7YxbBzcbRVq3fRUV3JeezhopGBZCiCZldnLTunVrdu7cSXh4eK3l//77L4GBgZaKS7RQFfsqTHUzBb8UYCip3U3bq5eX6e6MT18fnFwsO96MJdQYjSw4UTBcYTTieqJg+DEpGBZCiCZhdnJz8803M3nyZLy9vbnkkksA2LBhAw8++CA333yzxQMUjs1QYaBwwyndtFNP66bdygX/oSfuzlwZgGtw836Us6WoiElpaew8UTA86ETBcGcpGBZCiCZjdnIzZ84c9u/fz+DBg3F2Pr660Whk7NixPP/88xYPUDiWkxNPnkxmin4vwlhp/K+BFnz7+ZoSGu9e3k3aTbuxCmtqeDwzk0UnCoYDTxQMj5OCYSGEaHJmdQVXSnHgwAGCgoLIyspix44duLu707VrVzp06GDNOC1GuoI3PX2RnoJfCkw9m6oOnNZNu73ONIie32A/XPwsM09TU1BKseJEwfCREwXD49u04eWOHWklBcNCCGExVusKrpQiKiqKxMREoqKiiIqKOq9AAd5++21efvllsrOziY+PZ8GCBQwcOLDe9lVVVTzzzDN88sknHDlyhHbt2vHEE08wYcKE845FWIYyKkoTSk3JTNHmolrdtDU6DX6X+P3XTTvWtt20GyvzRMHwDycKhqNPFAxfJgXDQghhU2YlN05OTkRFRZGXl2eRxObzzz9nypQpvP322wwYMIB33nmHYcOGkZSURFhYWJ3rjBo1iqNHj/LBBx8QGRlJTk4Oer3+vGMR56c6p5qCn44XAuevzafmWO1u2u7R7qZkxm+QH1oP+y2srTEamXfoELNPKRh+vEMHHgsLQ+fU/AqchRCipTF7hOLvv/+eF154gYULF9KlS5fz2vlFF11Er169WLhwoWlZbGwsI0eOZO7cuWe0//HHH7n55pvZu3cvAQEBjdqnPJayDKPeSPGWU7ppbzutm7aXFr/BfqZu2u4R7jaK1LL+PFEwvOtEwfClfn4sio4mxsPDxpEJIYRjs+oIxbfddhvl5eV0794dV1dX3N1rf2nln7hFfy7V1dVs27aNxx57rNbyIUOGsHnz5jrXWb16NX369OGll17i448/xtPTk2uuuYZnn332jDiE5VUeqDQ9air4uQBD8WndtHsc76btP9Qf3/6+OLk6zl2MwpoaZmRm8s4pBcPzIiO5PTjYLh+pCSGEIzM7uVmwYIFFdpybm4vBYCA4OLjW8uDgYI4cOVLnOnv37uWPP/7Azc2NVatWkZuby7333kt+fj4ffvhhnetUVVVRVfVfAWtxcbFF4m8JDJUGin4vMt2dKU8ur/W+c6AzAUNOzNc0xB9diM5GkVqPUoovjh3jwfR0jtYcf9R2R5s2vCQFw0II0WyZndxYenqF0//qVUrV+5ew0WhEo9GwbNkyfH19AZg3bx433ngjb731Vp13b+bOncvs2bMtGrOjUkpRkVbx32zaGwoxVpzSTdsJfPr6mHo2eff2RqN13LsWeysquDctjbUFBQDEnCgYvlQKhoUQollr1NTIe/bsYfHixezZs4fXXnuN1q1b8+OPP9K+fXvi4+MbtI1WrVqh1WrPuEuTk5Nzxt2ck0JCQmjbtq0psYHjNTpKKQ4dOlRnkfOMGTOYOnWq6XVxcTHt27dvUIwtgb5YT8Gv/82mXbmvstb7rm1dTcmM/xX+uPjbTzftxqoxGnn14EFm799P5YmC4Sc6dOBRKRgWQgi7YHZys2HDBoYNG8aAAQP4/fffee6550xTMrz//vt8+eWXDdqOq6srvXv35qeffuK6664zLf/pp5+49tpr61xnwIABrFixgtLSUry8vABIS0vDycmJdu3a1bmOTqdDp3O8xyWNpYyK0n//66ZdvKkYpf+vplzjqsF34H+zaXvGe7aompLNJwqGd58oGL7Mz4+FUjAshBB2xezeUv369eOmm25i6tSpeHt78++//9KxY0f++ecfRo4cSVZWVoO39fnnn3P77bezaNEi+vXrx7vvvst7771HYmIiHTp0YMaMGWRlZbF06VIASktLiY2NpW/fvsyePZvc3FwmTpzIoEGDeO+99xq0z5bYW6o697Ru2kdP66YdeUo37Uv90HrabzftxjIoxdP79vHc/v1SMCyEEM2QVXtL7dq1i08//fSM5UFBQeTl5Zm1rdGjR5OXl8czzzxDdnY2Xbp0Yc2aNabRjrOzs2vNNO7l5cVPP/3EAw88QJ8+fQgMDGTUqFHMmTPH3MNwaEa9kZK/S0y1MyVbS+CUFNbJ0wn/y/+bTdu9U8vuaXa0upoxSUn8WlgIwLjgYF7p1EkKhoUQwk6Zndz4+fmRnZ1NREREreUJCQm0bdvW7ADuvfde7r333jrfW7JkyRnLOnfuzE8//WT2fhxd5aHK/2bT/rkAfWHtgQ09u3makhnfAb446aR2BGBjYSGjk5LIrq7G08mJd2NiGFNPzZcQQgj7YHZyM2bMGB599FFWrFiBRqPBaDSyadMmHn74YcaOHWuNGEUdjFVGCjcen027YG0BZbvLar3v7O+M/xD/48XAQwPQhUrd0amUUrx88CCP792LAYjz8ODL+HhiZfZuIYSwe2bX3NTU1DB+/Hg+++wzlFI4OztjMBgYM2YMS5YsQatt3vUa9lpzo5SiIqPCVDdT+FshxvJTumlrwPtCb1PtjM8FPg7dTft8FNTUMD4lhdUnHqPeFhzMouhoPJv5764QQrRk5nx/m53cnLRnzx4SEhIwGo307NnTInNNNQV7Sm70JXoKfys09Wyq3HtaN+2Q07ppBzp+N+3zta2khBsTE9lXWYmrRsMbUVH8X0iIFA0LIUQzZ9WC4pM6depEx44dgTMH4hONo5SibGfZf7Np/1GEqjmlm7aLBt+LT+mm3bVlddM+H0opFh0+zJSMDKqVIsLNjS/j4+nl7W3r0IQQQlhYo5KbDz74gPnz55Oeng5AVFQUU6ZMYeLEiRYNriWoyauh4Of/umlXZ1fXet8two2AYSe6aV/mh7NXo/PRFqtUr2dSWhqf5uQAcG1gIEs6d8bPRe50CSGEIzL7m3LmzJnMnz+fBx54gH79+gHw559/8tBDD7Fv3z7pln0OyqAo/ue/2bRL/j6tm7aHE36XnTKbdqS73J05D0llZdyYmEhyeTla4IWOHZnWvr2cUyGEcGBm19y0atWKN954g1tuuaXW8uXLl/PAAw+Qm5tr0QAtzRY1N1WHq/6bTfunAvQFp3XT7uL532zaF/uidZPCVktYdvQod6WmUm40EurqyudxcVzs52frsIQQQjSCVWtuDAYDffr0OWN579690ev1dazR8hirjBRt+m827bJdp3XT9nPG/8rj3bT9h/rj1s7NRpE6pkqDgSkZGbyTnQ3AYD8/Po2Lo7UMyieEEC2C2cnNbbfdxsKFC5k3b16t5e+++y633nqrxQKzNzV5NeR8lkP+2nwKfi3AWHZaN+0LvP+bTftCb5ycZRA9a9hbUcFNiYlsLy1FA8zs0IGnwsPRymMoIYRoMRpdULxu3Tr69u0LwJYtWzh48CBjx46tNQP36QmQI6vJqyH9/nTTa5dgl/+6aV/pj2sruWtgbd/k5jIuOZkig4FAZ2eWxcUxNCDA1mEJIYRoYmYnN7t376ZXr17A8bFu4Pi8UkFBQezevdvUrqUVbLpHuRM0OgivHl4EXBWAVzcvNE4t6xzYSo3RyBOZmbx88CAA/Xx8+DwujvZu8rhPCCFaokYP4mev7GkQP3FuWVVV3JyUxB9FRQA81K4dL3TsiKuTPPYTQghH0iSD+Alha78UFHBLUhLHamrw0Wr5sHNnbggKsnVYQgghbEySG2F3jErx3P79zNq3DwV09/Tky/h4Ij08bB2aEEKIZkCSG2FXcquruS05mbUFBQDc2aYNb0RF4S6TXgohhDhBkhthN/4sKmJUUhKHqqpwd3JiYXQ049q0sXVYQgghmhlJbkSzp5Ti9awsHt6zB71SRLu782V8PF29vGwdmhBCiGZIkhvRrBXp9dyZksLKE9N6jAoK4r2YGHyc5VdXCCFE3eQbQjRb/5aWcmNiIhkVFbhoNMzr1In72rZtcWMoCSGEMI8kN6JZ+jA7m/vS06k0GgnT6fgiPp6LZFwiIYQQDSDJjWhWyg0G7ktPZ8mRIwAMDwhgaWwsgS4uNo5MCCGEvZDkRjQbaeXl3JiYyK6yMpyAZyMieCwsDCd5DCWEEMIMktyIZmFFTg53pqZSYjAQ7OLC8rg4LvP3t3VYQggh7JAkN8Kmqo1GHt6zhzeysgC4xNeXz+LiCNHpbByZEEIIeyXJjbCZA5WVjEpM5K+SEgAeCwvj2fBwnGXSSyGEEOdBkhthEz/k5XFbcjL5ej3+zs4s7dyZ/7VqZeuwhBBCOABJbkSTMijFrMxMnjtwAIA+3t6siIsj3N3dxpEJIYRwFJLciCZzpKqKMcnJ/FZYCMC9oaHMi4xEJ4+hhBBCWJAkN6JJ/F5YyOikJI5UV+Pp5MR7MTHcEhxs67CEEEI4IEluhFUZleLlgwd5Yu9eDECchwdfxscT6+lp69CEEEI4KEluhNUU1NQwNiWF7/LyALg9OJiF0dF4arU2jkwIIYQjk+RGWMXW4mJuSkpiX2UlOo2GN6KimBgSIpNeCiGEsDpJboRFKaVYdPgwUzIyqFaKjm5urIiPp5e3t61DE0II0UJIciMsplSv5660NJbn5AAwslUrFsfE4CeTXgohhGhCktwIi0gsK+PGxERSysvRAi926sTUdu3kMZQQQogmJ8mNOG+fHDnCpLQ0yo1GQl1d+Twujov9/GwdlhBCiBZKkhvRaJUGAw9mZPBudjYAV/j7syw2ltaurjaOTAghREsmyY1olD0VFdyUmEhCaSka4KkOHZgZHo5WHkMJIYSwMUluhNm+PnaM8SkpFBkMtHJxYVlsLEMCAmwdlhBCCAFIciPMUGM0MmPvXl49dAiAfj4+fBEXRzs3NxtHJoQQQvxHkhvRIFlVVYxOTGRTcTEAD7Vrx4sdO+Iik14KIYRoZiS5Eef0c34+Y5KTOVZTg49Wy+LOnbk+KMjWYQkhhBB1kuRG1MuoFHP27+fpfftQQHdPT76MjyfSw8PWoQkhhBD1kuRG1Cm3uprbkpNZW1AAwMSQEF6PjMRdJr0UQgjRzElyI87wZ1ERo5KSOFRVhbuTEwujoxnXpo2twxJCCCEaRJIbYaKU4rVDh3hk7170ShHt7s6X8fF09fKydWhCCCFEg0lyIwAo0uu5MyWFlbm5AIwOCuK9mBi8neVXRAghhH2Rby7BjpISbkpKIqOiAheNhvmRkdwbGiqTXgohhLBLkty0YEopPjxyhPvT06k0GgnT6VgRH8+FPj62Dk0IIYRoNEluWqhyg4F709L46OhRAIYHBLA0NpZAFxcbRyaEEEKcH5sPL/v2228TERGBm5sbvXv3ZuPGjQ1ab9OmTTg7O9OjRw/rBuiAUsvL6bt9Ox8dPYoT8HxEBN927SqJjRBCCIdg0+Tm888/Z8qUKTzxxBMkJCQwcOBAhg0bxoEDB866XlFREWPHjmXw4MFNFKnj+CInhz7btrGrrIxgFxd+6d6dGR064CT1NUIIIRyERimlbLXziy66iF69erFw4ULTstjYWEaOHMncuXPrXe/mm28mKioKrVbL119/zY4dOxq8z+LiYnx9fSkqKsKnBdWWVBuNPLxnD29kZQEwyNeX5XFxhOh0No5MCCGEODdzvr9tduemurqabdu2MWTIkFrLhwwZwubNm+tdb/HixezZs4dZs2ZZO0SHsb+ykoEJCabEZkZYGD937y6JjRBCCIdks4Li3NxcDAYDwcHBtZYHBwdz5MiROtdJT0/nscceY+PGjTg3cPyVqqoqqqqqTK+LT8xq3VKsycvj9uRk8vV6/J2dWdq5M/9r1crWYQkhhBBWY/OC4tPHUlFK1Tm+isFgYMyYMcyePZvo6OgGb3/u3Ln4+vqaftq3b3/eMdsDvdHIE3v3cvWuXeTr9Vzg7c323r0lsRFCCOHwbFZzU11djYeHBytWrOC6664zLX/wwQfZsWMHGzZsqNW+sLAQf39/tKdM3Gg0GlFKodVqWbduHZdffvkZ+6nrzk379u0duubmSFUVtyQns76wEID7QkN5NTISnZPNc1khhBCiUcypubHZYylXV1d69+7NTz/9VCu5+emnn7j22mvPaO/j48OuXbtqLXv77bf59ddf+fLLL4mIiKhzPzqdDl0Lqi3ZUFjIzUlJHKmuxkur5f2YGEa3bm3rsIQQQogmY9NB/KZOncrtt99Onz596NevH++++y4HDhzg7rvvBmDGjBlkZWWxdOlSnJyc6NKlS631W7dujZub2xnLWyKjUrx04ABPZGZiBOI9PPgyPp7Onp62Dk0IIYRoUjZNbkaPHk1eXh7PPPMM2dnZdOnShTVr1tChQwcAsrOzzznmjYD8mhrGpaTwXV4eALcHB7MwOhrPUx7hCSGEEC2FTce5sQVHG+fmn+JibkpMZH9VFTqNhjeiopgYEiKTXgohhHAodlFzI86PUoqFhw/zUEYG1UrRyc2NFfHx9PT2tnVoQgghhE1JcmOHSvR67kpL47OcHABGtmrF4pgY/GRuKCGEEEKSG3uTWFbGDbt3k1pRgbNGw4sdO/JQu3byGEoIIYQ4QZIbO/LxkSPcnZZGudFIW1dXPo+PZ4Cvr63DEkIIIZoVSW7sQKXBwOSMDN7LzgbgSn9/lsXGEuTqauPIhBBCiOZHkptmbk9FBTclJpJQWooGmBUezpMdOqCVx1BCCCFEnSS5aca+PnaM8SkpFBkMtHJxYVlsLEMCAmwdlhBCCNGsSXLTDNUYjczYu5dXDx0CoL+PD5/HxdHOzc3GkQkhhBDNnyQ3zUxWVRWjExPZVFwMwNR27XihY0dcZNJLIYQQokEkuWlGfsrPZ0xyMrk1NfhotSzu3Jnrg4JsHZYQQghhVyS5aQYMSjFn/35m79uHAnp4ebEiLo5IDw9bhyaEEELYHUlubOxYdTW3JSezrqAAgP8LCeG1yEjcZdJLIYQQolEkubGhzUVFjEpMJKu6GncnJxZFRzO2TRtbhyWEEELYNUlubEApxfxDh3h07170ShHj7s6X8fF08fKydWhCCCGE3ZPkpokV6fXckZLCqtxcAEYHBfFeTAzeznIphBBCCEuQb9QmtKOkhBsTE9lTWYmLRsP8yEjuDQ2VSS+FEEIIC5LkpgkopfggO5v709OpUooOOh1fxMdzoY+PrUMTQgghHI4kN1ZWZjBwb1oaS48eBeDqgACWxsYS4OJi48iEEEIIxyTJjRWllpdzw+7dJJaX4wQ8FxHB9LAwnOQxlBBCCGE1ktxYyec5OUxMTaXUYCDYxYXP4uK41N/f1mEJIYQQDk+SGwurMhp5eM8e3szKAuBSPz+Wx8bSRqezcWRCCCFEyyDJjQXtr6zkpsRE/ikpAeDxsDBmh4fjLJNeCiGEEE1GkhsL+au4mGE7d1Kg1+Pv7MzHsbFcHRho67CEEEKIFkeSGwuJcXfH39mZSHd3VsTH08HNzdYhCSGEEC2SJDcW4ufiws/duxOq06GTx1BCCCGEzUhyY0ER7u62DkEIIYRo8eQWgxBCCCEciiQ3QgghhHAoktwIIYQQwqFIciOEEEIIhyLJjRBCCCEciiQ3QgghhHAoktwIIYQQwqFIciOEEEIIhyLJjRBCCCEciiQ3QgghhHAoktwIIYQQwqFIciOEEEIIhyLJjRBCCCEcSoubFVwpBUBxcbGNIxFCCCFEQ5383j75PX42LS65KSkpAaB9+/Y2jkQIIYQQ5iopKcHX1/esbTSqISmQAzEajRw+fBhvb280Go1Ft11cXEz79u05ePAgPj4+Ft12c+DoxweOf4xyfPbP0Y9Rjs/+WesYlVKUlJQQGhqKk9PZq2pa3J0bJycn2rVrZ9V9+Pj4OOwvLTj+8YHjH6Mcn/1z9GOU47N/1jjGc92xOUkKioUQQgjhUCS5EUIIIYRDkeTGgnQ6HbNmzUKn09k6FKtw9OMDxz9GOT775+jHKMdn/5rDMba4gmIhhBBCODa5cyOEEEIIhyLJjRBCCCEciiQ3QgghhHAoktwIIYQQwqFIctNAv//+OyNGjCA0NBSNRsPXX399znU2bNhA7969cXNzo2PHjixatMj6gZ4Hc49x/fr1aDSaM35SUlKaJmAzzZ07lwsuuABvb29at27NyJEjSU1NPed69nIdG3N89nQNFy5cSLdu3UwDg/Xr148ffvjhrOvYy7U7ydxjtKfrV5e5c+ei0WiYMmXKWdvZ23U8qSHHZ2/X8Omnnz4j1jZt2px1HVtcP0luGqisrIzu3bvz5ptvNqh9ZmYmw4cPZ+DAgSQkJPD4448zefJkVq5caeVIG8/cYzwpNTWV7Oxs009UVJSVIjw/GzZs4L777mPLli389NNP6PV6hgwZQllZWb3r2NN1bMzxnWQP17Bdu3a88MILbN26la1bt3L55Zdz7bXXkpiYWGd7e7p2J5l7jCfZw/U73T///MO7775Lt27dztrOHq8jNPz4TrKnaxgfH18r1l27dtXb1mbXTwmzAWrVqlVnbTN9+nTVuXPnWssmTZqk+vbta8XILKchx/jbb78pQBUUFDRJTJaWk5OjALVhw4Z629jzdWzI8dn7NfT391fvv/9+ne/Z87U71dmO0V6vX0lJiYqKilI//fSTGjRokHrwwQfrbWuP19Gc47O3azhr1izVvXv3Bre31fWTOzdW8ueffzJkyJBay4YOHcrWrVupqamxUVTW0bNnT0JCQhg8eDC//fabrcNpsKKiIgACAgLqbWPP17Ehx3eSvV1Dg8HAZ599RllZGf369auzjT1fO2jYMZ5kb9fvvvvu4+qrr+aKK644Z1t7vI7mHN9J9nQN09PTCQ0NJSIigptvvpm9e/fW29ZW16/FTZzZVI4cOUJwcHCtZcHBwej1enJzcwkJCbFRZJYTEhLCu+++S+/evamqquLjjz9m8ODBrF+/nksuucTW4Z2VUoqpU6dy8cUX06VLl3rb2et1bOjx2ds13LVrF/369aOyshIvLy9WrVpFXFxcnW3t9dqZc4z2dv0APvvsM7Zv384///zToPb2dh3NPT57u4YXXXQRS5cuJTo6mqNHjzJnzhz69+9PYmIigYGBZ7S31fWT5MaKNBpNrdfqxGDQpy+3VzExMcTExJhe9+vXj4MHD/LKK680y/+Up7r//vvZuXMnf/zxxznb2uN1bOjx2ds1jImJYceOHRQWFrJy5UrGjRvHhg0b6v3yt8drZ84x2tv1O3jwIA8++CDr1q3Dzc2twevZy3VszPHZ2zUcNmyY6d9du3alX79+dOrUiY8++oipU6fWuY4trp88lrKSNm3acOTIkVrLcnJycHZ2rjO7dRR9+/YlPT3d1mGc1QMPPMDq1av57bffaNeu3Vnb2uN1NOf46tKcr6GrqyuRkZH06dOHuXPn0r17d1577bU629rjtQPzjrEuzfn6bdu2jZycHHr37o2zszPOzs5s2LCB119/HWdnZwwGwxnr2NN1bMzx1aU5X8PTeXp60rVr13rjtdX1kzs3VtKvXz++/fbbWsvWrVtHnz59cHFxsVFU1peQkNDsbhOfpJTigQceYNWqVaxfv56IiIhzrmNP17Exx1eX5nwNT6eUoqqqqs737Onanc3ZjrEuzfn6DR48+IyeNXfccQedO3fm0UcfRavVnrGOPV3HxhxfXZrzNTxdVVUVycnJDBw4sM73bXb9rFqu7EBKSkpUQkKCSkhIUICaN2+eSkhIUPv371dKKfXYY4+p22+/3dR+7969ysPDQz300EMqKSlJffDBB8rFxUV9+eWXtjqEczL3GOfPn69WrVql0tLS1O7du9Vjjz2mALVy5UpbHcJZ3XPPPcrX11etX79eZWdnm37Ky8tNbez5Ojbm+OzpGs6YMUP9/vvvKjMzU+3cuVM9/vjjysnJSa1bt04pZd/X7iRzj9Gerl99Tu9N5AjX8VTnOj57u4bTpk1T69evV3v37lVbtmxR//vf/5S3t7fat2+fUqr5XD9JbhroZHe903/GjRunlFJq3LhxatCgQbXWWb9+verZs6dydXVV4eHhauHChU0fuBnMPcYXX3xRderUSbm5uSl/f3918cUXq++//942wTdAXccGqMWLF5va2PN1bMzx2dM1nDBhgurQoYNydXVVQUFBavDgwaYvfaXs+9qdZO4x2tP1q8/pX/6OcB1Pda7js7drOHr0aBUSEqJcXFxUaGiouv7661ViYqLp/eZy/TRKnajsEUIIIYRwAFJQLIQQQgiHIsmNEEIIIRyKJDdCCCGEcCiS3AghhBDCoUhyI4QQQgiHIsmNEEIIIRyKJDdCCCGEcCiS3AghrOLSSy9lypQptg7DRCnFXXfdRUBAABqNhh07dtTZ7uuvvyYyMhKtVtus4hdCNJwkN0KIFuHHH39kyZIlfPfdd2RnZ9OlS5c6202aNIkbb7yRgwcP8uyzz573fvft23fWZEoIYXkycaYQwm4YDAY0Gg1OTub/XbZnzx5CQkLo379/vW1KS0vJyclh6NChhIaGnk+oVlFTU9PsJosUojmSOzdCOLBLL72UyZMnM336dAICAmjTpg1PP/206f267ioUFhai0WhYv349AOvXr0ej0bB27Vp69uyJu7s7l19+OTk5Ofzwww/Exsbi4+PDLbfcQnl5ea396/V67r//fvz8/AgMDOTJJ5/k1BlfqqurmT59Om3btsXT05OLLrrItF+AJUuW4Ofnx3fffUdcXBw6nY79+/fXeawbNmzgwgsvRKfTERISwmOPPYZerwdg/PjxPPDAAxw4cACNRkN4ePgZ669fvx5vb28ALr/8ctM5yMvL45ZbbqFdu3Z4eHjQtWtXli9fXmtdo9HIiy++SGRkJDqdjrCwMJ577jkA0+zsPXv2RKPRcOmll5rWeeaZZ2jXrh06nY4ePXrw448/nnFtvvjiCy699FLc3Nz45JNP2L9/PyNGjMDf3x9PT0/i4+NZs2ZNnedEiBbL6rNXCSFsZtCgQcrHx0c9/fTTKi0tTX300UdKo9GYJmPMzMxUgEpISDCtU1BQoAD122+/KaX+m1C1b9++6o8//lDbt29XkZGRatCgQWrIkCFq+/bt6vfff1eBgYHqhRdeqLVvLy8v9eCDD6qUlBT1ySefKA8PD/Xuu++a2owZM0b1799f/f777yojI0O9/PLLSqfTqbS0NKWUUosXL1YuLi6qf//+atOmTSolJUWVlpaecZyHDh1SHh4e6t5771XJyclq1apVqlWrVmrWrFlKKaUKCwvVM888o9q1a6eys7NVTk7OGduoqqpSqampphmZs7OzVVVVlTp06JB6+eWXVUJCgtqzZ496/fXXlVarVVu2bDGtO336dOXv76+WLFmiMjIy1MaNG9V7772nlFLq77//VoD6+eefVXZ2tsrLy1NKKTVv3jzl4+Ojli9frlJSUtT06dOVi4uL6dhPXpvw8HC1cuVKtXfvXpWVlaWuvvpqdeWVV6qdO3eqPXv2qG+//VZt2LDB3F8NIRyaJDdCOLBBgwapiy++uNayCy64QD366KNKKfOSm59//tnUZu7cuQpQe/bsMS2bNGmSGjp0aK19x8bGKqPRaFr26KOPqtjYWKWUUhkZGUqj0aisrKxa8Q0ePFjNmDFDKXU8uQHUjh07znqcjz/+uIqJiam1r7feekt5eXkpg8GglFJq/vz5qkOHDmfdzunHXp/hw4eradOmKaWUKi4uVjqdzpTMnK6uc6yUUqGhoeq5556rteyCCy5Q9957b631FixYUKtN165d1dNPP33W+IRo6eSxlBAOrlu3brVeh4SEkJOTc17bCQ4OxsPDg44dO9Zadvp2+/bti0ajMb3u168f6enpGAwGtm/fjlKK6OhovLy8TD8bNmxgz549pnVcXV3POIbTJScn069fv1r7GjBgAKWlpRw6dMjsYz2VwWDgueeeo1u3bgQGBuLl5cW6des4cOCAad9VVVUMHjy4wdssLi7m8OHDDBgwoNbyAQMGkJycXGtZnz59ar2ePHkyc+bMYcCAAcyaNYudO3c28siEcFxSUCyEgzu9AFWj0WA0GgFMhbnqlDqYmpqac25Ho9GcdbsNYTQa0Wq1bNu2Da1WW+s9Ly8v07/d3d1rJS11UUqd0ebkMZ1r3XN59dVXmT9/PgsWLKBr1654enoyZcoUqqurTfE1Vl0xn77M09Oz1uuJEycydOhQvv/+e9atW8fcuXN59dVXeeCBBxodhxCORu7cCNGCBQUFAZCdnW1aZskuy1u2bDnjdVRUFFqtlp49e2IwGMjJySEyMrLWT5s2bczaT1xcHJs3b66VpG3evBlvb2/atm17XsewceNGrr32Wm677Ta6d+9Ox44dSU9PN70fFRWFu7s7v/zyS53ru7q6AsfvAJ3k4+NDaGgof/zxR622mzdvJjY29pwxtW/fnrvvvpuvvvqKadOm8d577zXm0IRwWHLnRogWzN3dnb59+/LCCy8QHh5Obm4uTz75pMW2f/DgQaZOncqkSZPYvn07b7zxBq+++ioA0dHR3HrrrYwdO5ZXX32Vnj17kpuby6+//krXrl0ZPnx4g/dz7733smDBAh544AHuv/9+UlNTmTVrFlOnTm1Ut/FTRUZGsnLlSjZv3oy/vz/z5s3jyJEjpiTEzc2NRx99lOnTp+Pq6sqAAQM4duwYiYmJ3HnnnbRu3Rp3d3d+/PFH2rVrh5ubG76+vjzyyCPMmjWLTp060aNHDxYvXsyOHTtYtmzZWeOZMmUKw4YNIzo6moKCAn799dcGJURCtCSS3AjRwn344YdMmDCBPn36EBMTw0svvcSQIUMssu2xY8dSUVHBhRdeiFar5YEHHuCuu+4yvb948WLmzJnDtGnTyMrKIjAwkH79+pmV2AC0bduWNWvW8Mgjj9C9e3cCAgK48847LZKozZw5k8zMTIYOHYqHhwd33XUXI0eOpKioqFYbZ2dnnnrqKQ4fPkxISAh33303AM7Ozrz++us888wzPPXUUwwcOJD169czefJkiouLmTZtGjk5OcTFxbF69WqioqLOGo/BYOC+++7j0KFD+Pj4cNVVVzF//vzzPk4hHIlGnXofVwghhBDCzknNjRBCCCEciiQ3QgghhHAoktwIIYQQwqFIciOEEEIIhyLJjRBCCCEciiQ3QgghhHAoktwIIYQQwqFIciOEEEIIhyLJjRBCCCEciiQ3QgghhHAoktwIIYQQwqFIciOEEEIIh/L/aIQU3UtnwB4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plotVar(varX1Ary, varY1Ary, varX1Total, varY1Total)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 321,
   "id": "510b8282",
   "metadata": {},
   "outputs": [],
   "source": [
    "def SIMPLS(X, Y, numFactor, tolerance = 1e-6):\n",
    "    X = np.matrix(X)\n",
    "    Y = np.matrix(Y)\n",
    "    \n",
    "    varXAry = []\n",
    "    varYAry = []\n",
    "    varSAry = []\n",
    "    \n",
    "    numSample = X.shape[0]\n",
    "    numHuman = X.shape[1]\n",
    "    numPathogen = Y.shape[1]\n",
    "    eig_ary = []\n",
    "    \n",
    "    Y0 = Y - getMean(Y) # Centralize Y\n",
    "    S = (X.T @ Y0) # Calculate S\n",
    "    for a in range(1, numFactor+1):\n",
    "#         eigenvalues, eigenvectors = np.linalg.eig(S.T @ S)\n",
    "#         dominant_index = np.argmax(np.abs(eigenvalues))\n",
    "#         q = np.matrix(eigenvectors[:, dominant_index])\n",
    "\n",
    "        eigenvalues, eigenvectors = sp.eigh(S.T @ S)\n",
    "#         dominant_eigenvalue = eigenvalues[-1]\n",
    "        q = np.matrix(eigenvectors[:, -1]).T\n",
    "        r = S @ q\n",
    "        t = X @ r\n",
    "        t = t - np.mean(t)\n",
    "        normt = np.sqrt(t.T @ t)\n",
    "        t = t/normt\n",
    "        r = r/np.sqrt(r.T @ r)\n",
    "        p = X.T @ t\n",
    "        q = Y0.T @ t\n",
    "        u = Y0 @ q\n",
    "        v = p\n",
    "        if a > 1:\n",
    "            v = v - V @ (V.T @ p)\n",
    "            u = u - T @ (T.T @ u)\n",
    "        v = v/np.sqrt(v.T @ v)\n",
    "        S = S - v @ (v.T @ S)\n",
    "        if a == 1:\n",
    "            R, T, P, Q, U, V = r, t, p, q, u, v\n",
    "        else:\n",
    "            R = np.hstack([R, r])\n",
    "            T = np.hstack([T, t])\n",
    "            P = np.hstack([P, p])\n",
    "            Q = np.hstack([Q, q])\n",
    "            U = np.hstack([U, u])\n",
    "            V = np.hstack([V, v])\n",
    "        varXAry.append(np.diag((P.T @ P)/(numSample - 1))[0])\n",
    "        varYAry.append(np.diag((Q.T @ Q)/(numSample - 1))[0])\n",
    "        # tolerance\n",
    "        if abs(np.sum(varYAry) - getVariance(Y)) < tolerance:\n",
    "            B = R @ Q.T\n",
    "            h = np.diag(T @ T.T) + 1/numSample\n",
    "            return B, T, R\n",
    "    B = R @ Q.T\n",
    "    h = np.diag(T @ T.T) + 1/numSample\n",
    "    return B, T, R"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 322,
   "id": "1b1a70d9",
   "metadata": {},
   "outputs": [],
   "source": [
    "test3X = np.array([[-4,2,1], [-4,-2,-1],[4,2,-1],[4,-2,1]])\n",
    "test3Y = np.array([[430, -94], [-436, 12], [-361, -22], [367, 104]])\n",
    "B3, T3, R3 = SIMPLS(test3X, test3Y, 3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 323,
   "id": "aee3e9f7",
   "metadata": {},
   "outputs": [],
   "source": [
    "def bootstrap_pls(X, Y, n_components=3, n_bootstraps=1000):\n",
    "    B,T,R = SIMPLS(X, Y, numFactor = n_components)\n",
    "    original_coefs = B\n",
    "    \n",
    "    bootstrap_coefs = np.zeros((n_bootstraps, X.shape[1], Y.shape[1]))\n",
    "    \n",
    "    for i in range(n_bootstraps):\n",
    "        indices = np.random.choice(np.arange(X.shape[0]), size=X.shape[0], replace=True)\n",
    "        X_boot = X[indices]\n",
    "        Y_boot = Y[indices]\n",
    "        B,T,R = SIMPLS(X_boot, Y_boot, numFactor = n_components)\n",
    "        bootstrap_coefs[i] = B\n",
    "\n",
    "    # Calculate means and standard errors\n",
    "    coef_means = np.mean(bootstrap_coefs, axis=0)\n",
    "    coef_stds = np.std(bootstrap_coefs, axis=0)\n",
    "    \n",
    "    return original_coefs, coef_means, coef_stds"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 324,
   "id": "7dc41188",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1, 0, 1, ..., 1, 0, 2],\n",
       "       [0, 2, 2, ..., 1, 2, 2],\n",
       "       [1, 2, 0, ..., 0, 0, 1],\n",
       "       ...,\n",
       "       [2, 1, 0, ..., 0, 1, 0],\n",
       "       [2, 2, 1, ..., 2, 2, 1],\n",
       "       [2, 1, 1, ..., 1, 0, 1]])"
      ]
     },
     "execution_count": 324,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "v.X"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 325,
   "id": "14c92131",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0, 0, 1, 1, 0, 0, 0, 0, 0, 1],\n",
       "       [1, 0, 1, 1, 1, 1, 1, 1, 0, 0],\n",
       "       [1, 0, 0, 0, 0, 0, 1, 0, 0, 1],\n",
       "       [1, 0, 0, 0, 0, 1, 1, 0, 1, 0],\n",
       "       [1, 0, 0, 0, 1, 0, 0, 1, 1, 0],\n",
       "       [1, 0, 0, 1, 1, 0, 0, 1, 0, 1],\n",
       "       [1, 0, 0, 0, 1, 1, 1, 1, 0, 1],\n",
       "       [0, 1, 1, 1, 1, 0, 0, 1, 0, 0],\n",
       "       [1, 0, 1, 0, 0, 0, 1, 0, 1, 0],\n",
       "       [1, 0, 0, 1, 0, 1, 1, 0, 0, 0],\n",
       "       [1, 0, 1, 1, 1, 1, 0, 1, 1, 0],\n",
       "       [1, 0, 1, 0, 0, 0, 0, 1, 1, 0],\n",
       "       [0, 0, 0, 1, 0, 1, 1, 1, 1, 0],\n",
       "       [1, 1, 1, 0, 0, 0, 1, 1, 1, 0],\n",
       "       [1, 1, 0, 1, 1, 1, 0, 0, 0, 0],\n",
       "       [0, 0, 0, 1, 1, 1, 1, 1, 0, 1],\n",
       "       [1, 1, 1, 1, 0, 1, 0, 1, 0, 1],\n",
       "       [1, 1, 1, 1, 1, 1, 1, 1, 0, 0],\n",
       "       [0, 1, 0, 1, 1, 0, 0, 0, 0, 1],\n",
       "       [0, 0, 0, 0, 1, 0, 0, 0, 1, 0]])"
      ]
     },
     "execution_count": 325,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "v.Y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 326,
   "id": "7e353880",
   "metadata": {},
   "outputs": [],
   "source": [
    "def jackknife_pls(X, Y, n_components=3):\n",
    "    n_samples = X.shape[0]\n",
    "    B,T,R = SIMPLS(X, Y, numFactor = n_components)\n",
    "    B_original = B\n",
    "    \n",
    "    jackknife_coefs = np.zeros((n_samples, X.shape[1], Y.shape[1]))\n",
    "    \n",
    "    for i in range(n_samples):\n",
    "        X_jack = np.delete(X, i, axis=0)\n",
    "        Y_jack = np.delete(Y, i, axis=0)\n",
    "        B,T,R = SIMPLS(X_jack, Y_jack, numFactor = n_components)\n",
    "        jackknife_coefs[i] = B\n",
    "    \n",
    "    coef_mean = np.mean(jackknife_coefs, axis=0)\n",
    "    coef_std = np.std(jackknife_coefs, axis=0, ddof=1)\n",
    "    \n",
    "    t_stats = B_original / coef_std\n",
    "    \n",
    "    return t_stats"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 327,
   "id": "67c45d84",
   "metadata": {},
   "outputs": [],
   "source": [
    "t_stats = jackknife_pls(v.X, v.Y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 328,
   "id": "0d1e69ba",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "matrix([[-2.30355099,  2.05127306,  0.69339066,  1.17657988,  1.19614326,\n",
       "         -1.77463335, -3.04351424, -0.10921281, -1.98523198,  2.68590327],\n",
       "        [ 0.78373183,  1.36586603, -2.91896636,  5.92020353,  5.65441327,\n",
       "          6.76136268,  1.60977375,  1.65317449, -6.33042535,  1.93720216],\n",
       "        [ 2.62458589, -3.36155281, -1.59883729, -5.15850889, -6.23584339,\n",
       "         -0.38143488,  5.61498215, -2.21145683,  6.57769161, -3.41070092],\n",
       "        [-0.52106451,  1.10123301,  0.73246698,  0.75285639,  0.84832412,\n",
       "          0.25579731, -0.22809666,  0.68134385, -0.19663109, -0.58696256],\n",
       "        [ 2.37166083, -2.31751364, -4.36879665, -2.65551163, -2.9714653 ,\n",
       "          1.17325067,  3.8600053 , -2.52577962,  1.10542453,  1.02105561],\n",
       "        [-1.21818416, -0.21375393,  0.25824819, -2.17343351, -1.67654569,\n",
       "         -3.42801373, -2.16972158, -1.28048766,  0.43481026,  1.6370281 ],\n",
       "        [-1.86584194,  1.9294331 ,  2.80790382,  1.00500425,  1.11432996,\n",
       "         -0.91765453, -1.68840501,  1.48377484,  0.54231242, -1.78809947],\n",
       "        [ 3.98429693, -3.59321654, -3.22867834, -4.70755023, -4.90890708,\n",
       "          3.12333294,  6.56784428, -2.17942318,  5.69382841, -2.94495073],\n",
       "        [ 2.06518207, -3.02654792,  1.19488996, -7.23532314, -8.63483762,\n",
       "         -2.86220348,  3.73304639, -2.42548427, 11.6831704 , -6.2230892 ],\n",
       "        [-0.46768427, -1.35829604, -0.87920397, -3.40018855, -3.40385346,\n",
       "         -3.98216117, -1.99124631, -3.7798472 ,  1.71394077,  2.15239343],\n",
       "        [-1.50385299,  0.36340418,  0.75271586, -0.4692519 , -0.39864619,\n",
       "         -2.14586103, -2.18865341, -0.53483271,  0.15028227,  0.73937417],\n",
       "        [-0.69328625, -1.18545785,  3.33015326, -4.43688723, -6.02174855,\n",
       "         -4.0195824 , -0.37688479, -0.95629451, 10.04209252, -5.26828482],\n",
       "        [ 1.9865474 , -2.82431293, -3.16225052, -3.06229435, -3.41900894,\n",
       "         -0.99973117,  1.40607271, -4.45649197,  0.90880212,  2.5826123 ],\n",
       "        [ 0.71600429, -1.35625385, -1.7415555 , -1.05722076, -0.91057177,\n",
       "         -0.44839234,  0.01226585, -1.52514804, -0.57477664,  1.92583098],\n",
       "        [ 2.445324  , -2.54261744, -0.36198753, -3.81038489, -4.45307594,\n",
       "          0.68015201,  5.67980611, -0.83454523,  5.85499042, -4.56926171],\n",
       "        [-0.90939248, -0.78944979,  0.42598095, -4.03013169, -3.55855772,\n",
       "         -4.80772458, -2.10084495, -1.88835249,  2.07326746,  1.2605138 ],\n",
       "        [ 1.33128561,  0.46154246, -2.90911515,  2.56304819,  2.64491717,\n",
       "          4.41658273,  1.4729321 ,  0.77818174, -4.0488956 ,  3.42622782],\n",
       "        [ 3.22869775, -1.37561168,  0.27752941, -1.00085099, -1.34917232,\n",
       "          4.70881579,  6.34907062,  1.5265128 ,  4.46933533, -7.18222512],\n",
       "        [-1.42135065,  0.47838696, -2.89417248,  0.63901118,  1.10313473,\n",
       "         -1.65789391, -3.52476815, -2.31080005, -3.6646011 ,  6.53581006],\n",
       "        [ 0.80624275, -1.35422979, -6.48829707, -1.83809268, -1.08761029,\n",
       "         -0.27782368, -0.17351776, -3.4935575 , -4.36753883,  7.86113419],\n",
       "        [ 0.43878103,  1.27212661,  0.30266362,  2.6685957 ,  2.29664012,\n",
       "          2.92379506,  1.10494984,  1.94136379, -1.47427582, -0.95007863],\n",
       "        [-2.53963237,  0.60193652, -0.77763582, -0.32689284,  0.05337663,\n",
       "         -2.89580299, -3.89381742, -1.70934922, -2.22972587,  3.79448599],\n",
       "        [ 3.18389865, -0.89655342, -4.0330994 ,  1.65350611,  1.68645852,\n",
       "          6.01889498,  4.99689168,  0.57456075, -2.85695281,  1.01093184],\n",
       "        [ 1.28038755, -0.18420737, -2.39940549,  0.88543892,  1.0639037 ,\n",
       "          2.30462464,  1.46222004,  0.15928955, -1.64693527,  1.22395344],\n",
       "        [ 2.88539206, -0.92795307, -3.29855434,  1.32000198,  1.29134779,\n",
       "          6.88048233,  5.47161886,  0.83733017, -1.46945123, -0.29515945],\n",
       "        [ 0.48303867, -0.94526959, -2.78946881, -1.22112705, -0.93507575,\n",
       "         -1.0483966 , -0.47169221, -2.79778084, -1.18596071,  3.36456893],\n",
       "        [-0.33230078, -1.38668515, -2.03065929, -3.37590353, -2.98465962,\n",
       "         -5.44194134, -3.22755975, -3.91257995, -0.18187226,  4.1065162 ],\n",
       "        [-0.49589271, -0.88700422, -0.99069957, -2.08527196, -1.74042876,\n",
       "         -3.04191204, -1.82072106, -2.40242201,  0.40065529,  2.29505167],\n",
       "        [-2.54737497,  3.52807829,  1.58070604,  6.20340824,  8.21842766,\n",
       "          2.17237841, -3.78725444,  4.68413674, -7.43515741,  2.40541572],\n",
       "        [ 0.63423898, -1.05922879, -0.3974686 , -0.6411455 , -0.62044451,\n",
       "          0.09358693,  0.49176996, -0.2853666 ,  0.41795111, -0.01854954],\n",
       "        [ 0.84153231, -0.206393  , -5.45895857,  1.55902707,  2.51420417,\n",
       "          2.28661821, -0.03252642, -1.64895009, -6.74063582,  6.70475253],\n",
       "        [ 0.96349875, -1.83724332,  0.17021622, -4.91471504, -5.48908252,\n",
       "         -2.54842243,  1.25862361, -1.85157044,  4.61899036, -1.65047264],\n",
       "        [-1.79756538,  2.67585813,  3.94618843,  4.0116757 ,  3.74639703,\n",
       "          1.39548425, -1.75246876,  4.94856671, -1.51282141, -2.3850719 ],\n",
       "        [ 3.40124181, -0.78435616, -2.90546309,  1.52821536,  1.78794914,\n",
       "          8.06704751,  7.74415758,  1.35775446, -1.56790471, -0.91833179],\n",
       "        [ 1.31825392, -1.69759263,  2.08397858, -2.54486789, -2.71132296,\n",
       "          0.23031056,  3.00342732,  0.75539217,  6.03817196, -6.39860619],\n",
       "        [ 2.22477941, -1.77149848, -8.32968599, -1.11684194, -0.63018556,\n",
       "          3.1676252 ,  2.92356058, -2.96561816, -4.29079013,  4.99379605],\n",
       "        [-1.26078912,  1.79351655,  4.08300076,  2.23719969,  2.15274612,\n",
       "          0.07923232, -1.31938875,  2.16872496, -0.04470825, -2.22590203],\n",
       "        [ 2.03217042, -2.21137408, -2.81144272, -1.70596337, -1.98651718,\n",
       "          0.55376407,  2.38204296, -1.75760693,  1.01937344,  0.16553956],\n",
       "        [-1.77137281,  0.61376485,  0.05275456, -0.5602674 , -0.19206371,\n",
       "         -3.31192015, -3.13860291, -1.18821053, -1.68503821,  3.50504316],\n",
       "        [-1.30298003,  1.20648073,  2.63582881,  0.19589616, -0.06086117,\n",
       "         -1.06457173, -1.00438122,  0.90775941,  1.91171297, -2.78452718],\n",
       "        [-1.5152142 ,  0.39655944, -3.12298557,  0.17142624,  0.84388477,\n",
       "         -2.69193618, -5.67248824, -2.9074805 , -7.13020989,  8.91569151],\n",
       "        [-1.75010863, -0.2997236 ,  4.05833262, -2.3562561 , -2.45263659,\n",
       "         -4.43483622, -2.78462107, -1.02187568,  3.5616982 , -1.84647904],\n",
       "        [-0.55525195,  0.35890423, -0.99190726,  0.43175754,  0.77924943,\n",
       "         -0.37858476, -1.30871744, -0.62085841, -1.75363411,  2.56788028],\n",
       "        [ 3.0300408 , -3.046334  , -0.75537263, -4.99949711, -6.28620739,\n",
       "          0.03635636,  5.33174863, -1.78663801,  6.80828085, -4.6875551 ],\n",
       "        [ 3.3177411 , -2.8412542 , -5.71827274, -3.16417981, -3.29460183,\n",
       "          3.83393186,  6.77876329, -2.67284325,  1.60593071,  0.37711238],\n",
       "        [ 3.07628287, -2.66131912, -6.72485324, -2.47166898, -2.73561199,\n",
       "          2.83866868,  4.91749802, -2.77401333,  0.28213097,  1.22557878],\n",
       "        [ 1.91900919, -2.98129148, -1.20477742, -6.50695998, -8.66734296,\n",
       "         -2.72883404,  3.63360788, -3.57713853,  6.61974911, -2.35004871],\n",
       "        [ 0.31038138, -1.55064301,  0.83748383, -2.70673897, -4.43325892,\n",
       "         -2.72354292,  0.19617127, -1.93064618,  4.02303742, -1.65143607],\n",
       "        [ 2.00736618, -0.6440276 , -4.77481505,  1.02777882,  1.68136132,\n",
       "          3.12658591,  2.37981631, -0.85049213, -5.02266651,  4.9650983 ],\n",
       "        [-0.45846915,  1.36838935,  0.11015565,  1.32957388,  1.49198018,\n",
       "          0.59938162, -0.33179318,  0.59042841, -1.35153137,  0.35925966],\n",
       "        [ 1.65141133, -0.65324642,  1.02672369, -0.07554698, -0.45903987,\n",
       "          2.9097094 ,  3.8199797 ,  1.32213629,  3.61601822, -6.28147486],\n",
       "        [-1.04091353, -0.72012884, -1.34263422, -2.60307524, -1.99850336,\n",
       "         -4.3944491 , -3.17541007, -3.17748914, -0.45136566,  3.62628295],\n",
       "        [ 2.2223858 , -1.82583354, -3.39065385, -0.64033779, -0.25734212,\n",
       "          1.39771579,  1.26940574, -1.49439206, -2.7163754 ,  3.6656039 ],\n",
       "        [ 4.00359267, -2.64613438, -1.57925439, -3.73198798, -3.90045309,\n",
       "          4.40175422,  8.85059464, -0.24235606,  7.41702656, -6.017202  ],\n",
       "        [ 0.87259496, -1.91532564,  0.15933083, -2.85341973, -2.81269862,\n",
       "         -1.40883742,  1.01316604, -1.32582042,  2.27262107, -1.2202158 ],\n",
       "        [ 2.89220205, -0.40405023, -1.20822204,  1.38619272,  1.27993964,\n",
       "          4.05416626,  4.88041722,  1.42695346, -0.54175347, -1.35367317],\n",
       "        [-1.4147678 , -0.78575948,  3.09583158, -4.67545323, -4.26986035,\n",
       "         -6.86678449, -2.16920592, -1.08478686,  5.81684671, -2.50104688],\n",
       "        [-1.56988427,  3.03400179,  1.70627969,  7.10327419,  8.74911903,\n",
       "          4.28245857, -1.68127851,  5.18520131, -8.13844719,  0.3639542 ],\n",
       "        [ 2.14102515, -1.35574826,  0.66526133, -1.04747372, -1.22715887,\n",
       "          1.52097162,  3.91980945,  0.85482087,  2.79362988, -3.79488983],\n",
       "        [-3.22756181,  2.4843863 ,  2.22703584,  3.21950638,  3.53227966,\n",
       "         -2.62358375, -5.3329826 ,  0.81864352, -4.20366802,  3.34223018],\n",
       "        [-1.88126153,  1.45548581,  3.48670582,  0.62289306,  0.47344654,\n",
       "         -1.01444061, -1.21626748,  1.57199905,  1.75799518, -2.6061817 ],\n",
       "        [ 1.48555509, -2.71443658, -3.26856117, -7.62809199, -7.31108283,\n",
       "         -2.48178753,  2.24030008, -3.54840356,  5.07521491,  0.69974332],\n",
       "        [-2.39775972,  1.15335833, -0.07887221,  0.63296576,  0.90221648,\n",
       "         -2.1648265 , -4.44635869, -1.14876106, -2.73861521,  5.10010189],\n",
       "        [ 1.7412274 , -1.12016249, -0.69830351, -0.80859348, -0.97580862,\n",
       "          0.44030349,  2.44314887, -0.42583697,  1.0165049 , -0.94572047],\n",
       "        [-1.03529587,  1.12278188,  1.69071519,  0.47460291,  0.30321609,\n",
       "         -0.41366284, -0.47581686,  0.86028798,  1.51818011, -2.29749308],\n",
       "        [-0.33419909, -1.83636893,  1.2798823 , -5.94577253, -6.28427262,\n",
       "         -6.37701879, -1.44931387, -3.02973481,  6.22594133, -1.36615694],\n",
       "        [ 0.04973559, -1.70433698,  2.90599921, -4.14870236, -4.90552599,\n",
       "         -3.00624753,  0.38880968, -1.04719152,  6.69255976, -4.29107989],\n",
       "        [-2.31771177,  2.51547448,  6.39177254,  2.84913819,  2.99807564,\n",
       "         -1.09477605, -3.89832896,  3.2770152 , -0.79224011, -1.56436524],\n",
       "        [ 3.73527238, -2.5664385 , -4.14109513, -1.72812564, -1.57649887,\n",
       "          4.53932253,  5.4575615 , -1.14941977,  0.48177089,  0.31088439],\n",
       "        [ 1.43964327,  0.13401074, -0.80190589,  2.14253585,  1.79470965,\n",
       "          5.5038661 ,  2.81724676,  1.21677221, -0.91233325, -1.33429446],\n",
       "        [ 3.74781714, -2.38424584, -2.35408559, -2.40997417, -2.43249935,\n",
       "          5.19627701,  8.43086022, -0.06588828,  3.67598866, -3.72901373],\n",
       "        [-2.50382801,  0.21087168,  1.14838739, -2.11710835, -1.88257668,\n",
       "         -6.29001802, -5.4081348 , -1.99972891,  0.36816915,  2.28281872],\n",
       "        [ 3.06975726, -1.43501024, -3.26037596, -0.21038862, -0.19119874,\n",
       "          5.06872378,  4.72311483,  0.01712087, -0.15128112, -0.37674284],\n",
       "        [ 1.81138547, -3.10836926, -0.83083524, -9.35893557, -8.74039161,\n",
       "         -2.06801059,  3.49923207, -2.33114915,  9.19509686, -2.57103825],\n",
       "        [-2.1883934 ,  3.46537302,  1.21932302,  5.70280679,  6.07671037,\n",
       "          2.41066391, -3.26405753,  3.43035618, -6.47128467,  2.1879437 ],\n",
       "        [ 2.87424607, -2.8891161 , -0.15761397, -2.80100136, -3.27306979,\n",
       "          1.26828889,  4.09963917, -0.17561849,  5.48713034, -5.76552283],\n",
       "        [-0.24298452,  2.0656475 , -1.9150518 ,  4.26992289,  7.98016596,\n",
       "          3.99376769, -0.24401766,  2.13570627, -7.94993732,  2.95677977],\n",
       "        [-0.991739  ,  0.55953405,  0.1538642 ,  0.20504269,  0.27900018,\n",
       "         -1.18821   , -1.85148475, -0.42677163, -0.6307075 ,  1.1463732 ],\n",
       "        [ 1.18290995, -1.29812384, -5.85749044, -1.07035742, -0.57774835,\n",
       "          0.57384386,  0.49751823, -4.04919191, -3.14154186,  5.45050956],\n",
       "        [ 1.88851417, -0.53411152, -0.26051734,  0.64313238,  0.45008407,\n",
       "          2.47584163,  3.48566596,  1.9217498 ,  0.97854273, -2.87512398],\n",
       "        [-0.46231043, -0.34339164, -3.92423376, -0.85362482, -0.28225532,\n",
       "         -1.67288315, -2.13936921, -2.60440379, -2.53118613,  5.44493924],\n",
       "        [ 2.9381779 , -1.46843142, -1.42038113, -0.43567231, -0.6766937 ,\n",
       "          5.19128859,  7.3989918 ,  1.16267437,  2.25978759, -3.40468798],\n",
       "        [-1.6785033 ,  1.48882821,  4.20468139,  0.94170563,  0.73975855,\n",
       "         -1.35125217, -1.97997941,  1.56254646,  1.08145721, -2.42806541],\n",
       "        [ 2.78625833, -3.07648354, -2.34404921, -3.32005947, -4.40658911,\n",
       "          0.13950638,  3.73178394, -2.51633534,  3.7986155 , -1.19347174],\n",
       "        [ 2.22802819, -0.65450553, -5.01322459,  1.46957641,  1.85389769,\n",
       "          5.02979376,  2.80375849, -0.17158393, -4.03561229,  3.0587612 ],\n",
       "        [ 1.31984744,  0.63384741,  0.39111531,  4.02164185,  3.35068622,\n",
       "          7.51535922,  5.05003795,  2.85007211, -0.71005689, -3.20296574],\n",
       "        [-1.48585272, -0.28124172,  0.32713368, -2.70913768, -2.21445267,\n",
       "         -5.77973172, -3.88189782, -2.02594471,  0.62921326,  2.35152668],\n",
       "        [-0.90508627, -0.58224034, -0.46739368, -3.77980384, -2.93971405,\n",
       "         -3.98729415, -2.45317653, -1.95230959,  0.27359804,  2.68540319],\n",
       "        [ 2.61214478, -1.25351853, -5.36255059,  0.56801563,  1.02376096,\n",
       "          4.56501917,  3.70364062, -0.76822979, -4.49326634,  2.95758081],\n",
       "        [-2.84645176,  3.16822878,  0.79967563,  3.36940062,  4.47562213,\n",
       "         -0.98951832, -3.89512063,  0.85832151, -4.59258404,  4.21997633],\n",
       "        [ 1.95405922, -0.60933627, -5.66145351,  0.95569496,  1.55692994,\n",
       "          3.0665917 ,  2.02128875, -0.5041575 , -4.73078264,  3.19179152],\n",
       "        [ 1.7760421 , -2.30599554, -2.7187896 , -2.33449481, -2.09761939,\n",
       "         -0.31644662,  1.74598112, -2.68919789,  0.746639  ,  1.86339169],\n",
       "        [ 2.40365294, -3.17150613, -4.58436255, -5.77383414, -6.3272696 ,\n",
       "         -0.93875896,  4.43573835, -4.91514697,  5.51938106,  0.18802795],\n",
       "        [-0.30603968,  0.97269057,  4.47165132,  1.26774958,  0.86714627,\n",
       "          0.8589394 ,  0.90955316,  3.48069053,  2.2877043 , -4.55835139],\n",
       "        [-0.8336466 ,  0.1727225 , -4.29541078,  0.47915858,  1.05008481,\n",
       "         -1.45412292, -3.03283559, -2.89905551, -7.36132712,  8.75557268],\n",
       "        [-1.76912054,  0.32605952,  0.41254326, -0.65414805, -0.53821697,\n",
       "         -2.65738561, -2.88647529, -1.04961684, -0.24018752,  1.65712037],\n",
       "        [-2.30540055,  3.22818446,  1.87957088,  4.55036348,  6.67753641,\n",
       "          1.25157312, -2.91772509,  4.03168814, -5.57360117,  1.53510037],\n",
       "        [ 1.86342522, -1.38746647, -4.33084595, -0.69744681, -0.21186719,\n",
       "          2.04522692,  2.15956479, -1.55591331, -2.25606072,  2.95616704],\n",
       "        [-1.08803391,  2.82701984,  1.56521851,  5.28796472,  6.23268597,\n",
       "          4.25152785, -0.59566818,  4.63135685, -5.67461242, -0.58776385],\n",
       "        [-3.17114617,  1.12514544,  0.45054738, -0.22594724,  0.13772768,\n",
       "         -4.37281973, -6.22479089, -1.55742176, -2.69245866,  3.66193017]])"
      ]
     },
     "execution_count": 328,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t_stats"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "feaa2259",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.9.15"
  },
  "latex_envs": {
   "LaTeX_envs_menu_present": true,
   "autoclose": false,
   "autocomplete": true,
   "bibliofile": "biblio.bib",
   "cite_by": "apalike",
   "current_citInitial": 1,
   "eqLabelWithNumbers": true,
   "eqNumInitial": 1,
   "hotkeys": {
    "equation": "Ctrl-E",
    "itemize": "Ctrl-I"
   },
   "labels_anchors": false,
   "latex_user_defs": false,
   "report_style_numbering": false,
   "user_envs_cfg": false
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  },
  "vscode": {
   "interpreter": {
    "hash": "610c699f0cd8c4f129acd9140687fff6866bed0eb8e82f249fc8848b827b628c"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}