PLSDA-cancer-classification-512.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# PLSDA Cancer Classification on GBM Cell Lines"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [],
   "source": [
    "import csv\n",
    "import time\n",
    "\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib\n",
    "\n",
    "import scipy\n",
    "from scipy import stats\n",
    "from sklearn.cross_decomposition import PLSRegression\n",
    "import scipy.io as sio\n",
    "from sklearn.preprocessing import scale, StandardScaler\n",
    "from sklearn.model_selection import LeaveOneGroupOut, LeaveOneOut\n",
    "from sklearn.metrics import r2_score\n",
    "\n",
    "LOG_FLAG = False  # Flag to print debugging/log text\n",
    "\n",
    "# Increase font sizes\n",
    "matplotlib.rcParams.update({'font.size': 16})"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Parsing Data\n",
    "After performing preprocessing on our images, the input data for this program will be a '.csv' file. The following code parses the data in that file and stores it into two arrays. X will be the array that contains the pixel color values for all of the pixels in all 100 of our images. Y will be the array which lists the label (0) mes or (1) pro for each of the cell images. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "csv_filename = '/Users/Ryan/Dropbox/Sharing/project data/data.csv'\n",
    "\n",
    "IMAGE_DIM = 512   # Images cut into sizes of IMAGE_DIM x IMAGE_DIM\n",
    "data = []    # Data matrix to hold csv read data\n",
    "\n",
    "# Open CSV, read row by row, parsing into data matrix\n",
    "with open(csv_filename, 'r') as f:\n",
    "    reader = csv.reader(f)\n",
    "    for row in reader:\n",
    "        data.append(row[1:])  \n",
    "data = np.array(data)\n",
    "\n",
    "# Separate data matrix into X inputs and Y labels. Y labels are indicated \n",
    "# by the last column\n",
    "X = data[:, 0:-1]\n",
    "Y = data[:, -1]\n",
    "X = X.astype(int)\n",
    "Y = Y.astype(int)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Verify and Normalize Parsed Data\n",
    "Read data from data matrix, ensure that the labels are what we expect, and if so, zscore"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def verify_y_vector(y):\n",
    "    '''Ensures that y vector consists of only labels 0 and 1'''\n",
    "    for label in y:\n",
    "        if label != 0 and label != 1:\n",
    "            return False\n",
    "    return True"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Test integrity of parsed Y vector\n",
    "assert(verify_y_vector(Y) == True)\n",
    "\n",
    "# Normalize data\n",
    "x_zscore = scipy.stats.mstats.zscore(X, axis=0)\n",
    "\n",
    "y_zscore = scipy.stats.zscore(Y)\n",
    "Y_mean = np.mean(Y)\n",
    "Y_std = np.std(Y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Inputs matrix size: (100, 262143)\n",
      "Outputs matrix size: (100,)\n"
     ]
    }
   ],
   "source": [
    "# Print data shape\n",
    "print(\"Inputs matrix size: {}\".format(X.shape))\n",
    "print(\"Outputs matrix size: {}\".format(Y.shape))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Preprocessing\n",
    "\n",
    "We characterize some properties of our data. Namely, we first compute R2X and R2Y for increasing number of components. We indicate that a 90% threshold of captured Y-variance is the determining factor in the number of components used to perform our PLS Regression. This process is illustrated below, accompanied by the scores plot to visualize the resulting relationships."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Calculation of R2X and R2Y\n",
    "This will be used to determin how many primary components to use"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def R2 (model, observed):\n",
    "    diff = np.squeeze(np.subtract(model, observed))\n",
    "    sum_sq_err = np.sum(np.square(diff))\n",
    "    sq_obs = np.sum(np.square(observed))\n",
    "    return 1 - (sum_sq_err/sq_obs)\n",
    "\n",
    "def compute_r2x (scores, loadings, observed):\n",
    "    loadings_transpose = np.ndarray.transpose(loadings)\n",
    "    model = np.squeeze(np.dot(scores, loadings_transpose))\n",
    "    return R2(model, observed)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "92.8152576224 % of variation for y accounted for using 2 PCs.\n"
     ]
    }
   ],
   "source": [
    "min_num_pcs = 0     # Stores the min num PCs required \n",
    "                    # to capture 90% of Y-variance\n",
    "max_pcs = 20    # Upper bound on PCs\n",
    "    \n",
    "R2X_range = []\n",
    "R2Y_range = []\n",
    "\n",
    "reach_90_x = False  # Flag set in the for loop below if we reach \n",
    "                    # 90% X variance captured\n",
    "reach_90_y = False  # Flag set in the for loop below if we reach \n",
    "                    # 90% Y variance captured\n",
    "\n",
    "# Compute captured variance for the PC range for both X and Y\n",
    "for i in range(max_pcs):\n",
    "    pls = PLSRegression(n_components=(i+1))\n",
    "    pls.fit(x_zscore,y_zscore)\n",
    "    \n",
    "    # Compute R2X from scores and loadings\n",
    "    R2X = compute_r2x(pls.x_scores_, pls.x_loadings_, x_zscore)\n",
    "    R2X_range.append(R2X)\n",
    "    if (R2X > 0.9) and not reach_90_x:\n",
    "        print((100 * R2X), \"% of variation for x accounted for using\",(i+1), \"PCs.\")\n",
    "        reach_90_x = True\n",
    "    \n",
    "    # Compute R2Y\n",
    "    Y_model = np.squeeze(pls.predict(x_zscore))\n",
    "    R2Y = R2(Y_model, y_zscore)\n",
    "    R2Y_range.append(R2Y)\n",
    "    if (R2Y > 0.9) and not reach_90_y:\n",
    "        print((100 * R2Y), \"% of variation for y accounted for using\", (i+1), \"PCs.\")\n",
    "        min_num_pcs = i+1   # Set Min num PCs here for use later \n",
    "        reach_90_y = True"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "def flatten_logo_output (logo_output):\n",
    "    '''Flatten data from LOGO crossvalidation in order to easier analyze'''\n",
    "    flat = []\n",
    "    for i in range(len(logo_output)):\n",
    "        for j in range(len(logo_output[i])):\n",
    "            flat.append(logo_output[i][j])\n",
    "    return flat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cv_predict_and_report_opt(component_range, X, Y):\n",
    "    # Threshold our results to classify our data into two categories:\n",
    "    # (0) mes and (1) pro. Then evaluate the prediction error for our\n",
    "    # model. We will do this for a range of threshold values from -1 to 1\n",
    "    # and plot the resulting ROC curve\n",
    "    q2y_results = []\n",
    "    opt_num_components = 0  # Corresponds to maximal Q2Y\n",
    "    opt_q2y_Y_test = []     # Holds the optimal Y_test\n",
    "    opt_q2y_Y_pred = []     # Holds the optimal Y_pred\n",
    "    \n",
    "    q2y_did_decrease = False  # Flag to indicate when Q2Y max found\n",
    "    q2y_max = 0\n",
    "    \n",
    "    for num_of_comp in component_range:\n",
    "        num_observations = X.shape[0]\n",
    "\n",
    "        logo = LeaveOneGroupOut()\n",
    "\n",
    "        # Assign a group to each observation\n",
    "        num_groups = 8\n",
    "        groups = []\n",
    "        for i in range(num_observations):\n",
    "            groups.append(i % num_groups)\n",
    "        groups = np.array(groups)\n",
    "\n",
    "        # Convert arrays to numpy arrays to be able to use indexing format\n",
    "        X = np.array(X)\n",
    "        Y = np.array(Y)\n",
    "        X = X.astype(float)  # To silence the warnings of performing int with\n",
    "        # float operations\n",
    "\n",
    "        pls_logo = PLSRegression(n_components=(num_of_comp))\n",
    "\n",
    "        Y_pred = []  # Stores predictions on the test set based on the\n",
    "        # trained model\n",
    "        Y_test = []  # Stores the true labels. This list is in the same\n",
    "        # order as Y_pred so we know the order in which we\n",
    "        # should compare\n",
    "\n",
    "        # Split into training and testing data and predict\n",
    "        for train_index, test_index in logo.split(X, groups=groups):\n",
    "            # Partition\n",
    "            X_train = X[train_index]\n",
    "            Y_train = Y[train_index]\n",
    "\n",
    "            # Normalize X\n",
    "            scaler = StandardScaler()\n",
    "            X_train_zs = scaler.fit_transform(X_train, Y_train)\n",
    "            X_test_zs = scaler.transform(X[test_index])\n",
    "\n",
    "            # Fit model\n",
    "            pls_logo.fit(X_train_zs, Y_train)\n",
    "\n",
    "            # Store data\n",
    "            Y_pred.append(np.squeeze(pls_logo.predict(X_test_zs)))\n",
    "            Y_test.append(Y[test_index])\n",
    "\n",
    "            # Log\n",
    "            if (LOG_FLAG):\n",
    "                print(np.squeeze(np.array(Y_test_pred)))\n",
    "                print(Y[test_index])\n",
    "\n",
    "        # Report the predicted sum of squared error\n",
    "        Y_test_flat = flatten_logo_output(Y_test)\n",
    "        Y_pred_flat = flatten_logo_output(Y_pred)\n",
    "        pred_sum_of_sqr_error = np.sum(np.square(np.subtract(Y_test_flat, \n",
    "                                                             Y_pred_flat)))\n",
    "        sum_of_sqr = np.sum(np.square(Y_test_flat))\n",
    "        q2y = 1 - (pred_sum_of_sqr_error / sum_of_sqr)\n",
    "        \n",
    "        # Check if Q2Y decreases from previous value. If so, mark \n",
    "        # previous component number as maximal q2y, and assign the \n",
    "        # optimal values to their arrays\n",
    "        if (len(q2y_results) > 0 and q2y > q2y_max):\n",
    "            opt_q2y_Y_pred = Y_pred_flat\n",
    "            opt_q2y_Y_test = Y_test_flat\n",
    "            opt_num_components = num_of_comp\n",
    "            print(\"Found new optimal. Number of components: {}\"\n",
    "                   .format(opt_num_components))\n",
    "            q2y_max = q2y\n",
    "            \n",
    "        q2y_results.append(q2y)\n",
    "        \n",
    "        print(\"{} component(s). Q2Y: {}\"\n",
    "              .format(num_of_comp, q2y))\n",
    "        \n",
    "    return q2y_results, opt_q2y_Y_pred, opt_q2y_Y_test, opt_num_components"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1 component(s). Q2Y: 0.7178551058029896\n",
      "Found new optimal. Number of components: 2\n",
      "2 component(s). Q2Y: 0.7142904988557763\n",
      "Found new optimal. Number of components: 3\n",
      "3 component(s). Q2Y: 0.7234547247224408\n",
      "Found new optimal. Number of components: 4\n",
      "4 component(s). Q2Y: 0.730054671840393\n",
      "5 component(s). Q2Y: 0.7269585384758526\n",
      "6 component(s). Q2Y: 0.7252545422994909\n",
      "7 component(s). Q2Y: 0.7232154366726098\n",
      "8 component(s). Q2Y: 0.7248002328938871\n",
      "9 component(s). Q2Y: 0.7253611161132534\n",
      "10 component(s). Q2Y: 0.7248991747728725\n"
     ]
    }
   ],
   "source": [
    "component_range = np.arange(1, 11)\n",
    "q2y_range, opt_Y_pred, opt_Y_test, opt_num_components = cv_predict_and_report_opt(component_range, X, Y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAswAAAJ4CAYAAABxirV5AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3XmYZFV5+PHvyzbADN3AMLKJrALjjuIWSVDBaFCSGPXnGjWLS6LRuEbcYjSuiZgYt7jHNXGLW1CjIOISV9xQFkUGcGRg2LqZGWaGYd7fH+eWfaem6nZ1963q6Z7v53nqqeWee865VbduvXXuuedEZiJJkiSpt13muwKSJEnSjsyAWZIkSWpgwCxJkiQ1MGCWJEmSGhgwS5IkSQ0MmCVJkqQGBsySJElSAwNmSZIkqYEBsyRJktRgxgFzRDw5IrJ2e8o06e9fS3vq7Ks6/2rb8eT5rsuOJCL2iIgXRMT3I+Km2vv0/hnkcW7XftW5bY6IqyLifyPiryJirwHzu19EnBkR50fE1VU+N0bEzyPiAxHxiIhYMuuN3glExBNqn8Ofz3Ddj9c+vxXDquM0dfh2VYd3zEf5i01ELImIG2r7xJPmu04arYjYs/b5P2YW6x/f5zjf67ZxGNvQpuo3b5BteV1DHgdGxB9GxKsi4ksRcX1tvecPWI97RsTLI+LLEbE6IjZFxLqIuDAi3hkRd29vq3uWf6cZfK7H98ljlyqfP4+Id1S/3bfU1jtggHrsHRF/EhFvrz6bG6o8rouIr0fEiyJi/9lu526zXbHmpRHxH5m5uYW8tDB9CHjUkPLeHTiouj0IeFZEnJaZl/VKHBGHA+8Efr/H4vHqthL4U+A3EfHczPyvodR84fsU8DZgH8r79d5BVoqIfYHTq6dfyMy1w6meRuyPgH1rz58I/Mc81WWnUAVafwdcnJk9Aw0teP8D3GO2K0fEe4E/67FoD+D46vaXEfFPmfl3sy1nBO4A/HS2K1eNadcAS3ss3h84qbr9bUQ8JjPPnWkZbXTJuB3wtBby0QIUEccyFSy/HTgaGKMEWbPZL66o1u3cDgB+F/h8tfx44LMRsWuPupwAfIepYPl/KT/qxwPLgSOABwL/DFwHHAL81SzquFPIzA3AJ6qnJ0fE7QZc9dFAp/XegGrx6G5RfsAM9gmp29+z7bG++7Z8/qo2Y++keVtePkAeNwHnAh+cYdlj1f2PgOcD9wRuQ/l9+3/ARUAAL4yIl84w79n4PZrfi4sHyOMKym/Pd2ZQ7q6UYHkz8DHg8cCxlGB5JWV/2wQcCHw+Io6bQd7A3APmS6v7MwY9Va5F5661xy/OzF9l5k2ZuS4zN80iv6zW7dyuy8xvZObplIMJwJ2Ah9dXiojlwGcoX4ZNwP/LzAdn5gcz8+LMvD4zL8/Mr2bmC4CjgDOBnEUddyadgDeAJwy4zhOr++uZ+qMzcpl5n8yMzHz6fNVhsYiIg4AHV08/03mZcuZBmo3NXcf67tv6+a7gDNwyzbY0nYF/KeU3bd/MfADwhhmW/X/A/TPzhMx8Y2Z+PzPXZuZVmflx4D7AJVXal0TEfjPeupm5eZr3ot9v7m8oZyYPzMzDM/NRwDdmUO6twL8Ch2fmozPzI5n5i8y8ITMvysxXMhU3LAX+YaYbNteA+ZXV/cHAM+aYlxamvTsPMvPGIZf1itrjB3Utez1wWPX4KdWBoq/MnMzM5wHPbK96i9J5wKrq8bTBUUQcDfxO9fQ/7aq1aDye0oID8GLge9XjJ/ZOLmkQmfnFzPxZZm6d5fpvzMyvNSyfAF5TPd0TeMBsyhm2qlHr85l5zSzXvzkz/zYz1zSk+QJTDW+nzbSMuQbM3wK+UD3+u4hYNtMMImJV1aH7FQ1pjqh1/L5/j+W/vRgvir+IiG9WHb4nIuJbEfHIrnWOrjqGXxoRG6NcWPaeiDhkwHrvGhF/ExHfjXIx2bqI+F5EPLNXd4Ee6x8dEf8SERdExGRE3BwRv6w66Pc9VRBTF8e9v3r+4Ij4dNXRf0tEnDtI/bvyjIh4bET8T0SsiXKh1tqIODsinhoR2/V1j4j3R0QC76+9tk3n/pnWYwA/qz2+ba3cg5kK5r6RmQOf0srMn3W/FhEHRLkA4wfV/nNLlAsHL4iID0a5GG7g/v8R8ZTa+3LUNGmPrqX9q65lSyPihdX+fH1Vr7VRLmT8eFVOr/5bs1a1Bnygenp8RJw4zSr1oHq77hgRsSIinlTV99Jqv98YEZdHxH9FxAObMq/2z4yIF1XPn1Dtp1dHxNaoXeAX01z0FxF3jYiXRsTXIuKa6v2ciIgfRcQbIuLQhno8pPY5HRQR+0TEP1SfxYbqmHBORPzRNO9XJ78/qrb/8uo9uSEifhoR766+4z2P1RGxX7UN345yYcvm6ljwieneyxnqdMc4PzN/ztQ+cWxE3HeQDKJcMPaMKBfwrolyYdKaKMfQ10bEXRvWPTrKRbw/qt6bm6v95wsR8bTocyFP9bmcERHfqb4zmyLiyuq9PqWhvG0+34Z0r6vSXNRj2dOjdvFalIu7zoxynN8YEddGxOcj4nf7lU/pvwxwXGx/8dR2ZVbrHhsRb672xclqf/xFlN+7Y/ptS239+0XEf1ffic77/Oam78OoVd+X/46I31T7fOeirmdHw8XcMXWR3luq56dXn8Fvovx+ztsZsSG7oPZ4m88xygWDm2OaC/Vj24vAnzesio5A573YJyL2mdGamTmjG/BkymnsBI4B7l57/pIe6e9fW35qj+WrqmWvaCjziFoe9++xvLPsL4HP1p53386o0j8ImOiT5grgkD716KR5CvDlhnLOBZY2bM+zKP1s+q1/C/BnfdY9t0rzfuAfe5U9w89zDDi7oS4JnA8c3LXe+6dZJ2dYj852rWpIs38t/7Nqrz+t9voTZ7pPd5VxB+Dq6bYNOGAGee4LbKzWe/k0af++SrcJ2L/2+sGU02rT1evEuWx/nzodXcv/zdOkvbRKd2Gf5RcOsA1vash/TZXmxZR+at3rvqOW9tvdr9WW3XuAetxAj+NNtf5D6u/5NJ/N8xu25wDgnAHqclCPdU+l9MVvWu/fgJjj539CLb+/rV5bztQx7O0D5rFqmrpe1GfdF1KOiU3rvqjHencFfj3Neu8Cdpvm893uva+le12/ugNPr5ZtrOqypk8dbgUe3VD+wO8X8IJp3qtNwOMbtue5wNY+617Dtt+bx8xiXzq+6TMbYP29KBcjN70vFwNH9Vn/+1WatwD/0mPdz8+wPr/Nr3q+G7DLHL9vd6rVp++xY4Z5/kEtzyf1+Y51lj+2x/KjgMlq+ZfoOqZ01fnE6rU9Wqj3P9fyHfg3d5o8/6vKb8tMP6vZFPbk2gYcU73W2YGvB8a70t+/ln7YAfOl1cHiHymBz/7AfSl9fDpv0KnAjcAPgIdROsffFngeUweaD/WpR6ecX9W+dHem/Hjci21/wPvl8fRams9S+gUeXOVxMnAWUwfRU3qsf261vPND8HlKJ/sDgCOBB8/w8zyrVp8PUH7896+261+YOnj+ANi9tt4SYBnbBqvL6rcZ1qOzXasa0pxUK+tdtdc/VHv9iDl+mb5e5bO+2idWVu/HbarP+LmUg+SMvrzAx6t8L54mXSfw+lTX6x+s7cOvpPwAr6hud6NcvPhV4B5z2f4B3pdr6BFg9Ph8ev4YUq4Ifw3lAH7nqv63oxwn3l/b357UZ/1O0HFlZz+o9tnl1Wf1u7W0TQHzPSkXhf4N5aLSY6rP+XhKK/kPq3XX9vqs2TagubSq19OAw6u6/D5Tfw42A0f3yGNPpn5wE/gw5XTpgZTv872AF1Eu2jmoa917UYKfpPyhfUxV9n7V/vCW2nt5xhw/+05gcQtwm9rrn2bquL+kYf1jKMfczvfq1ZQAen/Kse9UyvUE3+yx7gtq78+FlC4gR1bbeXtKv/rPAi/oWu82tX1lA+UP1jFMHWe/Usv3n6b5fOcaMN8CXA78Engc5WKsFcAjKf02k/LnbKy27q6U4+gbmQoCl3Xd9u4q729rdf4UpWHooGqbH8BUI88Wat+T2vr1oOoXwCOqffEwyr59HVO/fcn8BMwfrK3/35Tf9+VVvq9h6jf8Eno0WjH1fescPz4J3K/K42h6xCjT1Kee32rKd+7W6nP9OPCQWWzjMALm99TyPLbH8qjtHxPAkbVlu1EuvktKY1KvP+/1Ov8EWFc93lg9/2dm8dtMywFz9b3pNJZ+a8brz6LAJ9c2oBMw36naSRJ4VVf6+9fSDztgTuBxPZbvz9S/o1sowd9ePdK9pkqzCdhnmnL+vk9dP1BLc8+uZQcBN1fLXt1n/QD+s7Pj9Vh+bi3//2QOrUeUYaI6eZ3ZJ81za2me2bQ/zHFH7mzXqoY0X6zV5bG1179ZvXbzHOswVsv/WXPJq0fef1jL+9590tynlubhXcs6LYk9P6dh3yhnbzp1e1ifNO9k6s/ebWdZznNo+GPBtq10/zBNXn0D5gHqsQeln27SI+Bk24DqRnoHxEcwFdS+osfyl9fy+JuGugS1lpDq+c+r9b5I7Y9sn/dyA7P8saEM63hNlc//dC17eK3+j2rI45xaPe7VkG63rudHMRUAnUfzWbvudd9erbeVHkELJSA9q5bm2K7lbQbMCVwGLO+R5ndqaZ48k/y70t2OqRb/lzXsR/9dpfluj2UXV8t+Te2PUS3NCbX9OZl7wPxytv8TUL91t2Ler7buB/rk/8Raml7fufof1PfM5jvRkF+/26eYQQMSLQfMlD/QW6r8vtKQ7mBKA0FSjp27Va+/tlaf0waoc7/bzcBfzLDubQfMb6zl94QZrz+LAp9cK/CY2usfqV6brG8Yow2Y+/5joAxR0rceVZo719L0+gfeWbaa/j9SBzB1UHl717JXVK+vok8rXZXu8FpZd+ladm71+hb6dB2ZwWf5uSqvtcCefdIEpe9wAj9u2h/mWJfOdq3qen0Pyhf+o7X35Bds29rdacm7ao512K9WxiPnklePvHcHrq3y/rc+ad5SLb+OrtNZTP0rbqXFYRb1H6MEPAl8rMfyPSmtZAl8eQ7l1Lvd3K7H8k7AfHW/72At7awD5mr9zp/Fc3osqwdUfbvZUE5fJtsHm7sw9eN09gzr1Sn7VsoV4f3S7cLUmainzvI9qP/R6+42sEdtn+55Khu4S239noFcQ9mdUWy2MIPWqapenRauTzSkO5Kphp7XN3y+bQTMfbuKMRWovnUm+fdJdwkNp5mB42p1Orb2ev3sUN99BXhrLd1cA+bpbsd3rfu+6vUN1Lqr9SijczZsdY9lnQB3Y1MeM9iejwIvofzxOZSpeQMeA/y4ti0Dd/WgxYCZMu9A56zlzcCdp0n/sFrZr6UMw9r5jvxLw3orKWcPn0iJo/aj/CbcgfLHqPN93ErDn+se+bYWMFfb1jnrdl7T96Tfrc2psf++emP3ofSHmQ9faFjWGQJvE9DvitJf1h4f3JDX5zLzll4LMvPaWv7361rcmenwHGDPiFjW60YJmK6t0t6zTx1+lJm/aahjo4iIWv0+l5k9Z1XKsqd1Rpy4c5RJKYbp8Nj2wsFNlNPjnVmlrgD+qN/7PxeZeUOVP8Cre12MM4e8b6GcEQB4THRdNBgRu1PGLwb4r9x+dIkfVfcviIiHxgAXlrYpMycpp+ABTo+I8a4kpzM1qUXj2MsRcYcoF7yeH+XiuFtrn/d1taRN42R+ea77QBSPrC6Q+1VErO/a9944QD2g+bjTGXO0+8Kxu1H+XEPtwtkBdY4jPweuaziO7E350Yb+x5HpPKm6n2RqODkAqn20M+nPgyPiwIa6wuy382uZuWoG692dqckLPtYvUZbJjzqjfbT2Xe9jNvvITHTeq68AezfsE6sp4/3CtvvESbXHn2oo55NzqONcdT6jr2Tm9Q3pOvvkIdH/IutvT5PHQDLzsZn56sz8VmauzsxbMnNNZv4npdvU/1ZJHxoRfzzX8mai+k35GKXrEpTrDxonBsnMz1Oue4ASx32M8sf7x0xdgNprvQsz86GZ+YHM/GmWodw2ZubPswzndj/KMSSAN8eIhyGOiDtTum4GpaHi8TmLUUlaC5gz8xdMXTn9jKYri4eoKYC8ubpf2xDs3lx72vSBXjhNPX5e3R/R9XpnpqY/oxy0mm6dH9N+0wr/apo6TGeM8i+wXt9+OiNJBOXU36htoPShehGlxb27vp0DXxvB/POoWl+A86qrpz8aEX8dEbefZt3pdEbvOIDSglX3EKY+816jfLyIcsr1NpR+61dHxKci4rnRMLpAyzqB8J5sP7NjZ3ixm2j4wY2Iv6UcfJ9NOcU7Tv/jUHdQXjen/b+6Ovpsyp/BR1BaG/fuk7ypHtB83NlQ3XfnfXTt8Y+Ymc5x5E5MfxzpDJ004+nJo4w88bDq6cf7/KnuHPN3o/TP7dbZzusz88oZVqGz7kzfn8Nrjwc9th0xwzJmYks2z3bZbx+Zic4+8VdMv090Rgao7xNHVPfXVo0+/Uz3fs7EGVnGSe936x4FpPO5DvqZQv/Pda6/n9PKMg/Bn1POkMDg49jPWZRRdf6DqUm8XpuZ/z7g6i+g9DvehdK3ewOlC+Rs5lUAIDN/zNTQdgcBfUeoaVv1p+lLlOP4TZRuJTM9FgEtBsyVV1J+1PemXGQxare2lAZKcNjPumnW7SzvHmZvuh/eXvbs8/qGPq8Pqj6cyk19U22/fGbDsMxc90x/SzJzaZZJKF6fZUzJbpdV93tGmRp71jLzE5Qv89mU0zcHU1q33wpcEmVIt4GG0eqR93co3Ulg+4Nn5/kvM/P/eqz7f5Q+zp+m9OtcTulD+kbgR1GGvHvobOo1A19hKjj87fi7EbGCqT8An8gyQ+B2IuJk4E2U4OqHlB+Tu1D+BHRmh7xNbZWmofvmuv+/mXIhVFIuiDmNEqAtZ2rfe84A9YDBjindx5Ox2uPpvn/d2jyONHkspXsDwE8i4sTuG2XbO+OePqlHHp3tnNE2VmdQOgHkTN+f2Rzbhnlca+M3p/9K5b2azXCS9X2i81s16G/bSEXEnkx9D9v4vZrr8WMgmbma0g0EypmPUXkb5fsLZQSPgeOxKjCut0R/LzOnayQcxKdrj0fyXlRDIX6F8jt+M+X6m+83r9VfqwFzddrsPdXTp0bEbRuS/3a1AdIMPObtiEw33nS/g0/n+ZnT/LOu317Ras2n1A8qg25P93rDkDn4DEkdX689PrmFCnw1M0+lBE+nUUZd+UG1+L7AuRFxn1lm32k9/sOIGAOo7v+wa3mvev0wMx9OOTNwCmWGqPMo36E7Uqb7fMQs6zWtzLyVcloL4KSIOKJ6/DimvqNN3TE6kxtdAtw3M99Xnb5bm9XskIzgu151J+n8QXllZv5lZn4hyyyV13f2PZrPMs3VZO3xTIO1znHkWzM4jnSf0RhEPQD+V0r3hV63ztnEu/Y429HZzhltY7WvdYKamb4/szm2dR/XBvldgh3gt6l6rzpnR18zg33idbVs+jXydJvxXAttqM5udFpqd6Tfq0F0JuIYdndGACLijZRRTaB0g3rWDNd/DGWioo6To2tOgFmqT0gy9PciIm5DCZaPpDTkPjwzz5tLnm23MEMJLjZShh172QDpO6f5mn6cBppMZIRWTrP8DtX9qq7XO6eBjmb+TVIu0oISbDW5U3WfTPXx3ZF8jtLqCmWM7FZk5o1VIPWyzDyREqRupLS6nTHLbDtD4O1F6QoAZXipPavXP9RnvXq91mfmOVXfuZMp/WE7fX9fPst6DarXVNmd1uZVlAC+n7tV959qOL135znVbjB3ZCrQ+WhDumHWpX69xN36pupt6MeRiFjJ7Po9d7cyd7Zz/4g4rDvxNDrrzvT9WVV7POixbVXX6/XuJwvht6lzlm22+8Sq6v6AiDigId0dGpYN26rqftDPtL7OfOr8obyhMVULIuIfKRcrQ+l//JeZOeifP6pGkM5ET19l6rqFN0bEXD/7ejfdob4XVXeyL1O6Km2hXLD8pbnm23rAnOVCtLdXT/+MMjRQk6uq+6YLa2bTOjJMp1cd6rdTHWw6rZzf7Frc+cBOjYjlw6rcIKovUWee9of1mx2pujiwM0viT3P401/PWLXP1Vs+p53CuSMipjv41ss5h3LBJkz1GZyRLBcadfaLP+26/2ZmzrhvXWb+hKkLXWZVrxmU9XOmTjH+aXUQ7Zxe++A0B+fOPtZ0weLAn90c1Pf1nnWpLm49fYh1+DHl4hOY+fTSnePIgdFj5tOW1APfO07XWklpyQF4XNcFrV/uk+cgOuuePMOuVj+kjPcMU8eu7VTBQedPwde7Fl9Ve9zztyki9qCMIjBMnYaA6S7y7ewTD+lxQe4gvlF7/CcN6YZ2BmsAnc/o1GkuPu9cX/Gb2RxP21Ttt/eonp4/5LJeTBmxA8p1Lk+ozj4Muv6ulLHgxykNMH8K/AWlG95ewH/2ixMGVN+vhvZeVGdsv0Tp7reVMqb/p5vXGswwWpihDEeynjLEynR9Z75T3T8oyhTH24iI4ymTC+xIDqF/C+OZTPX7e2/XsrdSWi6WAu+bbuertn2Y3l3dr6CcGejlWUz9o3/nkOszFy+kXAEO8K6I6L4obRvVVeNvoAzl1nntgKY/MtUB5Yjq6XX90g2g0+3i5Ko/9Mldr3eXu3SA1rlOy9Jc6jWoTivzsUxdUQ1TF4D10/nx+oNefzgj4g/Y9lTgsNR/RLeburr6k/g2hnj6OcsV2m+unp4aEX/dL201mkf9WP0/TI2s8O99Rqeor39wp/vPIKqyOmcPfpzbX2Tby0eq+wMpkzEBkJkXUK4HADgjGqZWj+2nm38rpXVoV+A/mq6sr69bdePq7IuPjIhTe6TflfLd34Xyo/qeriS/ZOpi4n6B/sso2ztMnQvwDow+06NX3kwJrseB9/Rr0Ono8dvyTUpXKYCXV6ezu9c5gTIe+3zp/F7tRfmd3U5EPJ4ykRcM+fdqui6n1f76Pqb+7Ex79nAOdXkWZUIgKH9eH5kzH0Xo7ynD40FpmV6dmdcxNbb1nYF/6lP+dO/FiZQL16H8GT2nIfmsRcTelOPjiZQ6Py0zP9K81gzkzMeyezJT4+Id05DuNbV0nVuvcZg7TeZJaXV5IGUs1iOAv6a0wvyylsf9e+TRWfbkhvq8okqzaprt65tXbVl9pr87VfU9kakpF5P+M/09pZbmAsqFT8dQ+vQcRJl69FmUoelu6rH+udW675/pZ9enPvWZ/t5H+Te8f7Vdb2JqDMZtZvrrtT/MsR6d7Wr8fKbJ40SmJllIyr/MJ1CG1dmfqRnlXsfU9Nfn1ta/P6Xf5EcoF0wcX613aLVffr6Wd9+JJgaoZ32q7MuYGhd03z7pj6D8GH6GctbmzpQRNQ6iHODqk+W8scf6q+b63nblt5xtJzBI4BsDrFcfl/YsSn/wA6r3+ZXVe/CzWprtxnllahzmaWcJo3mmv86yjZQ/9cdWdbk/U2Mnd+qyscf6bYzTu4RtJz74IOXP022q9/gewPMpo/J0z/R3T6YmQbqa8ofxLtX+uqJ6/GTKMGCbgDvN4PP9/VqdXjDgOmO1+nysa1l9pr91wKso3Sz2r/bh+wNvAL7eI9/n1+pyAaXV63DKd+hoygW5n+quJ9vO9Lee8mPduajz95ia1SyBN/TZpvokB2+iHEf2oxxn3le93vltapwae5r3rjNR1Rd7LKtP1vHiah/drbrt2pX2b2ppf0QJ9I+u3quDKd+3v6W01K7tUVZ9pr9LKC2CnZlwn0YJ3nekmf4+QbkQen/KWYBXMTV5y3Qz/b1lpmX3yOvdlLMZL6w+p4Mpf1iOpvyuX1Sr62ca8jmw2o7O7XG19d7ctey4Hus/manxhb9XfWZNE8JsN2U1Zci+zm99r+Pl62t1emiP5RsoXUCeQGlgO4DyXbsn5RjYOTZsBR7R8F7cqWt7P1wr98Fdy/bvWncPtp3c7CXTvA/LmI+psfuk24+pg2Tn1m/CkDO60tVvP2fbGdDu32P9zrInN9TnFbQXMD+VbadW7b6dS/OsVH/O1CQQTbfre6x7Lu0GzGOUFqCmepwPHDzd/jDHenS2q/HzGSCfI6f5bOq3VdRm1GPbSXaabh+k68dqFvX8RFeeH29Ie8SA9fpyr/2OlgPmKs9PdZX9lAHW2Y2pYLTX7TLKwbbzfJgB812Ymmil1+0/KEN0JUMKmKvlB1D+HE/32faajvZ+TE1M0nTbSo8f2oY6f7i23sAzNtb26Y3Afl3L7s7UdMT9bv3eozOYalTpd9tuf6BMH796mvXeTf+p3scow2v1W/e1TZ8v7QTMAfzfoO9XVebGhjp3br/pU5fnMhV8dd+uZdvf4vkImPdiarbCfreLgaP6rN92wDzd+5yUluW9G/J55oD5JD0mQGGw2Qbrt7d0rb8f5dqkpMRbvWZB3p2pmU+vYfs/8NN9P5Ny3VTj7Hoz3JZHdq07yGyD3bcTZ/KZD6tLBlkmgeh52qRH2tdS/s2eS3lTb6a0qryC8g9lTb9158lmyg/mcygtrzdRWjF+QPmXf0pmru+3cma+lxLYvYJyMLyO8u9uPeXL/l+UluhjhrYFU3WZpAx6/3jK4PrXUFozr6N0+n86ZTrbq/pmsgPJzMuyjHDxu5Qr+39EOUuxhbJvXUhpkX04cPvM/O/a6t+ivBevpbTCrKLsi5uqx/9FmWb3T3MGfcP66O5+0Xd0DMrB7H6UU2ZnUybhWU/ZD38NfJYy6cnvN+13LauPhrGRhgkiOjJzC2Vc3844nxsp350LKF2CTmBqgqGhytLv+x6UlsLfUPaPayh/th6bmU+iHFCHXY9rKX/UHkUZduk3lM/1OsrQTu8CHkRpRe5e95uUlvFnUv6IrKF8dzcCl1Na8V9ACR4u7l6/l6rrxsOrp1/PzF/PYHM+XN0vYWoSnk5dz6e0Aj6H8gfhuqquVwHfpZyR7NnfuPp9uCPljN6FlFbqDZR95QuUY+V2Y8xmGfv1eEpL03cpDTibKUH0xykNOH9Z7Ze9yp2kHEdeRxkOclNV7/+ltLLN9sLfgWWJBB5CaYG/gGmGQ8vMd1CuG3oVpbvj9ZTflnWUFs//pDTY9LyAKzPPpLTAf4YSIG/tHzgMAAAgAElEQVSi/JF9G+X7OdMxsVuVmTdnGSno4ZQ6dvb5Gyj9sJ8D3DVH03f53ygjFX2e0qJ9PeU4MkH57v47ZTSgJ2SfoTZ3EO8EDqN81o/NbeejACBL947HUfajFZQuUvVhEB9H+b39FuX3qvP7dA0ljngxpYF1aN1SRiGqyFySJElSD0NrYZYkSZIWAwNmSZIkqYEBsyRJktTAgFmSJElqYMAsSZIkNTBgliRJkhp0T0WqllRjFB5CGWdWkiRpFPahTE7juMEtMmAenkMok0pIkiSN0m0pk/SoJQbMw3MTwJVXXsnY2Nh810WSJC1yk5OTHHbYYeDZ7dYZMA/Z2NiYAbMkSdIC5kV/kiRJUgMDZkmSJKmBAbMkSZLUwIBZkiRJamDALEmSJDUwYJYkSZIaGDBLkiRJDQyYJUmSpAYGzJIkSVIDA2ZJkiSpgQGzJEmS1MCAWZIkSWpgwCxJkiQ1MGCWJEmSGizIgDkijouIv4mI90fETyNiS0RkRLx0jvmeGhFnRcS1EXFzRFwUEa+OiGVt1V2SJEkLy27zXYFZ+ivg2W1mGBHPAc4EEvg6cDXwu8CLgUdExEmZeW2bZUqSJGnHtyBbmIELgH8GHg+sBD44l8wi4gTgjcCtwEMz8+TM/H/A0cDZwHHAO+ZUY0mSJC1IC7KFOTPfXX8eEVvnmOUZQADvy8wv1MrZEBF/AfyK0sp8fGZeNMeyJPWzdStMXAmb18Eey2D8MNhlCP/rF1M5i2lbRlWO27Jzl7OYtkUjsyAD5jZFxB7AQ6unH+lenpmXR8Q3Kd0zHg68doTVk6a3WA7+ay+GCz8H1/4CtmyE3faEA24PK0+HFcdZznyVsdjKcVt27nIW07ZopHb6gBk4Fti7evz9Pmm+TwmYTxhJjaRBLZaD/9qL4dvvgA3XwfihsPtSuGU9XPUTmFgN93m65SzmbRlVOW7Lzl3OYtoWjZznBuDI6v7GzLypT5oru9JK09u6FW64HK7+WbnfOteeQ106B+WrfgJ77w/Lb1/ur/pJeX3txQujnK1bSzC+4TpYcTwsGYNddi33K44vr1/0+bm/f4upnMW0LaMqx23ZuctZTNuieWHADPtU9+sb0qyr7sf6JYiIJREx1rnV8tXOaO3F8I0z4auvga+9odx/48z2gtjFdPCfuLK0XI8fChHbLouAsUNh7SUl3VwspnIW07aMqhy3ZecuZzFti+aFAXN7zgAmardfz291NG9G0fK7mA7+m9eVbh67L+29fI+9y/LN63ov3xnLWUzbMqpy3Jadu5zFtC2aFwbM0OmG0WfvBqAzcclkQ5rXAuO1223nXjUNzbC6S4yq5XcxHfz3WFb6RN/S5yTP5g1l+R5znD9oMZWzmLZlVOW4LTt3OYtpWzQvDJhhVXW/b0T060ZxWFfa7WTmpsyc7NyYCsS1oxlmd4lRtfwupoP/+GHlAsKJ1ZC57bJMmFwNK44t6eZiMZWzmLZlVOW4LTt3OYtpWzQvDJjhYmBD9fjEPmk6r58//OpoqIbdXWJULb+L6eC/yy5ltI29l8Pai2DjJGzdUu7XXgRLl8PxD5v7EHaLqZzFtC2jKsdt2bnLWUzbonkR2f0juABFxPuBJwEvy8x/nMX6HwMeBbw7M5/Stexw4FJgV2DloBOXVBf+TUxMTDA21vdaQY3S1q2lJfmqn5TuEfUW4MxyMDvkrnC/58z+YHbD5aXFeu/9SzeMbhsn4ebr4QEvhv0On10ZHfWhi8YOLcH45g0liF26HO49hCGShl1O99B1K44tPy7DHiJvoZazmLZlVOW4LTt3OYtpW3qYnJxkfHwcYLw6262W7DQBc0Q8E3gm8N3MfGLXsrtTxlreCjwsM79Yvb438FngFOCTmfnIGdTJgHlHM4pgdhRBed1iO/gvlklYRlnOYtqWUZXjtuzc5SymbeliwDw8C3LikirAfVvtpaOr+6dFxMNqrz88M6+qHh8AHAes6c4vM8+PiOcBZwJnRcTXgGsok5UcTOm28fR2t0IjN0h3iZt+M7fuEp3TcROrS3Dcq0W2zdNxK44r3UqGfVAeVTm77DL3lvedrZzFtC2jKsdt2bnLWUzbopFZkAEzZTzke/d4/bZsOzrFkkEzzMw3RcRPgecB96KMmnEFZfSL1zZMaqKFon4BW68W5rYulFtxXJnJqdMie9NvSr6H3HU4p+M8+EuSNFSLokvGjsguGTugUXeXmIfTcZKknZddMoZnobYwa7EaZpA56u4StshKkrQoGDBrx9HrwrIDbl+C3La6MYy6u4QkSVrwDJi1Y6gPXTZ+aLkw75b1pfvExOoS5LYZNI/iAjZJkrQoGDBr/nVPJ93pW7xkDFbsU7pPXPT5EuTaXUKSJI2YTWqaf6OaTlqSJGkWDJg1/0Y1nbQkSdIsGDBr/tXHR+6lrfGRJUmSZsGAWfNv/LAyGsbE6jIecl1mGfJtxbElnSRJ0ogZMGv+dcZH3nt5ucBv4yRs3VLu117U/vjIkiRJM2AEoh1DZ3zkg+8CN18P1/2y3B9yV7h3i0PKSZIkzZDDymnH4fjIkiRpB2TArB2L4yNLkqQdjE13kiRJUgMDZkmSJKmBAbMkSZLUwIBZkiRJamDALEmSJDUwYJYkSZIaGDBLkiRJDQyYJUmSpAYGzJIkSVIDZ/rTYLZudcpqSZK0UzJg1vTWXgwXfg6u/QVs2Qi77QkH3B5Wng4rjpvv2kmSJA2VAbOarb0Yvv0O2HAdjB8Kuy+FW9bDVT+BidVwn6cbNEuSpEXNc+rqb+vW0rK84TpYcTwsGYNddi33K44vr1/0+ZJOkiRpkTJgVn8TV5ZuGOOHQsS2yyJg7FBYe0lJJ0mStEgZMKu/zetKn+Xdl/ZevsfeZfnmdaOtlyRJ0ggZMKu/PZaVC/xuWd97+eYNZfkey0ZbL0mSpBEyYFZ/44eV0TAmVkPmtssyYXI1rDi2pJMkSVqkDJjV3y67lKHj9l4Oay+CjZOwdUu5X3sRLF0Oxz/M8ZglSdKiZqSjZiuOK0PHHXwXuPl6uO6X5f6Qu8K9HVJOkiQtfo7DrOmtOA6W396Z/iRJ0k7JgFmD2WUX2O/w+a6FJEnSyNlEKEmSJDUwYJYkSZIaGDBLkiRJDQyYJUmSpAYGzJIkSVIDA2ZJkiSpgQGzJEmS1MCAWZIkSWpgwCxJkiQ1MGCWJEmSGhgwS5IkSQ0MmCVJkqQGBsySJElSAwNmSZIkqYEBsyRJktTAgFmSJElqYMAsSZIkNTBgliRJkhoYMEuSJEkNDJglSZKkBgbMkiRJUgMDZkmSJKmBAbMkSZLUwIBZkiRJamDALEmSJDUwYJYkSZIaGDBLkiRJDQyYJUmSpAYGzJIkSVIDA2ZJkiSpgQGzJEmS1MCAWZIkSWpgwCxJkiQ1MGCWJEmSGhgwS5IkSQ0MmCVJkqQGBsySJElSAwNmSZIkqYEBsyRJktTAgFmSJElqYMAsSZIkNTBgliRJkhoYMEuSJEkNDJglSZKkBgbMkiRJUgMDZkmSJKmBAbMkSZLUwIBZkiRJarCgA+aIeFREnBsRN0TE+oj4cUS8MCJ2n0VeSyPijIj4fkRMRsQtEbEmIj4fEX84jPpLkiRpx7fbfFdgtiLiX4BnA1uAc4B1wAOB1wOnR8TvZ+bNA+a1HDgPuEOVz7eAG4FjgIcCD42IN2fms1vfEEmSJO3QFmQLc0T8MSVYXgfcOzMfnJmPAG4P/BQ4CXjVDLJ8OSVY/gFweJXfozPzHpSAeQvwrIi4T5vbIUmSpB3fggyYgRdX96/LzPM7L2bmtcBfV0+fGRHjA+b3wOr+9Zl5fX1BZp4FfLV6et9Z1leSJEkL1IILmCPiUOCe1dOPdC/PzG8AVwJLgNMGzHbjgOmuHTCdJEmSFokFFzADJ1T312fmZX3SfL8r7XS+UN3/XUTsX18QEacBDwDWAJ+dSUVHZutWuOFyuPpn5X7r1vmukSRJ0qKxEC/6O7K6v6IhzZVdaafzeuBewIOByyPim0xd9HcP4JvAX2TmxMyrO2RrL4YLPwfX/gK2bITd9oQDbg8rT4cVx8137SRJkha8hRgw71Pdr29Is666Hxskw8xcHxGnA68BnkcJnDuuA74CrG7KIyKWULqBdNdzeNZeDN9+B2y4DsYPhd2Xwi3r4aqfwMRquM/TDZolSZLmaCF2yWhdRBxMaUX+G+ClwFHAMkqr8w+Avwe+ERFNQfAZwETt9uth1pmtW0vL8obrYMXxsGQMdtm13K84vrx+0eftniFJkjRHCzFgvqm6X9qQZll1Pzlgnv9BuZDwZZn5msy8LDPXZ+b3gIdRhqq7K/D8hjxeC4zXbrcdsOzZmbiydMMYPxQitl0WAWOHwtpLSjpJkiTN2kIMmFdV94c1pOksW9WQBvjtqBsPqp5+tHt5Zt4CfKJ6emq/fDJzU2ZOdm5MBfbDsXld6bO8e5//DXvsXZZvXtd7uSRJkgayEAPmH1b3yyOi30V9J1b35/dZXne72uN+LdKdi/3277N89PZYVi7wu6VPV+7NG8ryPZb1Xi5JkqSBLLiAOTN/DXyvevq47uURcRKlhXkTcNYAWdYv5rt3nzSdGf76DWM3euOHldEwJlZD5rbLMmFyNaw4tqSTJEnSrC24gLnymur+RRFx986LEbEceFv19C31YeAi4uERcVFEnF3PKDOvYCoA/9eIOKK+PCKeADy6errdRCnzZpddytBxey+HtRfBxknYuqXcr70Ili6H4x9W0kmSJGnWIrtbJxeIiPhX4FnALcDZlGHmTgH2pYx48aDMvLmW/snA+4DLM/OIrrzuRJn++gDKrH/fpszqtxK4Y5XsQ8ATc8A3LCLGgImJiQnGxgYa3W52eo3DvOLYEiw7pJwkSTuNyclJxsfHAcar66nUkoU4DjMAmfnsaoKRZwC/A+wOXAq8DnhTZm6eQV4XVEHzc4A/oIyYsQS4AfgS8N7M/FjLm9COFcfB8tuX0TA2ryt9lscPs2VZkiSpJQu2hXlHN7IWZkmSJGxhHiabISVJkqQGrQbMEXFJRPxPm3lKkiRJ86ntFubbMjVmsSRJkrTgtR0wXw7YYVeSJEmLRtsB8yeB34uIA1rOV5IkSZoXbQfMrwEuBr4UEfdqOW9JkiRp5Noeh/nTlCmpfwf4v4j4DbAKuLlH2szMB7dcviRJktSqtgPmU2uPAzi0uvXiANCSJEna4bUdMD+o5fwkSZKkedVqwJyZZ7eZnyRJkjTfnOlPkiRJatB2l4zfiogDgd9jqg/zauC8zLx6WGVKkiRJbWs9YI6IfYA3A48Hdu1afGtEfAh4dmbe1HbZkiRJUttaDZgjYgnwFeBEyigZPwMurRYfBdwJeBJwh4g4OTM3tVm+JEmS1La2+zA/E7gn8GPgPpl558z84+p2F+DewI+qNM9ouWxJkiSpdW0HzI8BJoEHZ+Z3uxdm5veAhwDrgMe2XLYkSZLUurYD5uOAr2bm2n4JqmXnAMe3XLYkSZLUurYD5l2BzQOk2zyEsiVJkqTWtR20rgJ+LyL27JegWvZ7VVpJkiRph9Z2wPxZ4EDggxGxf/fCiFgOfBC4DfCZlsuWJEmSWheZ2V5mJUj+EWWykvXA/wCXVYuPAk4DlgJXAnfPzOtbK3wHExFjwMTExARjY2PzXR112ZpbuWr9Vay/ZT1Ld1/KwUsPZpewl5AkaeGanJxkfHwcYDwzJ+e7PotJq+MwZ+b1EfFA4KPAPYBHA52IPKr7HwCPWczBsnZsv7rxV5x9xdlcNnEZm27dxJJdl3Dk+JGccrtTOGrfo+a7epIkaQfT+kx/mflL4J4RcTJwf7adGvvczPxa22VKg/rVjb/iwxd+mBs23cBBex/EXrvtxc1bbubC6y9kzfo1PH7l41sLmm3FliRpcWg9YO6oAmODY83IMIPMrbmVs684mxs23cDR40cTUU56LNtjGUfvfjSXTlzK2VeczRHjR8y5zFG2Yi+mwHxU27LYypEkDdfQAmZppoYdZF61/ioum7iMg/Y+6LfBckdEcNDeB3HZxGVctf4qDl12aJ9cpjfKVuzFFJiPalsWWzmSpOGbU8AcEd8FPpeZr+qx7BBgfWZOzKUM7RxGEWSuv2U9m27dxF677dVz+V677cXVG65m/S3rZ13GqFuxF0tgPqptWWzlSJJGY67NQycCR/ZZdiVw5hzz1wC25lZWr1vNJTdcwup1q9maW+e7SjPSHWQu22MZu+6yawkyx4/mhk03cPYVZ895u5buvpQluy7h5i0391x+85abWbLrEpbuvnTWZcykFXsuRvWewVTwd+H1F7Lvkn05YuwI9l2yLxdefyEfvvDD/OrGXy2IbVls5UiSRmeYXTKCqZExNCSL4bTvqLpKHLz0YI4cP5ILr7+Qo3c/epuyMpM1G9awcv+VHLz04FmXMYpWbBjdezaKFvNRbctiK0eSNDpefbKADbvlr26YrdiDBJmbbt005yBzl9iFU253Cvst2Y9LJy5l3eZ13Lr1VtZtXselE5ey35L9OOV2p8ypq8QoWrFhdO/ZKFrMR7Uti60cSdLoeNHfArWYRnyoB5nL9li23fK2gkyAo/Y9isevfPxvt+fqDVezZNclrNx/ZSvbM4pWbBjdezaKFvNRbctiK0eSNDoGzAvUYhrxYVRBZsdR+x7FEeNHDGXEh04r9pr1a7h04tJt3rM1G9a00ooNiyswH9W2LLZyJEmjY5eMBWoUp31HdfHSKLpK9Crz0GWHcux+x3LoskNbzbvTir1y/5XcuOlGVk2u4sZNN7Jy/5WtjY4wqvesE/yt2bCGzNxmWSf4O3L8yDkFf6PalsVWjiRpdKL7R3BGK0dsBX4D/LzH4lMblgFkZj541oXv4CJiDJiYmJhgbGys9fxXr1vNW3/4VvZdsm/Plr91m9dx46YbecYJz5h1C/MoyqhbDBcw1o1i0opRvGf9zjJ0WszbHIptMY2PvNj2Z0k7vsnJScbHxwHGM3NyvuuzmLTRJeOQ6jbTZbOP1LWoRnzoGGZXifnQacUeplG8Z8Pu910vZxSf/2IrR5I0fHMNmF/dSi00Y6PoKzsfFy+NIshcbBZLYA6j+/wXWzmSpOGaU8CcmS9rqyKaucUy4oMWBoM/SdLOylEyFrjFMOKDJEnSjmxOF/2pv2Ff9DdKXrwkSdKOz4v+hscWZk3Li5ckSdLOzIBZA7H/qiRJ2lnZRChJkiQ1MGCWJEmSGhgwS5IkSQ0MmCVJkqQGBsySJElSgzmNkhERvzOX9TPzW3NZX5IkSRq2uQ4r9w1gtjOfZAvlS5IkSUM114D1W2wfMO8O3Kt6vA64vHp8OLCsSv894JY5li1JkiQN3ZwC5sw8qf48IpYAXwYuA16YmZ/sWv4nwOuBzcCD5lK2JEmSNAptX/T3EuDuwAO6g2WAzPwUcApwD+ClLZctSZIkta7tgPkxwFcz84p+Capl5wCPbrlsSZIkqXVtB8yHUfotT2c9cNuWy5YkSZJa13bAfCNwUkTs3i9Btex+wETLZUuSJEmtaztg/iJwCPDeiBjvXhgRY8C7qzRfbLlsSZIkqXWROdthlHtkFnEYcD6wP3ATcBZlxAyAI4DTgDHgOuDEpr7OC13152BiYmKCsbGx+a6OJEla5CYnJxkfHwcYz8zJ+a7PYtLqxCGZeWVEnAx8CLgb5SLATkQe1f2PgScs5mBZkiRJi0frM+1l5s+Bu0fE7wH3Z+rivtXA1zLz3LbLlCRJkoZlaFNTZ+Z5wHnDyl+SJEkahbYv+pMkSZIWlaEEzBFx34j4SESsioj1EfGu2rJTIuKVEXHgMMqWJEmS2tR6l4yIOAN4FdsG47vWHq+jTKF9FfD2tsuXJEmS2tRqC3NEPBh4NbAGeBxwaHeazPwOcC3wsDbLliRJkoah7RbmvwU2A3+QmT8FiIhe6X4M3L7lsiVJkqTWtd2H+Z7AdzvBcoNrgINaLluSJElqXdsB8zJK3+TpjLNtv2ZJkiRph9R2wHwNcMwA6Y4Dft1y2ZIkSVLr2g6YvwHcLSLu2y9BRJxGCaq/2nLZkiRJUuvaDpjPrO4/FREPi4ht8o+IU4F3A1uAt7RctiRJktS6VgPmzPw+8ELgQOAzwA1AAn8SEdcBX6Jc7PeCzLygzbIlSZKkYWh9pr/MfCNljOUfAvsAAYwB+wEXAg/PzH9tu1xJkiRpGFqf6Q8gM88CzoqI2wBHUQLzKzPzymGUJ0mSJA3LUALmjsy8hjJyhiRJkrQgtT019uaIeOcA6d4REZvaLFuSJEkahrb7MO/GYK3Wuw6YTpIkSZpXrV/0N6C9gVvmqWxJkiRpYCMPmCNiDLgfg02hLUmSJM2rOXeLiIhLul76k4g4qaG8Q4DdgXfNtWxJkiRp2NroR3xM7XFSxlwea0i/hTKpyRktlC1JkiQNVRsB8+2r+wAuAT4JvKhP2s3A1Zm5uYVyJUmSpKGbc8CcmZd2HkfEh4Gv1V+TJEmSFrJWh3bLzD9tMz9JkiRpvs3XsHKtiIhHRcS5EXFDRKyPiB9HxAsjYvc55PlHEfHZiFhTTcRyTUR8KyJe3mbdJUmStDBEZrabYcQy4GnAqcChwJ59kmZmHjeHcv4FeDblIsJzgHXAA4F9gW8Av5+ZN88gvz2ADwGPAm4G/g+4GjgIuCOwa2YeMIP8xoCJiYkJxsaaroGUJEmau8nJScbHxwHGM3NyvuuzmLTaJSMiDgbOA46iXATYZNaRekT8MSVYXgecnJnnV68fQAmeTwJeBTx/Btm+ixIsfxp4SmZeWytvF+Bes62vJEmSFq5WW5gj4r3Ak4GfAq8HLgT6/sOZ7cWBEfFd4J7ASzPz1V3LTgK+DmwCDszMiQHyOwX4CnABcPfMnPMshLYwS5KkUbKFeXhabWEG/gC4htLqe2PLeQMQEYdSgmWAj3Qvz8xvRMSVwGHAacBHB8j2b6r7f2kjWJYkSdLi0XbAvB/whWEFy5UTqvvrM/OyPmm+TwmYT2CagDkidgVOqZ6eFxEHAY8BjqO0Uv8Q+GRmrptrxSVJkrTwtB0w/7rl/Ho5srq/oiHNlV1pmxwFLKse3wd4W+15xz9FxGMy85yBaylJkqRFoe1h5T4J/G5E7N1yvnX7VPfrG9J0WoMH6Ty8vPb4PcAPKF0+9gHuBpwFrAA+ExG33371IiKWRMRY51arpyRJkhawtgPmVwJXAf8VEStazntY6qN5rAYenJnfz8x1mflj4A8pFwMuo/+U3wBnABO12yha2yVJkjRkbXfJeBNlZIxHAr+oRrO4AtjaI21m5tNmUcZN1f3ShjSdLhWDXCF6U+3x+zNzU31hZt4aEf8O/BtlbOl+XgucWXu+DwbNkiRJC17bAfNfMjW+8hjNAWZSJjiZqVXV/WENaTrLVjWkqeeXlJbmX/VJ03n94H6ZVIH2b4PtiOmGoZYkSdJC0HbA/JSW8+vlh9X98og4ss9IGSdW9+dPl1lmrouIi4HjgX4z+XVed6QMSZKknUyrAXNmvqfN/PqU8euI+B7lwrzHAb0mLjmM0tp71oDZfhx4GaVF/E09lj+ouv/ubOosSZKkhavti/5G5TXV/Ysi4u6dFyNiOWVYOIC31Gf5i4iHR8RFEXF2j/zeDNwAnBYR23QTiYjHAI+vpZMkSdJOZEEGzJn5aUrwugz4dkR8ISI+AfwSuDPwTUqLcd04ZTKSo3vkdy3waGAj8I6IuCAiPh4R51MmPgngVZk5aIu1JEmSFolWu2RExDtnkHy2o2R0Vn52RHwTeAbwO8DuwKXA64A3ZebmGeb35Yi4K/BiSteMP6KMsnEW8K+Z+b+zraskSZIWrsjM6VMNmllEr+Hj6jqFBSVg3rW1wncw1eQlExMTE4yNDTJ/iiRJ0uxNTk4yPj4OMJ6ZgwytqwGNapSMXYDDgdMos+f9K2UyEEmSJGmHNtJRMiLiZcA/A38GnNBm2ZIkSdIwjPSivyz9P/4OWA+8apRlS5IkSbMx8lEyMnMLZUKRB02XVpIkSZpv8zWs3BJgv3kqW5IkSRrYyAPmiDgWOAm4ctRlS5IkSTPV9jjMj2tYvA+wEngSsBfwn22WLUmSJA1D28PKfYipsZZ7ier+f/CiP0mSJC0AbQfMH6F/wLwZWA18JTPPa7lcSZIkaSjaHof5CW3mJ0mSJM23+RolQ5IkSVoQ2u6SsY2IWAEcWj1dnZlrh1meJEmS1LahtDBHxJMi4gJgDfCD6rYmIn4aEU8cRpmSJEnSMLQeMEfEO4D3AnegjIpxEzBZPb4j8L6IeHvb5UqSJEnD0GrAHBGPAJ5KCZBfBKzIzH0zcz/gAODvgBuBp0bEn7RZtiRJkjQMbbcwPw3YAjwoM9+Qmdd1FmTm9Zn5T8CDga1VWkmSJGmH1nbAfA/gvMz8fr8E1bKvASe2XLYkSZLUurYD5qXA1QOku7pKK0mSJO3Q2g6YrwbuOkC6uzJYYC1JkiTNq7YD5q8CKyPi+f0SRMQLKCNonNNy2ZIkSVLrIjPbyyxiJWXM5SXAt4EPAJdVi48CngTcC9gE3CMzL2yt8B1MRIwBExMTE4yNjc13dSRJ0iI3OTnJ+Pg4wHhmTs53fRaTVmf6y8wLI+IxwAeB+wL36UrSGZf5TxdzsCxJkqTFo/WpsTPzsxFxHPB04GRqU2MD5wLvzMyr2i5XkiRJGobWA2aAzFwDvGIYeUuSJEmj1PrU2JIkSdJi0vbU2Csj4sURcbeGNCdUaY5ts2xJkiRpGNpuYf5r4JXADQ1prgdeBfxVy2VLkiRJrWs7YH4A8JPMvLxfgmrZT4BTWi5bkiRJal3bAfNtgV8NkO5XVVpJkiRph9Z2wLw7sHWAdLcCe7VctiRJktS6tgPmXwMnDpDuROA3LZctSZIkta7tgPmrwOER8bR+CSLiqcARwDktly1JkiS1ru2A+UzgFuAtEfFP9aHjIuLYiPgn4K3AZuBNLXfoo9QAACAASURBVJctSZIkta7VgDkzLwGeCiTwXODCiNgSEVuAC6vXAJ6emT9vs2xJkiRpGFqf6S8zPwDcDzgL2FiVsQuwqXrtpMx8f9vlSpIkScOw2zAyzczvAadHxG7AiurltZm5ZRjlSZIkScMylIC5owqQrxpmGZIkSdIwtdolIyL2jYjfiYiDG9IcUqUZb7NsSZIkaRja7sP8bODrwGENaW5bpXlGy2VLkiRJrWs7YD4NuCwzv9svQbXsMuD0lsuWJEmSWtd2wHwkcPEA6S6mTF4iSZIk7dDaDpjHgBsHSHcjsF/LZUuSJEmtaztgvhY4boB0xwE3tFy2JEmS1Lq2A+ZvAydExP36JaiW3b1KK0mSJO3Q2g6Y3wEE8KmIeGj3wuq1T1Cmzv73lsuWJEmSWtfqxCWZ+ZWIeAfwdOCzEbGGqYsAjwMOogTU78rML7ZZtiRJkjQMrc/0l5l/HRG/AF4CHFzdOq4HXpuZb2y7XEmSJGkYhjI1dma+KSLeAtwTOLx6+Qrgu5l5yzDKlCRJkoZhKAEzQBUYf6u6bScidsvMLcMqX5IkSWpD2xf9TSsi7hARbwR+PeqyJUmSpJkaWgtzXUQsAx4L/AWlm0ZQRsqQJEmSdmhDDZgj4neBPwceCexNCZRvAb4AfGCYZUuSJEltaD1gjoiDgCdRAuVjKEEylBblZwEfzczr2i5XkiRJGoZWAuaI2BV4GKXLxUOAXSmB8g3AR4FTgGMz8y1tlCdJkiSNypwC5og4jtKS/ETgNkz1TT4HeA/w35m5KSK+Dhw7x7pKkiRJIzfXFuYLKQFyUMZZfh/wvsy8Yq4VkyRJknYEbfVh/jXwXOAzmXlrS3lKkiRJ826u4zB/kdLCfFvg48CvI+INEbFyzjWTJEmSdgBzCpgz8zTK1NcvB1YBBwLPAy6IiG9FxF9UYzBLkiRJC1Jktjd/SEQ8EPhL4I+BPSmtzxuArcCyzNy1tcJ2cBExBkxMTEwwNjY239WRJEmL3OTkJOPj4wDjmTk53/VZTFqdGjszz8nMxwGHUMZc/gmwFNgHICJWRcQrI+KYNsuVJEmShqXVgLkjM2/MzLdk5gnAicA7gEngdsBLgIsj4hvDKFuSJElq01AC5rrMPD8z/xo4mDID4HnVovsOu2xJkiRprlqfGrufzNwIfBD4YNUl489GVbYkSZI0WyMLmOsy85eUrhmSJEnSDm3oXTIkSZKkhcyAWZIkSWpgwCxJkiQ1MGCWJEmSGhgwS5IkSQ0MmCVJkqQGBsySJElSg6EEzBFxRES8JiLOjYifRcRra8vuGRF/HhFjwyhbkiRJalPrE5dExBOAdwJLgAAS+E4tyRjwLuBW4D/aLl+SJElqU6stzBFxL+B9lGD4JcD9KEFz3VeBSeD0NsuWJEmShqHtFuYXUgLkh2bmeQAR28bLmbk1In4E3LHlsiVJkqTWtd2H+STge51gucFVwMEtly1JkiS1ru2AeT/g8gHS7Qns0XLZkiRJUuvaDpivA243QLpjgDUtly1JkiS1ru2A+TvAiRGxsl+CiLgvcCfgmy2XLUmSJLWu7YD5bZQLCT8REXfqXhgRxwLvpQw19/aWy5YkSZJa12rAnJlfBt4MrAR+HBE/5/+zd+fxUVX3/8dfn+z7BAhrQNlUFHEHquICVFxRccG1ClatWq36RflVW5eqbbUVt9ZaXMEN+aqtorigUL9VUFFwg4oFDIIkEDBkQvbt/P6YmZidhNxkMsn7+XjM42buPXPOmUwm+eTOuZ9PIDieZGYrgdXAPsD9zrllbR3PzM4OFkfZYWZFZvaFmc00s1gP+j7JzFzw9m5b+xMRERGRyOR5pT/n3HXAL4FtwAgCaeYGAAcRyL98vXPuhraOY2YPAP9LINfzcuAtAuun7wGWmFliG/ruQaC4imvrPEVEREQksnle6Q/AOfeImT0KHAoMJRCYbwI+cs5VtLV/MzsduBYoBI5xzq0M7s8AlhBIb3cnsLuB+V+AvsDfgSvbOl8RERERiVyen2EOcc5VOeeWO+decM4975x734tgOejm4PbuULAcHHM7cFXw7tVm5mttx2Y2BbgAuI/AmWsRERER6cbaLWBuL2aWCYwO3n2+/nHn3AcEzmbHAye1su8MAmeVvwFubdtMRURERKQr8DRgNrMrzazczE5pps3kYJtLd3OYg4PbPOdcVhNtPq3XtqUeATKAnzvnSndnciIiIiLStXh9hnkKsAN4s5k2bwD5wJm7OcaQ4HZjM2021Wu7S2Z2LnAW8BfnnHJEi4iIiAjg/UV/+wKrnHNVTTVwzlWZ2VfAfrs5RmpwW9RMm8LgNq0lHZpZP+BhYD0/ro9uFTOLJ7AMJCS1qbYiIiIiEjm8PsPcG8hpQbstQB+Px26LR4EewKXOueLd7OMmwF/r9r1HcxMRERGRMPI6YN5JIOfyrvQDStowBkByM21SgtuCXXVmZhcDk4G/O+fe2805AfwR8NW6DWxDXyIiIiLSSXi9JONL4Agzy3TObW6sgZkNJFBs5OPdHGNDcDuomTahYxuaaRMyJbgdbWbv1TvWL7g9tNaxc51zW+p34pwrA8pC982sBUOLiIiISGfndcD8PDAeeNnMTnXO5dY+aGZ9gJeAWBpJCddCnwW3vcxsSBOZMg4Lblc2cqwphzVzLB04Jvh1Qiv6FBEREZEI5/WSjLnAR8AYYL2ZPW9mtwZvzxO4qG4MgYIgT+zOAM6574FPgnfPr3/czMYROMNcRiAjx676O905Z43dgOnBZotr7d+wO/MWERERkcjkacDsnKskUCzkNQJrjM8Fbgvezg3uWwicHGy7u/4Q3P7azA4J7TSzXsDfgnf/6pzz1zo2xczWmNniNowrIiIiIt2M10sycM7lA6cFA9kTgD2DhzYCbznnVngwxitm9hDwK+CjYBBcBEwksHxiKXBLvYf5gH3QkgoRERERaQXPA+YQ59xKWreGuLX9X2tmS4FfAkcQWBe9HrgbuN85V95eY4uIiIhI92HOuXDPoUsyszTA7/f7SUtrUf0UERERkd1WUFCAz+cD8DnndplaV1qu3c4wm1lfIJNmlkA455a11/giIiIiIl7wPGA2s1MIFPHYVelr1x7ji4iIiIh4ydOA1cyOB/4JRAOFQBYtqLYnIiIiItJZeX2G9zcEguW7gN8Hq9+JiIiIiEQsrwPmg4EvnHO3etyviIiIiEhYeF3pzwFfe9yniIiIiEjYeB0wf0kgM4aIiIiISJfgdcD8IHBk7XLVIiIiIiKRzNOA2Tn3IoFKe2+b2eVmNsDL/kVEREREOprXaeVC5aijgUeAR8ysmsDa5vqccy7ey/FFRERERLzmdZaMxvqL9ngMEREREZEO43XAHOtxfyIiIiIRp7rasTm/hKLySpLjYshMTyQqysI9LdlNngbMzrkqL/sTERERiTTrcnfy9qqtrN9WSGllFQkx0QzrncLx+/dleJ/UcE9PdoPXZ5hFREREuq11uTt5aukG8orK6e9LICkukeLySlZl+8n2lzD9yMEKmiOQ12nlRERERLql6mrH26u2kldUzl59UkhNiCU6ykhNiGWvPinkFZWzaPVWqqsby4UgnZnnAbOZxZjZtWb2gZltM7MyMytv5Fbm9dgiIiIi4bI5v4T12wrp70vArO56ZTOjvy+BdbmFbM4vCdMMZXd5nVYuDngHGAdoZbuIiIh0G0XllZRWVpEUl9jo8cS4aLYWlFJUXtnBM5O28voM86+Ao4DFwEjgGQI5mJOBg4BZQDlwJ8qoISIiIh2sutqxKa+YNVsK2JRX7OnyiOS4GBJioiluIiAuKa8iPiaa5DhdQhZpvH7FpgKFwFTnXH6waAnOuRLgS+BGM1sKvAx8AfzD4/FFREREGtXe2Ssy0xMZ1juFVdl+UuJj6izLcM6R4y9lVKaPzPTGz0BL5+X1GeZ9gI+cc/nB+w7AzGrGcc69AnwOXOPx2CIiIiKNCmWvWJXtJz0plqEZKaQnxbIq289TSzewLndnm8eIijKO378vPZPjWJtbyM7SCiqrq9lZWsHa3EJ6JscxaWRf5WOOQF4HzLHAtlr3S4NbX712a4ADPR5bREREpIGOzF4xvE8q048czP4DfOQXV7BhexH5xRWMyvQppVwE83pJxhagX6372cHtCODDWvv7o5LZIiIi0gFak71iUM+kNo83vE8qQ49NUaW/LsTrM8xrgOG17n9EIFvGDRb8CTWzccDRwH89HltERESkgR+zVzR+njAxLpqyyipPs1dERRmDeiYxol8ag3omKViOcF4HzG8Dg8zssOD9JQSC6NOBTWb2MYEMGgY84vHYIiIiIg0oe4W0ldcB8zxgOlAE4JyrBqYQCJoHAKMJLMX4m3PuSY/HFhEREWkglL0ix1+Kc3XXKYeyVwzvk6LsFdIkT/+Vcs7lAnPr7fsGGGlmI4GewDfBdiIiIiI1qqtdu6z7DWWvyPaXsDY3sJY5MS6akvIqcvylyl4hu2T1/9MSb5hZGuD3+/2kpaWFezoiIiKdWnvnSK4/RlllYBnG8D4pTBrp3RjhVFBQgM/nA/A55wrCPZ+uRIt1REREJKxCOZLzisrp70sgKS6R4vJKVmX7yfaXeJaOTdkrZHe1KWA2syOCX65wzpXVut8izrllbRlfREREIlv9HMmhtG+pCbGkxMewNreQRau3MjQjxbPlGV6kjpPupa1nmD8AqoH9CKSJ+4Bgdb8WcB6MLyIiIhGso3Mki+yOtgasywgEvsX17ouIiIjs0o85khvPUJEYF83WglJPcySLtFabAmbn3Ljm7ouIiIg0p3aO5NSE2AbHlSNZOgOv8zCLiIhIF1Nd7diUV8yaLQVsyiumutq7D5OVI1kigaf/rplZBfC5c260l/2KiIhIeLR3ujflSJZI4PXnG6UELv4TERGRCNeR6d6mHzm4JjDfWlBKfEw0ozJ9XSZHskQ2rwPmdUBvj/sUERGRDtbR6d6UI1k6M6/XMD8HHGVmQzzuV0RERDpQa9K9eSWUI3lEvzQG9UxSsCydhtcB8/3Au8ASM5tqZg0vdxUREZFO78d0b41/GJ0YF01ZZZXSvUm34PWSjK8JBOF7AvMAZ2Y5QGP/fjrn3D4ejy8iItJtVFe7dlvCoHRvIj/y+qd8eK2vLXjLbKKtCpyIiIjspvbOXhFK97Yq209KfEydZRmhdG+jMn1K9ybdgtcB814e9yciIiL1dET2CqV7E/mRpwGzc269l/2JiIhIXR2ZvULp3kQCtPBIREQkgrQme8WgnkltHk/p3kQUMIuIiESUH7NXNL52ODEumq0FpZ5mrwilexPprtolYDaz04CzgX2ANAIX/9WnLBkiIiKtpOwVIh3P03eTBT4beh6YSuNBMgSyYxjKkiEiItJqyl4h0vG8Llzyc+AcYDVwMvBPAoHxSOB04GUCwfJdwN4ejy0iItJpVFc7NuUVs2ZLAZvyiqmu9uY8USh7Rc/kONbmFrKztILK6mp2llawNrdQ2StE2oE5592JXjP7NzAGGOqcyzazp4CLnHPRtdpcBTwETHTO/Z9ng3cyZpYG+P1+P2lpaeGejoiIdKD2zpFcf4yyysAyjOF9UpS9ohsrKCjA5/MB+JxzBeGeT1fidcCcB3zhnBsfvP8kcDEQ42oNZGZfA9865072bPBORgGziEj31DBHcgzF5ZU1uYu9yJEc0p6V/iTyKGBuP15fEZAEZNe6XxbcpgK1X7jPgOM8HltERCSsOjJHMih7hUhH8XoN81Ygo9b93OB2eL12PQBdjSAiIl1Ka3Iki0jk8DpgXg8MrXX/UwIX+f0itMPM9gHGB9uKiIh0GT/mSG78A9zEuGjKKqs8zZEsIu3P64D5HWCome0bvP82sBm41Mw+NLP5wIdALPCMx2OLiIiEVe0cyY1RjmSRyOT1O3YeEAekADjnys3sHOAVYGzwBrAQeMDjsUVERMJKOZJFuiZPA2bn3Abgd/X2LTOzIcCxQE/ga+fcp16OKyIi0lLtmVkilCM521/C2tzAWubEuGhKyqtqsmQoR7JI5PE0rZz8SGnlREQ6n47Ij1x/HOVIlo6itHLtp01nmM1sFfA48Kxzbrs3UxIREfFew/zIiRSXV7Iq20+2v8TT/MjD+6Qy9NgU5UgW6SLaetHffsAs4Hsze9HMTrD6eXRERETCrH5+5NSEWKKjjNSEWPbqk0JeUTmLVm/1rHw1/JgjeUS/NAb1TFKwLBLB2howzwBWEbjQ70wCF/NtNLM7zGxos48UERHpIMqPLCJt0aaA2Tl3v3PuQGAM8HfAD2QCvwHWmtliMzvfzOLbPlUREZHdo/zIItIWnuRhds596py7CugPXAgsARyBAiXPADlm9rCZHerFeCIi0jVVVzs25RWzZksBm/KKPVsiofzIItIWXqeVKwOeB543s0HAJcDFwGDgCuAKM/uKwIWCzznndng5voiIRK72zGCh/Mgi0hZeV/qr4Zzb5Jz7nXNuKPBTAkVNSoEDgAcJVAAUERGpyWCxKttPelIsQzNSSE+KZVW2n6eWbmBd7s429R/Kj9wzOY61uYXsLK2gsrqanaUVrM0tVH5kEWlWuwXMtTnnljjnLgSmAtsAA7SuWUREOiyDxfA+qUw/cjD7D/CRX1zBhu1F5BdXMCrT52lKORHpetp9sZaZ9QcuAqYDexEIlgG+au+xRUSk82tNBotBPZPaNJbyI4vI7miXgNnMYoBTCaxhngREEwiUCwgszXhC5bFFRARqZ7BofP1wYlw0WwtKPctgEcqPLCLSUp4GzGa2P/Bz4AKgFz+eTX4feAJ40TmnJJciIlKjdgaL1ITYBseVwUJEwq3Nv33MzEcgQJ4OHBLaDWwB5gJPOufWtnUcERHpmpTBQkQ6uzYFzGb2PHA6gQv4DKgC3iBwNnmhc66qzTMUEZGwq6527bbuN5TBIttfwtrcwFrmxLhoSsqryPGXKoOFiISdObf7Vx2bWXXwy3XAk8Ac59wWLyYW6cwsDfD7/X7S0tLCPR0Rkd3WnvmRmxqnrDKwDGN4nxQmjfR2HJGuqqCgAJ/PB+BzzhWEez5dSVuXZDxD4AK+f3sxGRER6VxC+ZHzisrp70sgKS6R4vJKVmX7yfaXeJqOTRksRKSzalPA7Jy72KuJiIhI51I/P3JobXFqQiwp8TGszS1k0eqtDM1I8XR5hjJYiEhn0yGFS0REJPK0Jj+yiEhXpoBZREQa9WN+5MY/jEyMi6asssqz/MgiIp2VklqKiESw9sxeofzIIiIB+i0nIhKh2jt7hfIji4gEKGAWEYlAHZG9QvmRRUQCtIZZRCTC1M9ekZoQS3SUkZoQy159UsgrKmfR6q1UV+9+nv2Q4X1SmX7kYPYf4CO/uIIN24vIL65gVKbP05RyIiKdmc4wi4hEmNZkr/AiRZvyI4tId6eAWUQkwvyYvaLxtcOJcdFsLSj1NHuF8iOLSHemJRkiIhGmdvaKxih7hYiItyI6YDazs83sPTPbYWZFZvaFmc00s4b5j5rv52Azu8nMFpvZVjOrCPb5vpn9srX9iYi0p1D2ihx/Kc7VXaccyl4xvE+KsleIiHgkYk8/mNkDwLVAJbAEKAQmAPcAk81sknNul+WnzCwGWBm8Wwh8AmwFBgKHA+OAi8zseOdcvudPRESklZS9QkSkY1n9sxORwMxOB/5JIMA9xjm3Mrg/g0DwPAqY5Zy7oQV9xQAfEQi0FzjnymodGwW8DfQHnnLOXdKKOaYBfr/fT1paWoufm4h0He1ZVATq5mEuqwwswxjeJ4VJI73JwywikaWgoACfzwfgc84VhHs+XUmkBszLgdHAb51zv693bBzwPlAG9HXO+ds41oXAM0AJgR/AihY+TgGzSDfW3kVFQto7KBeRyKGAuf1E3JIMM8skECwDPF//uHPuAzPbBAwCTgLmtXHIz4LbRCADyGljfyLSxXVEUZEQZa8QEWl/kXjR38HBbZ5zLquJNp/Wa9sWewW35UCeB/2JSBfWkUVFRESkY0RiwDwkuN3YTJtN9druFgtUBJgZvPt67fXNjbSNN7O00A3QAkKRbqg1RUVERCQyRGLAHApEi5ppUxjctnXx8G0EMmUUAr/eRdubAH+t2/dtHFtEItCPRUUaX/GWGBdNWWWVp0VFRESkfUXcGuaOYmYXAbcC1cAlzrm1u3jIH4H7at1PRUGzSLdTu6hIakLDFO4qKrJrzjkqKiqorq4O91REwiIqKorY2NgGn1JJ+ETib+ydwW1yM21SgtvdukLUzM4Gngzevcw59+KuHhNcrlE7Jd3uDC0iES5UVGRVtp+U+Jg6vwtCRUVGZfpUVKQR5eXl5ObmUlxcTFVVVbinIxJW0dHRJCUl0adPH+Li4sI9nW4vEgPmDcHtoGbahI5taKZNo8zsDALZN6KAXzjnntzFQ0REaqioyO4pLi5m06ZNREdH06NHDxITE4mOjtbJB+l2nHNUVVVRUlKC3+9nw4YNDBw4kKQkZcMJp0gMmENp3nqZ2ZAmMmUcFtyubORYk4IFUV4AooErnXOP7f40RaS7Gt4nlelHDq7Jw7y1oJT4mGhGZfpUVKQJ27dvJzY2lj333JPo6OhwT0ck7FJSUujZsyffffcd27dvZ4899gj3lLq1iAuYnXPfm9knBHIxnw80VrhkEIHlEW+0tF8zmwz8L4HvyZXOudmeTVpEOpWOKPYxvE8qQ49NUVGRFqisrKSoqIj+/fsrWBapJTo6mp49e5KTk0NlZSUxMREXtnUZkfqd/wOB0ti/NrM3a5XG7gX8Ldjmr7Wr/JnZFAIX5m12zk2s3ZmZnQS8ROD7cYVz7tEOeA4iEgYdVYEPVFSkpSorAxlD4uPjwzwTkc4n9L5QwBxeEfmdd869YmYPAb8CPjKzxQTSzE0E0oGlwC31HuYD9gESau80sz7AP4A4AlktjjCzI5oY+gbn3HbPnoiIdKiOrMAnraf1yiIN6X3ROURkwAzgnLvWzJYCvwSOAGKB9cDdwP3OufIWdpUEhE5rDAQubqbt7YACZpEIVL8CX+iPUGpCLCnxMazNLWTR6q0MzUjRsgkREakjYgNmAOfc/xJYd9yStnOAOY3s3wDor6NIF9eaCnxaRiEiIrVFYqU/EZFWUwU+ERHZXQqYRaRbqF2BrzGqwCed3ZIlS4iKiiI5OZl169Y12e53v/sdZsbIkSMpKytrsl1n9tvf/hYz46677gr3VEQABcwi0olUVzs25RWzZksBm/KKqa52nvUdqsCX4y/Fubr9hirwDe+Togp80mlNmDCBq666iuLiYqZNm9Zo6fCVK1dy1113ERMTw9NPP63MIyIe0akUEekU2jvdmyrwSVdwzz338NZbb7F06VJmzZrFjTfeWHOsrKyMiy66iMrKSm699VYOPfTQMM60ba677jouvPBCevfuHe6piAA6wywinUAo3duqbD/pSbEMzUghPSmWVdl+nlq6gXW5Oz0ZJ1SBb/8BPvKLK9iwvYj84gpGZfqUUk4iQnJyMnPmzCEqKopbbrmF//znPzXHbrnlFlavXs0hhxzCb3/72zDOsu0yMjIYMWIEvXr1CvdURAAFzCISZvXTvaUmxBIdZaQmxLJXnxTyispZtHqrZ8szhvdJ5cpjh3H9cXtzzcS9uP64vbnimGEKliVijBs3juuvv56ysjIuvvhiKisrWbZsGbNmzSIuLo65c+cSGxvb4v5Gjx6NmfHSSy812eaBBx7AzDj77LNr9hUUFPDoo48yZcoUhg8fTnJyMikpKRxwwAHccsst+P3+RvsaOHAgZsb333/PP//5T8aPH0/Pnj0xMz744AOg6TXM5eXlPPPMM5x33nnss88+pKamkpiYyIgRI7juuuvIyclp8nsW6n/lypWcfvrpZGRkkJCQwMiRI3nggQcaLNWq7d133+Wss84iMzOT+Ph4evfuzZgxY7j99tvJy8tr0P6bb77hsssuY+jQoSQkJJCens4xxxzD888/3+QY0rkpYBaRsGpNujevhCrwjeiXxqCeSVqG0U2155r59nbXXXex77778umnn3LLLbfUrGn+3e9+x/7779+qvqZPnw7AnDlzmmzz1FNPAXDJJZfU7Fu5ciW/+MUvWLZsGf3792fy5MkceeSRbN68mbvuuouxY8eyY8eOJvu85557OOOMMygqKuKEE07g6KOP3mVp9OzsbC666CLeeustevbsyYknnsiECRMoKCjgwQcf5OCDDyYrK6vJx7/xxhv85Cc/Ye3atRx33HGMHTuWb775huuvv77O8pYQ5xxXXXUVxx13HC+//DJ9+/bljDPO4LDDDmP79u387ne/q3OWH+CFF17gwAMP5PHHHyc+Pp6TTz6ZQw45hE8//ZQLLriAyy+/vNnnKJ2Uc063drgBaYDz+/1ORJr2dY7fXfnsp+7Pb61x9y36psHtT2997a569lP3dY7eS11VSUmJ+89//uNKSko6bMy1WwvcXxevdde/8Jm78tlP3fUvfOb+unitW7u1oMPm0FbLly930dHRDnCAO/zww11lZWWr+9mxY4dLSEhwMTExLicnp8HxFStWOMANGDCgTv/fffedW7JkiauqqqrTvrCw0F1wwQUOcL/61a8a9JeZmekAFxMT415//fVG5/Sb3/zGAe7OO++ssz8/P9+99tprrry8vM7+8vJyN3PmTAe4U089tUF/Rx55ZM336fHHH69z7O2333Zm5mJiYlx2dnadY/fdd58DXEZGhvu///u/Bv1++OGH7vvvv6+5/9lnn7m4uDiXmJjoXnnllTpts7Ky3MiRIx3gnnvuuUafd2Na8/7w+/2h55nmOkEs1JVuOsMsImGldG/S0TpqzXx7Gz16NMcdd1zN/dmzZ+/yDG1j0tPTmTJlCpWVlTzzzDMNjofOLl900UV1+t9jjz0YP348UVF1Q4nk5GQeeeQRoqKiePHFF5sc95JLLuHk9phxhAAAIABJREFUk09u1Vx9Ph+nnHJKgyUnsbGx3H333fTp04c33niD4uLiRh8/depUfv7zn9fZN2nSJCZOnEhlZSXvvfdezf7y8vKaJSFPPPEERx99dIP+fvKTn5CZmVlz/6677qK8vJw//vGPnHbaaXXaDh48mMceewyAhx56qOVPWjoF/QUSkbAKpXtble0nJT6mzrIM5wLp3kZl+pTuTTzRlUqkL168mLfffrvm/rx58xg1atRu9TV9+nTmzZvHnDlzGmTeCK27DS3dqM05x4cffsi///1vNm3aRHFxcehTVuLj48nJyWHnzp2kpja8RuCss87arbkCfP755yxevJgNGzZQVFRUk2LPOUdlZSXr169v9HsxefLkRvvbd999effdd9m8eXPNvk8++YS8vDz69evHqaeeuss5VVVV1bweU6dObbTN2LFjSUxM5NNPP6WioqJVa80lvBQwi0iLVFc7NueXUFReSXJcDJnpiZ4EFEr3Jh2pq5RILygo4JJLLsE5x9VXX83s2bP505/+xJQpUxg9enSdto8++ijLli1r0MdvfvMb9tprLwAmTpzIHnvswX/+8x+WL1/OmDFjAFiwYAF5eXkceeSR7L333nUev2XLFs4888xG+64/18YC5sGDB7fmKQNQWFjIBRdcwIIFC3Y5ZmP22GOPRvenpaUBUFpaWrPvu+++A2CfffZp0dxyc3MpLCwEYMCAAbtsn5eXR9++fVvUt4SfAmYR2aX2zpEcSvcWGmNrQSnxMdGMyvQxaaQ3Y4hA7RLpjX9ikRgXzdaC0k5fIv36669n48aNTJw4kYceeoiMjAxuv/12pk2bxsqVK+sULPn3v//Nc88916CPSy+9tCZgjoqKYtq0adxxxx3MmTOnJmAOLcdo7OzyJZdcwrJlyxg3bhy33XYbBxxwAD169Kg5a9qnTx+2bdtWc8a5vsTE1n9qNHPmTBYsWMB+++3HH//4Rw477DAyMjKIi4sDYMyYMXzyySdNjll/+YiXQme5zYyLLrpol+1Dc5bIoIBZRJoVWu+ZV1ROf18CSXGJFJdXsirbT7a/xLP8xcP7pDL02JR2OYstElJ7zXxqQsOPwyNhzfzChQt58sknSUtL48knn8TMuPnmm3n11Vf57LPPuO2227j77rtr2j/77LM8++yzu+x32rRp3HnnncybN4/777+fvLw8Fi1aRHJyMuecc06dtgUFBbz99tvExMSwcOHCmjO0IX6/n23btnnzhGv53//9XwBefPFF9ttvvwbHmysZ3lqhs9HffPNNi9r36dOH+Ph4ysrKeOCBB0hPT/dsLhJ+uuhPRJrU0TmSle5N2lukl0jPy8vjsssuA+C+++6rCepiY2OZM2cOsbGx3HvvvXz88cet7nvIkCEce+yx5Ofn88orr/D0009TVVXFWWedRUpKSp22+fn5VFdXk56e3iBYBloUoO+OUJq6Pffcs8GxhQsXNpvGrrXGjBlDz5492bJlC6+99tou28fGxjJx4kTgx8Beug4FzCLSpHDkSBZpT6E18z2T41ibW8jO0goqq6vZWVrB2tzCTr9m/uqrryYnJ4cTTzyxQbaHUMGQqqoqpk2bVmc9bkuF8iw/9dRTNXmZG1uOMWDAAHw+H9u3b2fevHl1ji1durTdKg2OGDECgL/85S919n/99ddcddVVno4VFxfHTTfdBASWr4SKqtT28ccf17lQ8Pbbbyc2NpYZM2bw7LPP1izTqO2rr77ilVde8XSu0v4UMItIk35c79n4x9OJcdGUVVZ1+vWeIrVFaon0l19+mXnz5pGenl6Tnqy+m266iUMOOYQ1a9Zw6623tnqMM888k7S0NN5++23WrFnDsGHDGk2nFhMTUxMUn3/++Rx++OGcf/75HHnkkRx11FFMmTKlTro1r9x+++1A4HkedNBBnHfeeUyYMIEDDjiAffbZh7Fjx3o63owZM7j00kvJzc3lqKOO4tBDD+X888/npJNOYtiwYfzkJz+pUyhl9OjRNWfmf/aznzF48GBOOOEELrzwQk488UQGDhzIAQcc0GxVRemcFDCLSJOUI1m6qkgrkZ6bm8sVV1wBwIMPPthkMBoTE8PcuXOJi4tj1qxZfPTRR60aJzExkXPPPbfm/rRp0xp8uhRyww038I9//IPDDz+cNWvW8Nprr1FRUcEjjzzCE0880apxW+rss8/mvffeY8KECWzevJkFCxawbds27rzzThYuXEhMjLe/i8yMxx57jIULF3Lqqafy/fff89JLL7FixQoyMjK44447GDlyZJ3HnHvuuaxevZrrrruOtLQ03n//fV5++WXWrFnD3nvvzT333MMdd9zh6Tyl/VlTV5JK25hZGuD3+/2Nru8SiQTV1Y5H3lvPqmx/nZy1EFjvuTa3kFGZPq44Zlin/QhbOr/S0lKysrIYMmQICQkJ4Z6OSKfSmvdHQUEBPp8PwOecazy3nuwWnWEWkSZF+npPERERLyhgFpFmRep6TxEREa9o4aFIhGuvCny1KUeyiIh0ZwqYRSJYe1fgqy2UI1lERKS7UcAsEqE6qgKfiIhId6c1zCIRqKMr8ImIiHRnCphFIpAq8ImIiHQcBcwiEUgV+ERERDqOAmaRCKQKfCIiIh1HAbNIBMpMT2RY7xRy/KXUr9bpnCPHX8rwPilkpieGaYYiIiJdhwJmkQikCnwiIiIdRwGzSIRSBT4REZGOoQWOIhFMFfhERETanwJmkXbSESWrQRX4RERE2puWZIi0g3W5O3nkvfXc/85/eWjxWu5/57888t561uXuDPfURCTCXHjhhZgZJ510UrPtduzYQWZmJmbG448/3kGz805paSlmRkJCQrinItKAAmYRj4VKVq/K9pOeFMvQjBTSk2JZle3nqaUbFDSLSKv89a9/JTMzkzfffJPZs2c32e6Xv/wl2dnZnHzyyVx66aUdOEORrk8Bs4iHVLJaRLyWnp7OE088AcANN9zAt99+26DNSy+9xLx58+jVq1dEnl0GiI+P5+uvv+bLL78M91REGlDALOIhlawWkfZw/PHHc8UVV1BYWMjFF19MdXV1zbGtW7dy5ZVXAvC3v/2Nfv36hWuabWJmjBgxgr333jvcUxFpQAGziIdUslpE2su9997LsGHD+OCDD7j33ntr9l9++eVs376d8847j6lTp7aor1dffRUz46CDDmqyTVlZGb169cLM+Oabb2r2L1u2jBtuuIHDDjuMvn37EhcXR79+/TjttNN47733Gu3r73//O2bGFVdcwbZt27jmmmsYMmQIcXFxnHDCCUDza5i/+uorbrnlFo444ggGDBhAXFwcvXv3ZtKkSfzzn/9sdMy33noLM+OEE06gvLyc3//+9+y7774kJCSQkZHB2Wefzdq1a5t8/j/88AO33XYbhxxyCGlpaSQlJTFs2DDOPfdcFi1a1KB9RUUFs2fP5uijj6ZHjx7Ex8czdOhQrr76arKzs5scRyKDsmSIeKh2yerUhNgGx1WyWqQTqa4G/yYoL4S4FPANgqjOex4pOTmZuXPncvTRR3Prrbdy4oknsmLFChYsWMCAAQN4+OGHW9zXySefTJ8+ffjiiy/4/PPPGw2cFyxYQF5eHkcccQT77LNPzf6ZM2fy0UcfMXLkSEaPHk1iYiLr1q1jwYIFvPbaazzyyCP84he/aHTcLVu2cOihh1JUVMRRRx3FYYcdRmrqrnPG33PPPTz//POMGDGCAw88kLS0NDZu3Mi7777LO++8w0033cQf/vCHRh9bVlbGcccdx4oVKzjmmGPYb7/9+Pjjj3nppZf417/+xRdffEFmZmadx6xYsYJTTjmFLVu20KNHD44++miSk5PZuHEjr776Kvn5+UyaNKmmfX5+PqeccgpLly4lLS2NQw45hF69evHll1/y8MMP8+KLL7J48WL233//XT5X6aScc7q1ww1IA5zf73fSfVRVVbu/Ll7rrnjmUzfr7TXuvkXf1Nxmvb3GXfHMp+7hJWtdVVV1uKcq0mmUlJS4//znP66kpKTjBs1d49z//dm5ly93bv5Fge3//Tmwv5ObOXOmA9zIkSOdz+dzgHvjjTda3c+MGTMc4K699tpGj5944okOcI899lid/a+//rrbsmVLg/bvvfeeS05OdgkJCW7r1q11jj3yyCMOcIA78cQT3c6dOxs8vqSkxAEuPj6+wbHFixe7rKysBvu/+uor169fPwe4L774os6xN998s2bMMWPGuNzc3JpjRUVFbvz48Q5wv/rVr+o8Lj8/3/Xv398B7tJLL3VFRUV1ju/YscMtXry4zr4zzjjDAW7KlClu+/btNfurq6vd3XffXfN6VVVVNXgOu9Ka94ff7w895zTXCWKhrnTrvP9Ki0QglawWiQDbvoGP/g45X0JST+i1V2Cb82Vg/7Zvdt1HGN1xxx2MGjWK1atX4/f7ufzyyznxxBNb3c8ll1wCwHPPPUdFRUWdY9nZ2SxatIikpCTOOeecOsdOPvlk+vbt26C/Y445hssvv5zS0lJef/31RseMj49n9uzZpKSktGquEyZMYPDgwQ3277///tx0001A4MLHxkRHRzNnzhx69+5dsy8pKYnbbrsNgHfffbdO+7///e/k5OQwduxYZs+eTVJS3Tz36enpTJgwoeb+559/zj/+8Q/23HNPnn32WXr16lVzzMz4f//v/zFhwgRWr17NkiVLWvW8pfPQ58LSLbVnUZFQyeq3V21l/bZCthaUEh8TzahMH5NG9lXJapFwqq6Gr1+D4h+g9wgIXZwbnwa9U2HbGljzeiCI7qTLM+Lj4/nDH/7A5MmTAZg1a9Zu9bPffvsxZswYli9fzuuvv86UKVNqjj399NNUVVVx1llnNbpkYvv27SxcuJBVq1axY8cOKisD12WsWbMGoM6a59rGjh3LoEGDdmu+BQUFvPHGG3zxxRf88MMPlJeXA/D99983O+bw4cPZd999G+wP7du8eXOd/W+99RYAl156KVEt+BlYuHAhAKecckqD4Drk2GOPZcmSJSxbtoyf/vSnu+xTOh8FzNLtrMvdWRPMllZWkRATzbDeKRy/v3fBrEpWi3RS/k2wfS34Mn8MlkPMIC0Ttv030K7HnuGZYwvUPkPb2NnanJycmjOvtY0aNYoZM2bU3L/kkktYvnw5Tz31VJ2Aec6cOQBMnz69QR8PP/wwN954IyUlTWf7KSgoaHR/Y2eJW+Lll1/msssuY8eOHa0ec4899mh0f1paGhC42LC27777DoARI0a0aG6hNH8PP/zwLteRb9u2rUV9SuejgFm6lVBRkbyicvr7EkiKS6S4vJJV2X6y/SVMP3KwZ0GzSlaLdELlhVBZCrHJjR+PS4Kd2YF2Eczv9zN37twG+48//vg6AfO5557L9ddfz5tvvsnWrVvp27cvy5Yt45tvvmHo0KEcc8wxdR6/dOlSrr76auLi4pg1axYnnXQSgwYNIikpCTPjoYce4tprrw1dy9NAYmJiq59LVlYW559/PhUVFfzmN7/hnHPOYfDgwSQnJxMVFcWCBQs47bTTmhyzJWeJ2yKU4u/QQw/d5UV9hx12WLvORdqPAmbpNuoXFQnlSU5NiCUlPoa1uYUsWr2VoRkpOhMs0lXFpUBMAlQUBZZh1FdeHDge17o1tp3NiBEjmgwga/P5fJxxxhk899xzPPvss8yYMaPm7PK0adMa5JN/8cUXAZgxYwb/8z//06C/5tK07a5XX32V8vJyzjvvPO666652H3OPPfYgKyuLNWvWMG7cuF22Dy0xOfbYY+uk+5OupXMu0BJpByoqIiL4BkHGXuDfDPUDSuegYDP03jvQrpsIXfw3Z84ciouLmT9/PlFRUVx88cUN2ubl5QGw554Nl6sUFRXxyiuveD6/5sasrq5m3rx5no4Xygv9xBNPtOifjtAFl6+88kqDiyel61DALN2GioqICFFRsO9kSOoVuMCvtACqKwPbbWsguReMOKXTXvDXHsaPH8/gwYNZtWoVN998MwUFBUyYMKHRtb+hC+WefPJJioqKavaXlJRw+eWX11yA56XQmPPnzyc3N7dmf2VlJb/+9a9ZsWKFp+NdccUV9O3bl48++ogrr7yS4uLiOsfz8/PrZLs4/PDDOfnkk1m/fj1nnXUWmzZtatBnYWEhTz/9ND/88IOnc5WO031+I0i3V7uoSGNUVESkm+i9D/zkCuh/AJTkwQ/rAtsBB8LYKwLHuxEzY9q0aQA8+OCDwI9nneu77LLLyMzMZPny5QwZMoQpU6Zw5plnsueee7Jw4UKuueYaz+d35plncsABB5CVlcXw4cOZPHkyU6dOZciQITz44IPMnDnT0/HS09NZsGABvXv3Zvbs2QwaNIhTTz2Vc889lyOOOIL+/fvzpz/9qc5jnn32WY466igWLFjAXnvtxdixYznnnHM4++yzGT16NL169eLiiy9m586dns5VOo4CZuk2MtMTGdY7hRx/aYOP2Zxz5PhLGd4nhcz01l+UIiIRpvc+MO5/YPzNcMzMwPbI67tdsBwybdq0movj0tPT62TMqC0jI4MVK1Zw+eWXk5KSwptvvsnHH3/MSSedxOeff85+++3n+dzi4uJ4//33mTlzJv379+edd97h3//+N2PHjmX58uWMHz/e8zHHjBnDqlWr+PWvf82AAQNYvHgxr732Grm5uZxxxhnceOONddqnp6fzr3/9i6effprx48eTlZXFP/7xD9577z3Kysq46KKLWLBgAQMHDvR8rtIxrCXrc6T1zCwN8Pv9/prUNRJ+9bNkJMZFU1JeRY6/lJ7JcZ5myRCRliktLSUrK4shQ4aQkJAQ7umIdCqteX8UFBTg8/kAfM65xvPsyW7RGWbpVkJFRfYf4CO/uIIN24vIL65gVKZPwbKIiIg0Sos1pVNpzwp8ISoqIiIiIq2hgFk6jY6owBeioiIiIiLSUgqYpVPoyAp8IiIiIq2hNcwSdvUr8KUmxBIdZaQmxLJXnxTyispZtHor1dW6QFVEREQ6ngJmCTtV4BMREZHOTAGzhJ0q8ImIiEhnpoBZwk4V+ERERKQzU8AsYacKfCIiItKZKWCWsIuKMo7fvy89k+NYm1vIztIKKqur2VlawdrcQnomxzFpZF/lSRYREZGwUMAsnYIq8ImIiEhnpUWh0iKqwCciIiLdlQJm2SVV4BMREZHuTAGzNEsV+ERERKS70xpmaZIq8ImIdE7OOebPn88ZZ5zBoEGDSEhIoEePHhx00EHMnDmTjRs3Nvq4jRs3Mnv2bM444wz23HNP4uPjSUlJ4cADD+Tmm29m27Ztddpv2bKFjIwMzIz58+c3O6f58+djZvTu3ZstW7Z49lxFOgOrn8ZLvGFmaYDf7/eTlpYW7unslk15xdz/zn9JT4olNSG2wfGdpRXkF1dw/XF7axmFiOy20tJSsrKyGDJkCAkJCeGeTqeXnZ3NlClTWL58OWbGoYceyrBhwyguLubDDz9k+/btxMfHM2vWLH75y1/Weey4ceNYunQpMTExHHzwwQwdOpS8vDw+/vhjCgoK6N27N4sWLeKggw6qecyLL77I1KlT6dmzJ6tWraJ///4N5pSTk8P+++9PXl4eL774ImeddVa7fx+6i9a8PwoKCvD5fAA+51xBh0ywm9AZZmmSKvCJiHQuO3bs4KijjmL58uUcfPDBfPXVV3zyySe88MILLFiwgJycHO69914qKyu5+uqreeihh+o8PjMzk/vvv58tW7awfPlyXnjhBRYtWsS6des49thj2bZtG1OnTqWqqqrmMWeffTbnnXceeXl5/PznP290Xpdddhl5eXmcf/75CpalS1LALE1SBT4Rkc7l6quv5ttvv2XIkCEsWbKEkSNH1jkeExPDjBkzePDBBwG44YYbWLNmTc3x+fPnc91119GrV686j+vduzfPPPMMAGvXruXDDz+sc/zhhx9mwIABvPnmmzz66KN1jj3xxBMsXLiQzMxM/vrXv3r2XEU6EwXM0iRV4BMR6Ty+/fZbXnjhBQDuvfde0tPTm2x71VVXceCBB1JRUcGf//znFvU/cOBAMjIyANi0aVOdYz169OCJJ54AYMaMGXz77bcAfPfdd1x//fVAIHDu0aNH656USIRQwCxNUgU+EenKql01mws3898d/2Vz4WaqXXW4p9Ss1157jerqatLT0zn11FObbWtm/OxnPwNgwYIFDU56NGb79u3s2LEDoNF1yieccAK/+MUvKCws5OKLL6ayspLp06ezc+dOrrjiCo4//vjdeFYikUGfpUuzQhX4QnmYtxaUEh8TzahMH5NGep+HWUSkI3yb/y2LNy4my59FWVUZ8dHxDPENYeIeExmaPjTc02vUihUrADj44IOJidn1n+/Ro0cDgUD4u+++Y/Dgwc22v/fee6mqqqJ///4cccQRTbZ55513+OCDDzjmmGNYtmwZw4YN4957723dkxGJMAqYZZdUgU9EupJv87/lua+fY0fZDvol9SMxJpGSyhK+zvuaLUVbuGDfCzpl0BxK+da3b98Wta/dbtu2bc0GzO+++25N0Dtr1izi4uIabZeSksLcuXNrguWoqCjmzp1LcnJyC5+FSGTSkowuoLrasSmvmDVbCtiUV9wueZFDFfhG9EtjUM8kBcsiEpGqXTWLNy5mR9kOhvmGkRKXQnRUNClxKQzzDWNH2Q4Wb1zc6ZdntETtZRi1s17U99VXX3H22WdTVVXFNddcw3nnnddsv+PGjWPKlCkATJkyhSOPPNKbCYt0YjrDHOE6smy1iEikyynKIcufRb+kfpjV/cffzOiX1I8sfxY5RTlkpmSGaZaNC12Qt3Xr1ha1z83Nrfm6d+/ejbZZs2YNP/3pT8nPz2f69Ok12TV2JSUlpc5WpKvTGeYIFipbvSrbT3pSLEMzUkhPimVVtp+nlm5gXe7OcE9RRKRTKaoooqyqjMSYxrP7JMYkUlZVRlFFUQfPbNcOPfRQAFauXEll5a7z3y9fvhwAn8/HkCFDGhz/73//y4QJE8jNzeWiiy7i8ccfb/BPhIgEKGCOUCpbLSLSesmxycRHx1NSWdLo8ZLKEuKj40mO7XxrcidPnkxUVBR+v59XX3212bbOuZq8yqeddhpRUXX/3K9du5bx48eTk5PDhRdeyFNPPdWgjYj8SO+OCLU5v4T12wrp70to9GPF/r4E1uUWsjm/8T8KIiLdUf/k/gzxDWFL8ZZG88tvKd7CEN8Q+ic3TKsWbsOGDWPq1KkA3HjjjeTn5zfZ9m9/+xtffvklcXFxzJw5s86x9evXM378eLKzs7nwwguZO3eugmWRXdA7JEKpbLWISOtFWRQT95hIj/gerPevp7C8kKrqKgrLC1nvX0+P+B5M3GMiUdY5/zw+/PDDDB48mKysLCZMmMDq1avrHK+srOS+++7j2muvBeDRRx+tUw0wKyuL8ePHs3nzZn72s58pWBZpIV30F6Fql61OTYhtcFxlq0VEGjc0fSgX7HtBTR7mrcVbiY+OZ9+e+3bqPMwAPXv25P333+f0009nxYoVjBo1isMOO4xhw4ZRXFzMhx9+yLZt20hLS+PPf/4zF198cZ3Hn3nmmWzatIn4+HgALrnkkkbHufTSSxk3bly7Px+RSKFoKkKFylavyvaTEh9TZ1lGqGz1qEyfylaLiDRiaPpQBvsGk1OUQ1FFEcmxyfRP7t9pzyzXNnDgQJYvX878+fOZP38+n3zyCZ9//jkVFRUAJCUlsXLlSoYNG9bgsXl5eQCUlZXVrHFuzLHHHquAWaQWa0m5TGk9M0sD/H6/n7S0tHYZI5QlI6+onP6+BBLjoikpryLHX0rP5DimHzlYqeVEpNMrLS0lKyuLIUOGkJCQEO7pRCy/38/48eP57LPPmDRpEgsWLKg5kyyRqzXvj4KCAnw+H4DPOVfQIRPsJjr/v9JNMLOzzew9M9thZkVm9oWZzTSzhusTWtbfoWb2opltNbNSM8sys7+YWR+v5+6VUNnq/Qf4yC+uYMP2IvKLKxiV6VOwLCLSzfh8Pt5++2323XdfFi1axDnnnNOi9HMismsRuSTDzB4ArgUqgSVAITABuAeYbGaTnHMtTg9hZmcB8wh8Pz4BsoDDgKuBs81snHNunbfPwhsqWy0iIiG9e/fm3Xff5bHHHsM5x4oVKxg7dmy4pyUS8SIuYDaz0wkEy4XAMc65lcH9GQSC53HAncANLexvADCXwPfiF865R4P7o4E5wIXA82Y21nXS9SuhstUiIiIDBgzgtttuC/c0RLqUSFyScXNwe3coWAZwzm0HrgrevdrMfC3s7zogCXg3FCwH+6sCrgT8wGhgUlsnLiIiIiKRJ6ICZjPLJBC8Ajxf/7hz7gNgExAPnNTCbqc0018hsCB494xWTVZEREREuoSICpiBg4PbPOdcVhNtPq3XtklmlgoMr/e43e5PRERERLqeSAuYhwS3G5tps6le2+YMrvV1U322pj8RERER6WIi7aK/UJ60ombaFAa3LUl+XDvvWlN9tqg/M4snsBSksb5FRGQXOul11SJhpfdF5xBpZ5g7s5sIXCAYun0f3umIiESGqKjAn6Kqqqowz0Sk8wm9L0LvEwmPSPvu7wxuk5tpkxLctqTCzc5aXzfVZ0v7+yPgq3Ub2ILxRUS6vdjYWGJjYyksLNx1Y5FuZufOnTXvEQmfSAuYNwS3g5ppEzq2oZk2Id/V+nqPtvTnnCtzzhWEbtQNxkVEpAlmRmpqKn6/n5KSFtecEunySkpKKCgoIDU1FTMVJAunSFvD/Flw28vMhjSRKeOw4HZlI8fqcM4VmNk6ApkyDgO+akt/IiKyezIyMigpKWHjxo2kpaWRmppKdHS0ggTpdpxzVFVVsXPnTgoKCoiPjycjIyPc0+r2Iipgds59b2afEMjFfD7w+9rHzWwcgTPCZcAbLez2n8CNwf6eqtdfCjA5ePcfuz9zERFpTnR0NINgvcxhAAAQ/ElEQVQGDWL79u3s3LmT/Pz8cE9JJKxiY2NJT08nIyOD6OjocE+n27NIu/oyWBr7nzQsjd0L+BcwCpjlnLuh1mOmEFhjvNk5N7FefwOAtQSq/V3unHssuD+aQAD9M+AToFWlsc0sDfD7/X7S0lqSsENERCBwhq2iooLq6upwT0UkLKKiooiNjW31JywFBQX4fD4AX3B5qHgkos4wAzjnXjGzh4BfAR+Z2WICKeEmAunAUuCWeg/zAfsACY30l21m04B5wKNm9nMC65VHA0OBrcD5rQmWRURk95kZcXFx4Z6GiEiNSLvoDwDn3LXAOcCHwBEEymB/D/wamOCca9VVI865F4GxBJZdDCVQLjsaeBg40Dm3zrvZi4iIiEgkibglGZFCSzJERESkI2lJRvuJyDPMIiIiIiIdRQGziIiIiEgzFDCLiIiIiDRDAbOIiIiISDMUMIuIiIiINCPi8jBHmoICXaQqIiIi7U8xR/tRWrl2YmaZBHJDi4iIiHSkgc65zeGeRFeigLmdWKCe5YDg3e+BgcDO8M1IwiQVvf7dnX4Guje9/t1bOF7/VCBbFYq9pSUZ7ST4g7o5WMAEYKeSiHc/gf+bAL3+3ZZ+Bro3vf7dW5hef/2ctQNd9CciIiIi0gwFzCIiIiIizVDA3P7KgN8Ft9L96PUX/Qx0b3r9uze9/l2ELvoTEREREWmGzjCLiIiIiDRDAbOIiIiISDMUMIuIiIiINEMBczsys7PN7D0z22FmRWb2hZnNNLPYcM9N2o+ZzTEzt4tbQrjnKbvPzPYxs2uCr/VXZlYZfF1/24LH/tTM3jCz7WZWYmZrzOz3ZpbSEXOXttud19/Mbm/B74URHfk8pPXMLNbMJprZn83sEzPLN7MKM9tiZgvM7ORdPF7v/wilwiXtxMweAK4FKoElQCEwAbgHmGxmk5xzJWGcorS/pcC6Jo5VdeRExHNXEnh/t4qZXQ/cBzjgfWArcBRwM3CmmY1zzm33cqLSLnbr9Q/6Avi8iWP+3exTOs4xwDvBr7cAHwBFwH7AZAJ/3x8FrqhfaU/v/8imgLkdmNnpBH6ZFgLHOOdWBvdnEAiexwF3AjeEbZLSER53zs0J9ySkXawC7gU+A1YS+IP3s+YeYGYHA7MI/LM02Tn3ZnB/ErAAmAj8HTir/aYtHmn161/LK86529tpXtL+qoGXgQedc+/XPmBm5wDPAZcTOGHydK1jev9HOAXM7ePm4PbuULAM4JzbbmZXEfjP8mozu9M5pzMKIhHGOfd47ftmVt2Ch90EGPBU6I9lsK9iM/s58C2Bs0wjnHNrPJ2weGo3X3/pApxzSwic+Grs2HwzOw74OXARtQJm9P6PeFrD7DEzywRGB+8+X/+4c+4DYBMQD5zUgVMTkTAxszggtLaxsd8L3xE4IwUwpaPmJSKe+yy4HRTaofd/16AzzN47OLjNc85lNdHmUwJvpoOBeR0yKwmH8WY2CkgFfgCWA28451TxqfvZG0gKfv1pE20+JbCe8eAmjkvXcIiZ3Q30JLBm+TPgNefczvBOSzyyV3CbU2uf3v9dgAJm7w0Jbjc202ZTvbbSNV3UyL4cM7vEOfdWh89Gwin0Xs9vJjDS74XuYXLwVpvfzH7lnHu6sQdIZDCzfsC04N2Xax3S+78L0JIM76UGt0XNtCkMbtPaeS4SHl8QuOhzfwKvcV9gErAM6A8sMLNjwzY7CQf9XpD1BK5vOZjA2eWeBC4Afx3wAXPN7ILwTU/awsxigGcJvJZfAbNrHdb7vwvQGWYRjznn7q+3ayfwjpm9C/wTOA14ADioo+cmIuHhnHumkd1LCaQhewi4BrjfzF50zpV37OzEA38nkOniB+AsvYZdj84wey/0cUtyM21CCcoL2nku0okEc3LeFrx7oJkNaq69dCn6vSDNuZ1AurHewNjwTkVay8weJJAZYwdwnHPuv/Wa6P3fBShg9t6G4La5YCh0bEMzbaRr+rrW1wPDNgvpaBuC23QzS22ijX4vdFPOuTwgN3hXvxciiJnNAn4F5AOTnHOfNdJsQ3Cr938E+//t3XnwlVUdx/H3BzS1cZvIFpgQohqzMgTcUUy0UZpQM5dJJJhsk5ECERQdw1HBDR3c0hkLtJi0UYFccisw3BBEsnGUDEOUQMU1ZHX89sc5V66/7n1+2733x/J5zTDnPud5zvOc3/1x7+97z/2e8zhgrr3Si6WLpGrJ+/1yubDKftt6dSl77Fnx247FwJr8uF+VY/y+sI2S1JmU+wp+X9hiSLocGE1a7eTbEVFtBQy//rcCDphrLCJeBebnzR803S+pP+mT5HrgvgZ2zTYPp+TyPdKbqG0Dcj7jvXmz0vvCnsDBeXNGo/plm43BpGXHgurLjtlmJC8NeDYpWD4qIuZXO9av/62DA+b6mJjLcyT1KVVK6gLckDev813+tj6SeksanGdMl9d3yndzKv3fuCYiNja+h9aBLiUFRMMlHV2qzLfG/Q3QGbjTd/na+kjqLmmIpB0r7DsOKN05cHpErGxs76y1JF0MjCOlYRQGy2X8+t/CKc1DslrLkwBGAhuBv5CWkxkI7E6aGX1URKztuB5aPeQ/fjNIkz8WAq+RfudfB7rnw/4ADI2IDzqkk9Zu+YPwDWVVvYBPA68Cy8vqj4+IFWXtRgFXkf5wPkLKWz2UtNzgYqB/RKyqb++tvVr7+5fUm5SutzqXy4GdgL3ZdKOL2cDgiFiNbbYkDQZm5c0FwHNVDl0VEWOatPXrfwvmgLmOJJ0EjCAtH7Y9aR3O3wNXe8mZrVPOWx9JykfrScpZFilwfgqYGhFOxdnC5XW0Z7fg0J4RsbRJ2yOBs4D9SbPmlwF3AJN8t7ctQ2t///nbxbOB/YAvkd4XPgGsAp4m3S759oj4sD49tlqRNAyY2oJDX46IHhXa+/W/hXLAbGZmZmZWwDnMZmZmZmYFHDCbmZmZmRVwwGxmZmZmVsABs5mZmZlZAQfMZmZmZmYFHDCbmZmZmRVwwGxmZmZmVsABs5mZmZlZAQfMZlZzkpZKivzv+wXHPZyPGdbA7rWapDm5n4d3dF/qTdJ3Jc2V9F7Z7/DwFrQr/52X/q2XtEzS7ZIObab9FyRdLOlJSW9I2ijpHUkLJU2RtF/Nfkgzs1barqM7YGZbvUskzYyIDzq6I1ZMUm/gTtJgyl+BFUAAK1txmseAf+XHu5NuE38ScKKkMRFxVYXrjgUuIt0uejUwD3gd2AX4Bul28yMlXRERY9vwo5mZtYsDZjOrpzXAV4DTgRs7uC/WvOOA7YGJEXFeG89xc0RMK21I2hG4CRgKXC7pnoj4Z9n+S4FxwEZgDHBdRKwvP6GkA4FLSP+XzMwazikZZlZPU3J5gaRPdmhPrCW65/LFWp0wItYBI4D3gc7A90r7JA0kBcsAJ0fE5KbBcj7Hk8CRwORa9cvMrDUcMJtZPd0HPAJ8HhjV0kaSphXlNksalvdPq1YvaTdJV+Xc2nWSXpQ0TlKnfGw3STdJeiXn2i6WdGYL+jZA0oOS3pK0RtJTkk5rps1ASXdJWiFpg6TXJc2QdFCV40NS5MfDJT0h6d1c36O5PuZ220n6maTHc9vSc3CNpG5Njp2Qrzc8V00ty0Oe05LrFYmI1cDivFne//Nz+aeImNHMOSIi5jbpd9+cH/1qfl7fk/SSpDslHdvefpuZlThgNrN6K40gjpXUpUHX3B14AjgVWEAK2rsBlwJTJPXK9ccAj5PybnsB10gaV/GMyfGk3N5uwAPAfKAvcKukiqOfkq4EHgaOBZYBM4GX8vZcScMrtcttrwVuBj4A7iXl9kZzP7ykHYA/A78G9s0/30xgB+BMYJGkPmVNFgG3AEvy9mN5+xbg/uau10K75nJ97uPuwGG57pbWniyPTj9Byo9eBcwiPc9vAN9hU/BvZtZuzmE2s7qKiHmS7iJ9FX8eMLoBlz0WuBvoFxFrAHKAOA84AziCFECeWZqMmEckZwLjJV1batfESGB8REwqVUgaQApOR0t6MCIeKNv3Y+As0iS4EyLi2bJ9hwH3ADdKejQiKqVBDAUOySkJrXEhKYVhCXBkRCzN19yeFET/CLhD0l4RsSEiZgIz84h9L5rkIbeXpH2AL+bNRbnsw6ZBm/ltOO15pHzrIRExvcn1dgO+2oZzmplV5BFmM2uE8aRR0jMk7dmA660GTi8PeiNiISlFpBOwMzCqfOWOiJgF/IM0EtqvynmfKQ+Wc7tHgBvy5lml+pz6MSFvnlIeLOd2f2PTyhA/rXK9K1sbLOdJdiPy5qhSsJyvuZEU9L8G9ASqLvlXCzktZhBwF+l5/w/wx7x7j7JDX2/D6T+by/ua7oiId9vwIcPMrCoHzGZWdxGxGPgtKSXgogZc8umIqBSElUZxZ+fJaNX2d61y3lur1JdSCvpL6pwf75vPsyQinq7Sbk4uD66y/44q9UX6kT4QvBURdzfdmT9E3JY3v9WG8zfno/xn4B1SKkkv0mj3oIh4v0bXeSqX0yX1l+RvTM2sbvwGY2aNMgEYApwq6cqmI641tqxK/epm9v83lztW2f/vZup3ArqQRkxLKQi9ShP4CuxRpX5pM+0qKU3oq9ZX2JSr3K3gmLYqX4d5A+m5eBK4v8la3G+UPf4M8Eorr3MusA8pD/0YYK2khaQPIdMj4vnWd93MrDIHzGbWEBGxQtIUUqAziTQxq62a+3bsw3bubw/lstTHlaQJgkVWVaqMiLW16lQDtTT/+RnS76ETsB+tDJgjYqWkfsAAUr72IcABuRwv6dyIuKw15zQzq8YBs5k10mXAT4BBedJbNRtyuUuV/Y3Ig66kZ5X6HrlcB7yZH5cCwDcjYlgd+9TU8lxW6ytsGv1eXnBMXUXE25LmkgLeH5LynFt7jiCNKM+Bj/K3hwHXAxMl3RERS6q1NzNrKecwm1nDRMS7wMS8eXnBoaVA7v9WOpAk0lfwHWFIlfqhuXy0LO1gPmnkeG9JX6t7zzZZQEo9+ZSkwU13StoJOCVvzm5gvyq5JJeDJR1fdKCS/kXHRMS6iLgReJb0922f2nTTzLZ1DpjNrNGuJ+UQHwBUvHEHaT1dgNMk7V2qzMuiXUb6Cr8j9JU0trwiB3GlVSmuLtXnFSkuJKVozKgU7EnqLOkIpVs/10SezHh93pxcvipJfv6mAJ8j5Ti3ZVJhzUTEQ2y6e99tkkbnNaQ/RlJfUlrLmLK6MZK6Vzh2L+DLefPl2vfazLZFTskws4aKiPWSLgCmARVvlx0Rj0maRVpPeYGkR4G1pLV7dyUFfb9oTI8/5hpgkqShpFHMrsChpMGHKRHxsSXOIuK6HNSdTbpJyXOkCXFrSUFrb9JNVn5OmhhXK78irZYxEHhe0mzShMaDSLe/fhM4MSI2VD9FY0TEGElvkSaFTgYmSJpHmiy4M2mUuEc+vDwn+XzgCkkvAM+TntOuQH/S37Zb81KCZmbt5hFmM+sIvyOteVzkZOBiYAVwOHAgMJcUNC+q3qyuZgBHkSbyDQL2BxYCwyLil5UaRMRY0kS06aQA8GjShMeupNzb04Hba9nJiFifr3MG8HdSUH88sBG4FvhmwVJ3DRcRE0mjwpOAF0hL8p1Eym9+m/QBqU9EnFPWbAQwlbS+9wDgBFLe9kOkn3VYg7pvZtsApTkTZmZmZmZWiUeYzczMzMwKOGA2MzMzMyvggNnMzMzMrIADZjMzMzOzAg6YzczMzMwKOGA2MzMzMyvggNnMzMzMrIADZjMzMzOzAg6YzczMzMwKOGA2MzMzMyvggNnMzMzMrIADZjMzMzOzAg6YzczMzMwK/A8k5eC3ZblWQAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 700x700 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot results\n",
    "plt.figure(figsize =(7,7), dpi=100)\n",
    "plt.title(\"Number of PCs vs. Variance Accounted For {}x{}\"\n",
    "         .format(IMAGE_DIM, IMAGE_DIM))\n",
    "plt.xlabel(\"Number of PCs\")\n",
    "plt.ylabel(\"Variance Accounted For\")\n",
    "plt.scatter(np.arange(1,max_pcs+1), \n",
    "            R2X_range, \n",
    "            label=\"X-variance\",\n",
    "            alpha=0.5)\n",
    "plt.scatter(np.arange(1,max_pcs+1), \n",
    "            R2Y_range, \n",
    "            label=\"Y-variance\",\n",
    "            alpha=0.5)\n",
    "plt.scatter(component_range, \n",
    "            q2y_range, \n",
    "            label=\"Q2Y\",\n",
    "            alpha=0.5)\n",
    "plt.legend()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Verifying X_scores shape: (100, 2)\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAd4AAAHICAYAAADgEYtFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzs3Xl8VPW9//HXB4gxLMoOikoQWhe0WhvqUnFBr6iouFTtFVBal9b2KiLXVnupImq1F+WWaq3YVlFww1ZRat2qFKhKa1D8VVwqEBCpCJi4ICEm8P39cc4Jk8lsmcycOZO8n49HHiFnvnPme5KQ93yX8/2acw4REREJR4dCV0BERKQ9UfCKiIiESMErIiISIgWviIhIiBS8IiIiIVLwioiIhEjBK9JOmNkxZubMbEqh6yLSnil4pSiYWVczu9bMlpnZF2ZWa2bvm9mLZjbFzPoXuo6FYGZ/9cM0+NhuZh+b2QtmdlYeXm+8/zrjs3xeso/eceU7m9kkM3vQzP7lX1ezcjHlB5jZFWb2nP978aWZrTezP5rZoa245KzEvMlJ9lGR4DmXmdksM/unmTUkK+eX7WVml5jZk2a2yszqzGyTmT1tZiPzf4XSGp0KXQGRdMxsF+Al4ADgX8Bs4GOgH3AEcB2wBHimUHWMgGnAFqAE+CpwOjDCzH7inPvfgtasqSeAZQmOb4n7ui9wq//vKuBToHuK814G/ARYCTwHbAS+gvd9ON3MznPOPdKKemdrIfDXBMf/neDYr2Ie24T3+53M2cBv/LIvAOuAPYCzgBPN7MfOuWlZ1lnyTMErxeAKvNC9G/iBi1tuzcyGAp8XomIR8r/OuU3BF2Z2HPA8MMXM7nDOxQdbocxzzs3KoNwm4D+Apc65GjP7K3B0ivL/AI5xzi2MPWhmw/GC6TdmNs85V5ddtbP2V+fclAzLnoJ3vevNbBZwQYqy/wJOA55yzm0PDprZjcDfgZvM7AHnXKKAlwJTV7MUg8P8z7+JD10A59xy59z78cfN7BAze8TMPvS74j4ws8fN7Mi4cn3N7A4zW+N3UX5oZveZ2aAE51ztf/Qwszv9c24zs9Njygw2s3v9x770P99pZn0SnO94v3s0qON6M1tgZt/J6ju143vyAvAOUAYMTVfezM4ws8Vm9rnflV9pZhfFlZkF3Ot/eW9s12lr6pqk/pudc39xztVkWP6x+ND1jy8GFgA9gAPTncfMfuZfU7PWopnd7j92VSZ1ainn3FPOufUZln3ROTc/NnT94+8Cj+D1fByRh2pKDqjFK8Wg2v/8VRJ3UzZjZufidUlvx+veXAX0B4YD3wb+5pfri9dCKMdrIT4I7AOMA0aZ2ZHOuXfiTl8KvAh0BubF1tHMDsfr8t4ZeBKvm3Rf4AfASDMb5pwLyp7il1nvf/7Yr2OFX8eHM7nWDKQMRjP7MfALvO7Z+4Av8bosf2tmX3fO/cgvOg+vu3c0ybuM0znEH6ftAKwAnnfO5bu3ot7/3JBB2Z8DJwCTzOwZ/w0MZjYK+C+8n/utKZ4f76tmNhHvd2YN3vVuSvOc1mrJ9UohOOf0oY9If+D9oXfAZuD/gJFAjxTl+wNf4I0LDo17zIDdY76e5Z/72rhy3/OPL4g7vto//mdg57jHdsL741oN7Bf32Nn+8+6IOfYYUAf0TXANvTL83vzVP2/vuOPHAtv870Nn/9gxftkpMeWG4P2BXgf0jzneDXjTL390zPHx/rHxLfwZBs+L//gEGJvtdWbwvL2ArXhjoR0zfE65X68PgJ54Y60f+T/XPTI8xzFJrncLcFUGzw9+LytaeL274L2Rq830d0gf4X+oq1kizzn3BN7EGcMb730GqDazd8zsVjPbI+4pF+C1Rqc555bHncs5f9zLzEqBc4ENQPwEpHvxgucYM9srQbV+7JzbGnfsFLw/9Lc4596Oe91HgaX+68WqZ0cLJbb8xwleM5Ufmze7+0Yzm4s3wagD3huKVOO7Y4COwK0uppvTea3QKf6XqcYaM7UK+CFe0HfGC7cf4fVI3G9mJ+bgNZowsxK8Xo9S4CfOuW2ZPM85t9qv2wDgd3i/C32BS5xzH2T48huBq4D9gC7+ucbi9Wr8r5ldkvmVtMhdeG8Ufp7F75CERF3NUhScc/9rZjOBk/HGrobhdclOAi4ys5Occ6/4xYf5n59Lc9p98LqEn48PUeecM7OFeJO6DgJix5C3OufeTHC+4LaVoZb4XtkyoLeZ9XZed+MjwBnAm2b2EF435kvOuU/T1DuRYNwxaEUuBH7tnHs8zfMO8j83Gx9lx2zcgxI81iLOuUXAophDa4A7zWwF8CzezPSczUo3sw54rcajgLudc7Nb8nzn3ANmdjJwnn9olnPuDy14/nIg9k3fFuABM3sD7w3YFDP7rXMuZ+PjZnYz8J94vTE/z9V5JfcUvFI0/EB6yP8IxmfvxBuPnAl8zS+6q/853YzOXfzPHyV5fH1cucCGJOV7+p/PT/O6XYBNzrlHzKweuBKYiPcmosHMngEmOudWpDlPrD4uu7HDpN8D59wmM2ug+fXnjHPuOTNbC3zTzDo551o9LumH7j14oXkfcGmWp5rHjuD9dWvrBeCce9PM/o4312AQXk9Aq5nZDcDVePMUzsq0dS+Foa5mKVrOuQ14k6DqgAPNrJf/0Cf+593TnOIz/3Oy+yX7xZVrfOk05zvJOWcpPtbEXMNjzrkj8UL7FLzJXacAT5vZTmnqnwtJvwf+97MTza8/1zbh/S3aubUn8kP3Xrzu8QeA77m4mb8ZnqcvcAfe79I24C6/6zoXgjdInXNxMj90J+P1mIxOMAQiEaPglWJXx44xUvM/v+p/PiHNc9/Fm3jzTX+8N95R/uc3MqzLP/zPh6UslYBz7lPn3U5yAd6M4SF444P5FsxMPirBY0fHlQEvhMAbF241M+uGN+v7Y+fc5laeKwjd8/FmhF+QTej67sEb1/0eXrftN4AbWlM/v46dgEPwxrbX5uB8QeguBE51ztW29pySfwpeiTzzlsY7JMnD/wV0Bd6N6Wq9H3/2qL+4Ruy5zMx2A3DeYgqP4LX2JsWVuwCv6/qvLsE9wknMw/tjepWZNQtfMyuzmOULzWy4mXWMK9MBCJZFDKPl8iBemP6338oL6tGVHZOr7o8pH9zatWdLXsTMvp7g2M54i6KU0cpbp2K6l88HHsWbKZ1Vd6uZ/RAYBfzeHyOfivem6iozS7WIR+w5El1vR+AWYCDwdJZj+bHnm4oXuouBUWkm0UmEWA7H9kXywszm4d1S9C7wMvAh3jjuYXgtka14f3hejHnOucAcvFCZhzeW1hevZfdn59wVfrnY+3ifBV7Dm3R1Bl7INLmP18xWAzjnypPU9VDgabz7XZ8D3sLrri3Ha0G+4pw70S+7DO/Wp5fwblMyvNuADgb+5Jw7NYPvzV/986Yd4zWzY/AWk7jexaymFHMf70d4oVXvX385cKfbcR8vZtYT7zabL/HGTzf6348b07x2jf+8ZXi3LvUBjsebBf4G3qpTn8Q951Z2vAk5Ee8N0sN4vRzgzR5/xy87BW+C1mZgBonvYZ3lz1hOVc/98CY/rQO+HrTCzWyIX/dq4GvxdU1wnmV4t2T9A+/N2K54v3v74k3UOyp2yMF/ztX+4wBHAoOBP+HNhAb4nXMuuP98PF7rvsG/3kS9BfOcc9ncay35Vuj7mfShj3QfeEH4Y7yJI1V4QbsFL4hnAvsked4w4I944VCH9wfwj8C34soF43nv4wXKerxW3qAE51wNrE5T37388630X7cG+CdwOzAspty5wFy/3Ba8P+qVeOsO75Th9+avZHh/Kwnu44157Ey8RUU2+3VZClyc5Dwn4wXKFv98LoPXnoo3q3m9/z3+3L/Wq4GyJM9ZDQnvhQ0+jokpOytN2Sblk7zeTnhvvOqBQxM8fqF/noczuN7Lgb/gvdnYinc/9T/xuq0T3oMe87NM9jE+puyUDK53fLp66qMwH2rxioiIhEhjvCIiIiFS8IqIiIRIwSsiIhIiBa+IiEiItGRkFnr37u3Ky8sLXQ0REYmQpUuXbnLONdt3O56CNwvl5eVUVlYWuhoiIhIhZrYmfSl1NYuIiIRKwSsiIhIiBa+IiEiIFLwiIiIhUvCKiIiESMErIiISIgWviIhIiBS8IiIiIVLwioiIhEjBKyIiEiIFr4iISIgUvCKSXnVVoWsg0mYoeEUktcW3wa8O9j6LSKspeEUkucW3waJp3r8XTVP4iuSAgldEEgtCt77W+7q+VuErkgMKXhFpLj50A2GFr8aUpQ1T8IpIU9VV8MLU5qEbqK/1Hs9XOGpMWdo4Ba+INNVzEBx3LZSUJX68pMx7vOeg3L+2xpSlHVDwikhzwyfBUVc1D9+SMu/48Em5f02NKUs7oeAVkcTiwzeT0M22+7nQY8oiIVLwikhyQfhC+tDNdmy20GPKIiFT8IpIasMnweXL0odutmOzhRxTFikABa+IpJcq9HIxNluIMWWRAlHwikj2cjk2m82YskgRUvCKSHbyMTbbkjFlkSKl4BWR7ORrbDaTMWWRIqbgFZHs5WtsVhOppA1T8IpI62hsVqRFFLwi0noamxXJWKdCV0BE2ojhk2DomeomFklDLV4RyR2FrkhaCl6RMGnZQ5F2T8ErEhbtMysiFEnwmtkeZna7mb1iZlvMzJlZeYJyO5vZNDP70Mxq/fJHJSjXwcyuMbPVZrbVzN4ws7PCuBZpp/K1z6xa0CJFpyiCFxgCnAPUAItTlPs9cDFwLXAK8CHwrJkdHFfuBmAKcAdwErAEeNTMTs5ttUXI3z6zakGLFKVimdW8yDnXD8DMLgJOiC9gZgcB5wHfc87d6x9bCCwHpgKn+cf6Av8N3OKcu9V/+gIzGwLcAvw5z9ci7Um6tYwhu1tv4lvQ2Z6n2FVXaUKXFJ2iaPE657ZnUOw0oB54JOZ5DcDDwEgzK/UPjwR2AubEPX8OcKCZ6X+x5Ea+9pnNVwu62KjFL0WqKII3Q0OBKufclrjjy/GCdkhMuTpgRYJyAPvnrYbSvuRjLeNc7gZUzPI1Zi4SgrYUvD3xxoDjVcc8Hnz+xDnn0pRrwswuMbNKM6vcuHFjqysr7UQu1zLOVwu62KjFL0WuLQVvXjnn7nbOVTjnKvr06VPo6kgxydVaxvnaDaiYqMUvbUBbCt4aoEeC40ELtjqmXHczszTlRHInV2sZ52s3oGKgFr+0EW0peJcDg8ysc9zx/YEv2TGmuxwoBQYnKAfwVt5qKO1brvaZba+7AanFL21EWwre+UAJcHZwwMw6AecCzznn6vzDz+DNfh4T9/yxwJvOOb1dlvzJVShEeTegfLY4k7X4AfY6PFrfB5EkiuU+Xszs2/4/v+F/PsnMNgIbnXMLnXOvm9kjwC/NrASoAi4FBhETss65DWY2HbjGzD4HXsML5xH49/qKFIUo7ga0+Davu/eIy+CEG/PzGsMnweq/wcoXmx5//xXv9RW+EnFFE7zAo3Ff3+l/Xggc4//7u8BNwI1Ad+AN4ETn3Gtxz/0fYDMwAegPvAuc45z7U+6rLZJHUQvd4Bafl2+Hj5bDuMfz8zrvv9L8eGsXJREJiTW/q0bSqaiocJWVlYWuhkh0JJttPHhEbsO3uspbNCOdy5dF602JtAtmttQ5V5GuXFsa4xWRQkgWuuB1B88+I3evpQlW0gYoeEUke+lu8QEvfJ+dnLvXa8+3VEmboOAVkez1HORNpErnldtbP9s5dm3m9npLlbQJxTS5SkSi6IQbvYlU8bOMA0Eotqb7N9VuTC9MVehKUVHwikjrjXvcG8uND99ctESTrc0M0bylSiQNdTWLSG6Me9ybxRzIR+gGYtdmVuhKkVHwikjujHscDvfHfFsbulqbWdoodTWLSG6NvBGGXZiblujgEfkdOxYpAAWviORea8MwWHqyQ0nyMlqbWYqUuppFJFoW3wYLf+H9e3t98nIrX1Q3sxQltXhFJDqCyVQNdanLqZtZiphavCISDamWnoylxTKkyKnFKyKFF8xgTqfTzgpdKXpq8YpI4aXb/KCD30Y4+scKXSl6avGKSDQEgRrf3Rx0LWuFKmkj1OIVkcILZien2vxAoStthIJXRAordtch2BG+kH48V7cTSRFS8IpIy+Qy7J6d3HTXodjwvXxZ6tCND2yRIqHgFZHM5TLsZp/h7dMbv+tQcO5UXcvx2wRmWx+1mKUAFLwikplchR0k3kIQmodvqnokC+xMqcUsBaLgFZH0chV2AM9NTr7xQXDuZLsOZbJNYCZy+SZCpIUUvCKSWq7CLjjXq79LXaakzLunN76rOVfbBObyTYRIFhS8IpJcLvfETXeuQMVFiSdVpVtkI1lgx8rlmwiRLCl4RSS5dGEH3p65mdxjm+m5Rt6Y/PH4+3wDydZvjn1DkMs3ESKtoOAVkdSShV3g/VcybymmOtfgETDu8ZafI1noxk+eykWLWSQHFLwiUReFFtjwSd7G84m0tJs2UXAecVlmoRt/DtixnGSsZJOnWtpiFskDBa9IlEXllpfqquxnIicSH5wnpOheTnWOy5d5/479HqWbPJVpi1kkT7RJgkhUxbfaoHDhEHTTJtsvN5uN6YdPav3GB8sfa/o9Wv03r+s72eSp4HWD7+MLUxW6EjpzzhW6DkWnoqLCVVZWFroa0pYlmn0bhZZZlOqVbIZyOpcv2xH21VUa05WcMbOlzrmKdOXU4hWJmnS3vEDhwjd+675iCt1ErXKFrhSAxnhFoiRqt7wkep2W7B6Urzplcj9wrCj0Foj4FLwiURKVW16qq1JP7Mpk96B8Cb5HHZJ02JWUebcmFWLyVBRmoEvkKXhFoqbQt7wEgbvgZu/rZLcKFbyb1pof6tDJ+x6Nezz8VnlUZqBL5GmMVySKCjWWuvg2WPBz79/b673PhRpbTjbxKRjfDerXREwY52LWdKaiNANdIk8tXpGoCnssNQjd7Q3NHwt7LeNkrcd047vb65uOgYcZutp0QTKk4BWJsrDGUhffBgt/kTh0A2FN7Eq1ZV9UxsAD2nRBsqDgFYm6MCZSvTAVGupSl+u0c/5DLZPWY6HHwANRm4EuRUPBK9LeZbJrUIcSOPrH+Q21lrQeo7DsY9Ra31I0FLwiknrXoA4lcOw1+Q21bFqPYYyBp2utBnXoVNr0eC7eCKil3GYpeEXEEx++wX2y+Q5dyL71mM8x8JbcHtRQt+P7lYvQ1a1JbZpuJxKRHWI3Dzj2p+HdjhP72i1dCzq2frlaeznT24NiywW3MuUidHVrUpumFq+INBXbigx7fLI1Y7e5aiVmentQfLnt9d4EtDBeW4qagldEmivkhKBsxm5T3YLUEplO8EpWrmFr9q+vW5PaDW0LmAVtCygSgky7jXO1VWF1lddiTueCP8F9p6QvF7v9YK5euyXnlNBlui2gWrwiEk3Zhi5k10rMdILXoOG5v41Itya1KwpeESlO+VjAItPFOfKxiEdUFgaRvFPwikhx6jnI2/4vmWxbiZlO8MrHIh5RWBhE8k7BKyLFafFt8P4ryR/f6/DsAyvTCV75WMQj7M0xJHSaXJUFTa4SKbCwJiNlOsErV/cP5/uckleaXCUibVdYk5EyfX4+AlKh22YpeEWkOGkykhQpBa+IFC9NRpIipOAVkeIWPxlp6JmtO592BZI8U/CKSPEL1peG1q3XrF2BJAQK3iK0prau0FUQiYbY1unyx1q3XnOu1nsWSUPBW2RmrFnPoUveZsaa9YWuikhhxbZOW7urj3YFkhBpP94iMmPNeqav/giAX/qfJwzsX8gqiRRGbOt0wc2Ag+0NTcsE4QleV3Sy+2LTrfccPF8kR9TiLRIz1qxnWtV66rZ7C57Ubnf8cvVHavlK+5NoH9z40A0E6zU/Nznx2G0+1nsWSUPBWwSC0G2IW2RM4SvtTrLWaTIlZd56zq/+zvs6vvtYuwJJASh4I25NbR03r2oeuoHa7Y6bV63XhCtp+9K1TuOVlHnrNb//SuqxWy3EISFT8EbcvA01KQfiOxlcs3d/BpaVhlYnkYJI1zrtUAId/P8tiUI3kEn4KnQljxS8EdbY2k1RpsHB6X17hFYnkYJK1To99ho49qfe1xUXwcoX04/dVi1ufm5Q6EpeaVZzhA0sK+Wavfvzy9UfUbu9eV9zJ+AqtXalvQkCMRjrjW+dDj3Tax137pF8PDhoEd93iteKjt3gPni+SJ6oxRtxEwb254ryfpR1sCbHg9DV7UTSLqVqnQahmap1HHRDQ+IJVyJ5pBZvEQjCNWj5lhpcOUihK+1cJq3TRK3jZBOuYsuL5JFavEUiaPmCQlekUSat09jWcUsmXInkiTmX5D4VSaqiosJVVlYW5LXX1NZpTFckG1WLvTHddC5fpu5myYqZLXXOVaQrpxZvkVHoimRp0HAtliGRoOAVkfZDi2VIBCh426GorXIVtfpIG6fFMqTA2m3wmtmeZvYHM/vUzD4zs8fMbK9C1Ses8InatoLJ6qMwlrzSYhlSQO0yeM2sM/AisC9wATAO+AqwwMy6hF2fXIdhstCasWZ943aCUdhcIVl9ovbmQNqo4ZO8iVQKXQlZuwxe4GJgb+B059w859wTwGnAQOD7YVYk12GYLLSC16mNyLaCyepz7rIVTK/y6pRJ/dQyllbRRCopgPYavKcBS5xzK4IDzrkq4CVgdFiVyHUYpmpBJlp2slDhm6o+C2s2U+cyq59axiJSjNrrylVDgScSHF8OnB1GBdKFIdCiRTKShXhNfQN3rd2U9HnBtoKn9+0Ryq1KwcYPmUr2/Yh/kxH/uIhIVLXXFm9PoCbB8Wog4VY/ZnaJmVWaWeXGjRtb9eJB+CTa+ABavsduqhC/b93HHN2ja7O1ngNlHSxn2wpmUt9g44dk9Ukk/vsRtW5zEZGWaK/B22LOubudcxXOuYo+ffq06lzpwqclYZhJiC+s2cz4Ab2bvV5ZB+OK8n45aSm2pNs32cYPyXRix57DUes2FxFpqfYavDUkbtkmawnnXLLwaWkYZhri1w0Z0OT1ch26LZ0gFlx/aQbZ24C353CuewpERAqhvQbvcrxx3nj7A2+FVYn48M02DDMN8diNFnIdutl0+04Y2J/v7tE7o9e5d93GnPYUiIgUSnudXPUkcKuZ7e2cWwVgZuXAt4Crw6xIEH43r1rfqjCM3zowWYhPGNg/ZxOpcjFBbMqQPTCMmWs3sj1Fud+t3USPkk7NrjOQ6ZsWbTIhIoXWXlu8vwVWA0+Y2WgzOw1vlvNaYGbYlZkwsD9/P2y/VrdAM23R5moiVUu7fZN1AY8f0Dtl6ILX3Ty9aj0z1qzPuqdAtx9JQVVXFboGEhHtMnidc18AI4B/AbOBB4AqYIRzbnMh6pSrVliuQjydlnb7pgq94Fzpul/q3I4x5JZ2m0dt1S5pZxbfBr86WPv9CtB+u5pxzr0PnFXoeuRDWF2pp/f15qdNq1pPQ0zDN74Fmsk9t8HXt6xaT6odomPvO8602zzZOHSieojk3OLbYNE079/BZy1T2a61yxavtF7Qgn25ZjOxSdnJSBi6mU6+6phmlnN8S7qloRvQ7UcSiiB062u9r+trva/V8m3XFLzSYrEt2IU1m2mIfdA1L5cq9IJx32DMuCFVcxfo1aljxq1U3X4kBRUfugGFb7un4JUWSRamgQa87uQpKz7IKPSCcd95G2romMHrf/BlA9evWJdRXXX7kRRMdRW8MLV56Abqa73HNeGqXVLwFplCts7ShW6gdrvjrrWbuHTPPklDrxM0Bu20qvVMW7WebRnW4zdrN2b8fcjVQiUiLdJzEBx3LZSUJX68pMx7XLsjtUsK3iJSyNth0nXbxkq2WlagE7DN/wBocDTtrk7j6B5dW1A6dwuViLTI8Elw1FXNw7ekzDuuCVbtloK3SBT6dphMNzdItlpW8LwgdNPHd3Iv1Wxu8RuQfKzaJZJWfPgqdIV2fDtRMYnK7TDJVo0KpFotC7zVuVrSsk0mOMf0Fn4Pcrlql0jGgpB9YapCVwAFb+Tlet/e1kq0NOU3d+3CwprNzUI3dnnGCQP780n9Nn77wca0M5czVZfF90ChKwUxfBIMPVNjugKoqznSono7THy37SMHD2lcLSt2z9zY7uA1tXX8Zm3uQjeg+3GlaCh0xafgjbAo3w4TvzRlsFfuoUve5txlKxKOR2cyRhyvR8f0v6K6H1dEiomCN+KifDtMbODHL6oROx49rcoLZCDhtaRSsy3d9gnN34AogEUkyhS8RSDqt8OkXVTDPzxtldfyvaK8X04nF3xz1y5NlqjUDkQiEmUK3iIR1dthMl1UA3Zs7fdyzWbIYE3mo3t0pTSDxvHCms2sqa1r1S1XaiWLSFgUvEUkrC3/MtWSRTUCdc5f3znFU4IW/SMHD+HKQam3Cwy6medtqMloM4ZEAatWsrR5WpoyUhS8RSZKt8PM21ADeDsS5dIFA3o1WYDjW0lWqgoCGhLfWxwfvokCttALk4jknfYCjhzdx9sOxN5PmyuxgYXzwre1twqVdTDGD+jNdUMGNHmdf3z6RePX5r1cY+gO26ULZy5bmfScwYznT+q3MWvdJmDHnsDBvwu9MIlI3mgv4EhSi7eNy0c3avy4bgM0hi/sGJ+F5u/sgseSzdKOD934lmywX2/Q0j1z2cqE54t/vVnrNjWbaT0tQTe57guWNkN7AUeWgrcNy0c3arLJVEH4QtNFNa6KuXc3duw23SztpK/joLSD8XLNZqZVedfz0ieb+eauXRKG+Td37cI/Pv0i4XmSLV+p+4Kl6Gkv4EhT8LZRydZ3bk34pptMFQTZ6X17AN54dLLZ2Klmaad7nbrtrskErQbXPHzLOhgXDOjV5J7iTGmfXilq2gs48hS8bVC69Z2zDd9sV9JKNhs72fFMd0KKFRu+4IX5lCF7pDxPJ0vcFR6l27VEWkx7AUeegreNyff6ztmupJWs9ZjseLLXSaXBv1XpsYMHJ92WMLa+Vw3qn7ArXKErRU97AUeagreNCWN951yupJXqDUD865RmGMIDdt4p4/pGdWESkVbTXsCRpeBtg8JY3zkXgZXJjOvY17myvB/X7J16QY2je3RN+KYiVX2jtjD+Y+yFAAAgAElEQVSJSM4E4QsK3Qgx53K8T1s7UFFR4SorKwtdjbRix3rz1Y2azT3Ca2rrmqw0FVu3ZOeLPX7ushUsrNmc8NzprjMf9zSLRF51lcZ0Q2BmS51zFWnLKXhbrliCF7zwvXnVeq7Zu38kWnRBfeIX3Ahu/VlYszllXdfU1jXudJTK3w/bTwErIqHKNHjV1dzGRakbdcaa9Y07FMWvclXr3yIEqe85bs0Ytu7LFZEoUPC2A1Fo+c1Ys57pVeuTLloRK91tT9mMYWsjBBGJCgWv5F0w1lzXglGNloZvutDVRggSedksaKFFMIqSglfyKputAwPp7jnOZGZ1PlbwEsm5bHYQKtSuQwr7VlPwSl5lswpVIJN7jlONYedrBS+RnIrfQSiTIM3mObmgLQZzQsEreZdsTDZdFH9z1y4ZTQpLNpEqnyt4ieRENjsIZfqcXLdMCxX2bZCCV0LRbBUqa9zMKKmFNZuzDsYwVvASaZVsdhDK9Dm5bplqi8GcSrUIkEhOBa3Xm1et58pB3r8TdQXDjslSrV3aMtFraE1mKbhgB6Fkgh2Ehp65Y+GLTJ9TvRrefNQ7FrRQW7NiVbqwb+352yG1eCVUsWOyYS5tqY0QJFKy2UEok+f03Btevz93LVNtMZgXavFK6GJbsfGt0nwEY2xL+4ryfny3R2dWrFhBfX19zl5DpMV6nwwnHwlbP6Gk7mP6Lv89u2ysTL2ZQXAsvgVaUgbddoPqVc2f05qWaRD2iVq8wesedZWWo2whLRmZhWJaMrJYhLG05ZraOnrU1/HRRx8xYMAAysrKMGv5bGuRXHKffUjtJx+x7sOP6Pf6DHYZelz6gIzt/i0pg90OhvdfSf9ily/LLiQTdTdrt6NmMl0yUi1eyZnWbEAwYWB/Tu/bI6+TnQaWlbJi3VoGDBhA586d8/Y6Ii1hu+xGZzMGuO38+5tXs8shR6V/UhB2L0yFvQ5PH7qtbZnGt7QVuq2iMV7JiVwsyRjGDOP6+nrKypKMkYkUSrf+lO2+P/VlfTJ/zvBJcMGfYOWLycdgAxUXtT4ktcVgzih4pdWKbUlGdS9LFFnJzi1/0qDhqSdcAQweASNvzL5isYZP8rqrFbqtouCVVtGSjCIFFrREE4Xv4BEw7vHkz81mNrImUrWagleypiUZRSIiPnxLyuCIy1KHrpZ/LBgFr2RFSzJGx6xZszAzzIx//etfzR5fuHBh4+N/+ctfClDD3Bo/fjx77LFHoavRIsHPaPXq1fl7kfgx2BNSdC9r+ceCUvBKVrQkY/R069aN2bNnNzt+33330a1btwLUSEKXyRisln8sOAWvZC2MlaeKSoFX7znzzDOZM2cOsffm19bW8oc//IGzzjqrgDWTUKUag81mfWjJOQWvtIqWZPRFYLxs3LhxrFmzhr/97W+Nxx5//HG2b9+eMHgXLlzIcccdR7du3ejSpQsjR47kzTffbFLm2Wef5YgjjmDXXXela9eu7LPPPkyd2nS94DfeeIPTTjuNHj16UFZWxre+9S0WL17cpEzQPfz6668zfPhwOnfuzFe+8hXuuuuuZvWqqqpi3Lhx9O/fn9LSUvbee28mTJjQrFyqcy1duhQz44knnmj2vKAu27ZtA6C8vJyxY8cye/Zs9tlnH8rKyhg+fDjvvfceX3zxBd///vfp1asX/fr1Y9KkSTQ0NDSea+vWrUycOJEDDjiArl270r9/f0499VTeeeedZq9bcFr+MTIUvNJqmWxI36ZFZLxs4MCBHHXUUU26m++//37OOOMMunbt2qTsU089xXHHHUfXrl2ZM2cODz74IJ9//jnDhw9n7dq1AKxatYrTTjuNQYMG8cgjj/Dkk09y5ZVX8sUXXzSe57XXXuOII46gurqa3/72t/zxj3+kV69eHH/88SxdurTJa3722Wecd955jB07lieeeIJhw4Zx6aWXsmDBgsYyVVVVfPOb32TRokVMnTqVZ555huuuu45Nmza16Fzf+MY3GDZsGDNnzmzyvE8++YS5c+dy0UUX0bFjx8bjixYt4s477+QXv/gF9913HytXruSss85izJgxdOvWjYcffphLLrmE6dOnc/fddzc+r66ujs8//5zJkyfz1FNP8Zvf/IatW7dy+OGHs359xCYXZrM+tOSHc04fLfz4xje+4aS51Vu2FroKab311lu5PeGiW527sZ9z1+2y4+PGft7xkNx7770OcO+99577/e9/77p37+5qa2vdv//9b9exY0f33HPPuQULFjjAPf/888455wYPHuxGjBjR5Dyffvqp69Wrl5swYYJzzrlHH33UAe7TTz9N+tojRoxw++67r6urq2s81tDQ4Pbdd183evToxmMXXHCBA9yLL77YeGzr1q2uZ8+e7uKLL248Nm7cONelSxe3bt26pK+Z6bnuvfde16FDB7d69erGYzNmzHAdO3Z0a9eubTw2cOBA16NHD/fJJ580KQe4Cy+8sMlrf/3rX3fHHHNM0ro1NDS4L774wnXt2tVNnz69SV0AV1VVlfS5zuXh9zORCPzOtlVApcsgQ9TilZxpdxOpIjhedvbZZ1NXV8f8+fN54IEH6N+/P8cdd1yTMu+99x4rV65kzJgxNDQ0NH507tyZww8/nEWLFgFw8MEHU1JSwne+8x3+8Ic/sGHDhibnqa2tZeHChZx99tl06NCh8TzOOY4//vjG8wQ6d+7Mscce2/h1aWkpX/3qV3n//fcbjz333HOccsop7L777imvM5Nzfec736F79+789re/bTw2c+ZMRo0a1WxW9OGHH86uu+7a+PW+++4LwMiRI5uU23fffRt7BAJz587l0EMPpXv37nTq1IkuXbqwefNm3n333ZTXUDCJbj3SSlShUvCKZCOi42XdunXj9NNPZ/bs2dx///2MGTOGDh2a/jcPAvTCCy+kpKSkycef/vQnPv74YwCGDBnCs88+y/bt2xvHXA877DAWLlwIQHV1Ndu2beOGG25odp477riDmpoatm/f3vi6PXr0aFbf0tJStm7d2vj1xx9/nNGtQpmca+edd+a73/0u99xzDw0NDSxevJi33nqLH/zgB2nPt9NOOyU9Hvsa8+fP59xzz2W//fbjwQcf5O9//zuvvvoqffr0aVIucrT8Y0FlvEmCmZUCvYEPnXPbEzzeDfi6c25RsyeLtDUR3i7t/PPPZ9SoUWzfvp2HHnqo2eO9evUC4Oabb+b4449v9ngQOgDHHnssxx57LHV1dbz00ktce+21jBo1itWrV9O9e3c6dOjAj370I84///yEdYkP/XR69+7NunXrWvScVC699FKmT5/OE088weOPP055eXmzVmxrPPzwwwwZMoRZs2Y1Hquvr6e6ujpnr5E3wyfB0DM1plsAaYPXvIVtbwEuA0qBT8xsGvC/cQG8P7AA6Nj8LCJtUKq9UQvYiviP//gPzjnnHLp3787QoUObPb7PPvtQXl7O8uXLufrqqzM6Z2lpKSNGjGDz5s2MHj2aqqoqhg0bxvDhw3njjTc45JBDWhyyiZxwwgk89thjfPjhh+y2226tPt/gwYM54YQTmDZtGsuWLePaa6/NST0DW7ZsoVOnpn9GZ8+e3ThjOrKqq7zAVegWRCYt3u8C/w3MAl4FjgFuAkaY2ZnOuc15q51I1EVwu7SOHTsmbOkGzIxf//rXjB49mi+//JJzzjmH3r1789FHH/Hyyy+z1157ceWVV3LXXXexaNEiTj75ZPbcc082bdrEzTffzO67784BBxwAwPTp0znqqKMYOXIkF154IbvtthubNm3itddeY9u2bdxyyy0tqvv111/Pn//8Z4444gh++tOfMmTIENatW8czzzzDnDlzsvp+/PCHP2T06NGUlJRw4YUXZnWOZE488UTmzZvHxIkTOeWUU6isrOT222+ne/fuOX2dnFp8mzcMcty16mIukEyC90fAdOecPyDAXWZ2L/Ag8KKZneicK4J+FZE8id0btUjGy04++WQWLVrETTfdxEUXXURtbW3jGO65554LwEEHHcTTTz/NNddcw4YNG+jZsydHHnkkDzzwQOPWiocccgivvvoq119/PZdffjmffvopffr04ZBDDkk4lppOeXk5S5YsYfLkyVxzzTVs3ryZAQMGMHr06KyvddSoUZSVlTFq1Cj69euX9XkSufjii1m7di333HMPM2fOZNiwYcyfP58zzjgjp6+TM/G3vkFR/L62NeZc4rV2GwuYfQaMds4tiDu+P/A8UAP8B7AX8LJzrs13NVdUVLjKyspCV0Oy8Pbbb7Pffvvl5+RB951EyvPPP88JJ5zAX/7yl2YzvKMmr7+fiWbhR6CHpi0xs6XOuYp05TIZ7NgCdI0/6Jx7CxgOdAEWA/qLI+2bQjdSVq5cyfPPP8/EiRM55JBDIh+6eRXBW9/as0yCdzkwItEDzrlVeOHbANyTw3qJiLTKDTfcwEknnURpaSn3339/oatTOBG99a09yyR4nwbGm1nC2QLOuQ+AI4EILk4qIu3VrFmzaGhoYOnSpQlnd7cbWioycjIJ3unAHsBnyQo45zYBRwB756heIiKSK/GrVQU0xlsQaWc1+/fqfpFBua3AmlxUSkREciyCt761V2lbvGb2dTP72MySzuc3s9F+mQNyWz0REckZLRUZCZncx/tfwBvOueYbW/qcc0+Y2VLgcuCSXFVORERyTEtFFlwmY7zHArPTloIHgHY8X19EpEgodAsqk+DdHViRQblVwIDWVUdERKRtyyR4a0mwgEYCXYEI74MlIiItpvt7cy7TBTQy6UI+HnizddURkZaaNWsWZtb40a1bNw466CDuuOMOGhoaCl29ggu+P6tXry50VYrP4tvgVwdrZascy2Ry1YPANDOb65z7R6ICZnYY8H3gqkSPi0j+Pfroo+yxxx589tlnPProo1x22WVs2LCBqVOnFrpqUoy0oULeZBK8dwPfARaa2d3AfHbcrzsQOBVvJvPf/bIi7dKa2joGlpUW7PUPPvhghgwZAnj72q5YsYIZM2YkDF7nHPX19U02vS8W27ZtwznXbB9cyaH4tZ2DNZ1B4ZsDabuanXMNwIl4s5YvBZ7DWx7yHeBZ4IfAHGCUcy7iuz+L5MeMNes5dMnbzFizvtBVaTRs2DA+++wzNmzYQHl5OWPHjuWee+5h3333ZaedduKpp54C4MMPP+T888+nd+/elJaW8rWvfa3Z3rdBd+2SJUsYM2YMu+yyC7vvvjuXX345W7c2ndqxZcsWfvKTnzBo0CB22mknBg0axE033cT27dubnS+++3fKlCmYWZNjZsb//M//cMsttzSe85///Cdbt25l4sSJHHDAAXTt2pX+/ftz6qmn8s47Wr22VbShQt5l9JbRObcFuMjMpgLfAvb0H1oL/NU592Ge6icSeTPWrOeXqz8CaPw8YWD/QlYJgKqqKjp27EjXrt7cyAULFrBs2TKuu+46+vbtS3l5OV988QVHH300NTU1/PznP2fPPfdkzpw5jBs3ji1btnDJJU1vyx83bhz/+Z//yWOPPcYrr7zClClT6NGjB9dffz0ADQ0NjBw5krfeeouf/exnHHjggSxZsoQbbriB6upqbrstuz/as2bNYu+99+bWW2+lS5cu7L777tTV1fH5558zefJkdtttN6qrq7nzzjs5/PDDefvtt+nfv/A/g6ITbKiQTLChgu4DbpW0wWtmHYGfAROAXYBteN3NFzrnPslv9USiLQjd2u3evta1213Bwnfbtm00NDTw+eefM3fuXB577DFOPfVUOnfuDEBNTQ1Lly5tEkh33HEH7733HgsWLOCYY44B4KSTTuKjjz5i8uTJXHjhhXTsuGOL7fPOO68xZI8//nj+/ve/89BDDzUee+ihh/jb3/7GwoULOeqoowAat+O7/vrr+clPfkLfvn1bfG3OOZ577jnKypquNfy73/2uyfWPHDmSfv368dBDDzFx4sQWv067F2yokKjFCzuWmVTotkoms5p/AFwLvA7cCjwBjAb+L4/1asLMrjSz+Wb2oZk5M5uSouzpZva6mW01szVmNtl/8xBf7kgze9nMas1svZlNN7Mk23eINBcfuoEgfMPudt53330pKSmhZ8+e/PCHP2TMmDHcc8+O3ToPO+ywZq3ARYsWMWDAgMbQDYwdO5aNGzfy1ltvNTk+atSoJl8feOCBvP/++41fP/PMMwwcOJAjjjiChoaGxo8TTjiB+vp6lixZktW1nXjiic1CF2Du3LkceuihdO/enU6dOtGlSxc2b97Mu+++m9XrCNpQIQSZdDVfDPzWOff94ICZfR+4w8y+75z7Mm+1a1qHz4B5eG8EEjKzkcAfgd8DVwJfB34OdAN+ElPua8DzeGPUpwCDgGl4C4Ccm5crkDZlTW0dN69KHqy12x03r1rP6X17hDbh6vHHH2ePPfagW7duDBw4kJ133rnJ47vttluz51RXVyc8HgR0dXV1k+M9e/Zs8nVpaSl1dXWNX2/YsIE1a9ZQUlKSsI4ff/xxZhcTJ1Ed58+fz7nnnssFF1zAddddR+/evenQoQMnn3xys3FnaSFtqJBXmQTv3sB/xx17BPgN3qzm93JdqQSGOue2m1knUgQvcAvwN+dcMDC1wMy6ApPN7P+cc8FfyuuBD4CznXP1AGb2JXCfmf3COfdanq5D2oiBZaVcs3f/hC1egLIOxhXl/UKd5XzAAQc0zmpOJH7SEnhBmqh1uH79+sbHW6JXr14MGjSIuXPnJny8vLwcoPFNwZdfNn3fniyYE9X94YcfZsiQIcyaNavxWH19fbM3C5KlIGRfmKrQzbFMupq70nwv3s/9z91yW53E/K0JUzKzPYGD8WZYx5oNlAAn+eVK8GZpzw1C1zcX+BKvG10krQkD+3NFeT/KOjQNhSB0ozDBKp2jjz6aDz74gJdeeqnJ8QcffJC+ffuy//77t+h8J554ImvXrqVr165UVFQ0++jduzcAAwcOBODNN3esudPQ0MBzzz2X8Wtt2bKl2S1Fs2fPZts23VyRM8MnweXLFLo5lumNcAPMLHaT+44xx5tMsHLOrcpJzVpuqP+5yepZzrkqM9sCBH9BBgM7Jyi31cxWxpQTSSsI16DlW0yhCzB+/HhmzJjBmWeeyU033cQee+zBAw88wPPPP8/MmTObTKzKxJgxY7j33ns57rjjmDRpEgcddBBffvklK1eu5Mknn2TevHl07tyZYcOGMXjwYK666iq2b99OaWkpd955Z5Nu63ROPPFE5s2bx8SJEznllFOorKzk9ttvp3v37i39NkgqmkiVc5kG7x+SHJ+X4FjL/qfmTtAnVpPgsZqYx1OVq455vAkzuwR/y8O99tor+1pKmxOE7M2r1hdV6AJ06dKFhQsX8uMf/5irr76azz//nH322YfZs2czduzYFp+vpKSEZ599lltuuYW7776bqqoqunTpwuDBgxk1alTjgh2dOnXiiSee4Ec/+hHjx4+nZ8+eXHHFFRx66KGNM6TTufjii1m7di333HMPM2fOZNiwYcyfP58zzjijxfUWCZM513x8qkkBswtackLn3H1pznc83sSmdBY6546Je24noB643jk3Je6x8/AW+djPOfdO3GMfAM865y40syOAl4CTnHPPxJX7G1DnnEu5NnVFRYWrrKzM4BIkat5++23222+/vJy70CtXSfHL5++n5J+ZLXXOVaQrl7bFmy5Is/AykMlv1pYWnjdowfZI8FgPvNZsunI98TaFEGkxha6IZCL0xU79VbDysaZbEJhDgVeCg2ZWDnQGghsSVwJ17BgTDsrtjDeD+9E81E1ERATIbFZzUXDOvQ+8AYyJe2gsXvf00365L4FngHP8ruvAt4FS4Mn811ZERNqrotjew8wqgHJ2vFHY38y+7f/7z34rGuCnwJ/MbCbwEN4CGpOBGTH38AJMAZYAc83s1/65pwF/cM4tzeOliIjkTnWVZh0XoWJp8f4XXhfwI/7XZ/tfPwo0LvzqnPszXsv1MLxVqSbirVx1dezJnHPLgBOA3YCn/DL3Ay2aSCYiUjDapL5oFUWL1zk3HhifYdnHgMcyKLcIOLxVFZOi5JxLuBKSSCGlu8OkCW1SX9SKpcUrkhMlJSXU1ibYdUWkwGpra5Oucd1Esk3q1fItGgpeaVf69u3LunXr2LJlS8taGCJ54pxjy5YtrFu3Lv2Widqkvk0oiq5mkVzZZZddAPj3v/9NfX19mtIi4SgpKaFfv36Nv58JaZP6NkPBK+3OLrvskvoPnEgUaZP6NkNdzSIixUKb1LcJCl4RkWISH74K3aKj4BURKTZB+IJCtwhpjFdEpBgNn6SJVEVKLV4RkWKl0C1KCl4REZEQKXhFRERCpOAVEREJkYJXREQkRApeERGRECl4RUREQqTgFRERCZGCV0REJEQKXhERkRApeEVEREKk4BUREQmRgldERCRECl4REZEQKXhFRERCpOAVEREJkYJXREQkRApeERGRECl4RUREQqTgFRERCZGCV0REJEQKXhERkRApeEVEREKk4BUREQmRgldERCRECl4REZEQKXhFRERCpOAVEREJkYJXREQkRApeERGRECl4RUREQqTgFRERCZGCV0REJEQKXhERkRApeEVEREKk4BUREQmRgldERCRECl4REZEQKXhFRERCpOAVEREJkYJXREQkRApeERGRECl4RUREQqTgFRERCZGCV0REJEQKXhERkRApeEVEREKk4BUREQmRgldERCRECl4REZEQKXhFRERCpOAVEREJkYJXREQkRApeERGRECl4RUREQqTgFRERCZGCV0REJEQKXhERkRApeEVEREKk4BUREQmRgldERCRECl4REZEQRT54zeyrZjbDzP6fmW02sw/N7EkzOyhJ+YvN7B0zqzOzd83sB0nKnW5mr5vZVjNbY2aTzaxjfq9GRETau8gHL3ACcCxwH3Aq8EOgD7DEzL4RW9DMLgZmAn8ETgQeBe40s0vjyo30y7wKnATMACYDP8/rlYiISLtnzrlC1yElM+sNfOxiKmpmuwKrgfnOufP9Y52AfwNPO+cuiCl7D3AasJtzrt4/9jrwmXPu6Jhy1+KF717OufWp6lRRUeEqKytzdIUiItIWmNlS51xFunKRb/E65za5uHcHzrlPgX8BA2IOH47XEp4Td4rZQC/gSAAz2xM4OEm5ErwWsIiItBfVVaG+XOSDNxEz6wkcALwdc3io//nNuOLL/c/7pyrnnKsCtsSUExGRtm7xbfCrg73PISnK4AVuBwz4Zcyxnv7nmriy1XGPJysXHOuZ4DhmdomZVZpZ5caNG1teYxERiZbFt8Giad6/F00LLXxDD14zO97MXAYff03y/GuA84D/cs6tCKvezrm7nXMVzrmKPn36hPWyIiKSD0Ho1td6X9fXhha+nfL+Cs29DOyXQbkt8Qf8W4N+Dkx2zt0T93DQgu0BfBhzPGjBVicoF69HTDkREWmL4kM3EIQvwPBJeXv50IPXObcFeKelzzOzccCdwG3OuZsSFAnGcofSNHiDMdu3EpR7Jeb85UDnmHIiItLWVFfBC1OTP15f6z0+9EzoOSgvVSiKMV4zOwO4F/idc+6/kxR7BdgEjIk7PhavFfsSgHPufeCNJOXqgadzVG0REYmanoPguGuhpCzx4yVl3uN5Cl0oTFdzi5jZUcBDeGE5y8wOi3m4zjn3OoBzrt7Mfoa3YMY64C/ACOB7wGXOuS9jnvdT4E9mNtM/99fx7uGdke4eXhERKXJBN3J8d3NJGRx1VV67maEIghcvPEuBQ/BbrTHWAOXBF865u8zMAZOAq4D38SZh3Rn7JOfcn83s28B1wHjgI7yx40Rd2CIi0tbEh29IoQtFsHJVFGnlKhGRNmLxbd6Y7nHXtjp0M125qhhavCIiIvkxfFJeJ1IlUhSTq0RERPImxNAFBa+IiEioFLwiIiIhUvCKiIiESMErIiISIgWviIhIiBS8IiIiIVLwioiIhEjBKyIiEiIFr4iISIgUvCIiIiFS8IqIiIRIwSsiIhIiBa+IiEiIFLwiIiIhUvCKiIiESMErIiISIgWviIhIiBS8IiIiIVLwioiIhEjBKyIiEiIFr4iISIgUvCIiIiFS8IqIiIRIwSsiIhIiBa+IiEiIFLwiIiIhUvCKiIiESMErIiISIgWviIhIiBS8IiIiIVLwioiIhEjBKyIixaO6qtA1aDUFr4iIFIfFt8GvDvY+FzEFr4iIRN/i22DRNO/fi6YVdfgqeEVEJNqC0K2v9b6ury3q8FXwiohIdMWHbqCIw1fBKyIi0VRdBS9MbR66gfpa7/Eim3Cl4BURkWjqOQiOuxZKyhI/XlLmPd5zULj1aiUFr4iIRNfwSXDUVc3Dt6TMOz58UmHq1QoKXhERibb48C3i0AUFr4iIFIMgfKGoQxegU6ErICIikpHhk2DomUU3phtPLV4RESkeRR66oOAVEREJlYJXREQkRApeERGRECl4RUREQqTgFRERCZGCV0REJEQKXhERkRApeEVEREKk4BUREQmRgldERCRECl4REZEQKXhFRERCpOAVEREJkYJXJELW1NYVugoikmcKXpGImLFmPYcueZsZa9YXuioikkcKXpEImLFmPb9c/REAv1z9kcJXpA1T8IoUWBC6tdsdALXbncJXpA1T8IoUUHzoBhS+Im2XglekQNbU1nHzqvXNQjdQu91x86r1mnAl0sYoeEUKZGBZKdfs3Z+yDpbw8bIOxjV792dgWWnINRORfFLwihTQhIH9uaK8X7PwLetgXFHejwkD+xeoZiKSLwpekQKLD1+FrkjbpuAViYAgfAGFrkgb16nQFRARz4SB/Tm9bw+N6Yq0cWrxikSIQlek7Yt88JpZNzOba2YrzOwLM/vEzP5hZmMTlO1gZteY2Woz22pmb5jZWUnOe7GZvWNmdWb2rpn9IP9XIyIi7V3kgxfYCWgAbi/I4NoAAAw+SURBVAZOA84D3gZmm9nEuLI3AFOAO4CTgCXAo2Z2cmwhM7sYmAn8ETgReBS408wuzd9liIiIgDmX+Ob9qDOzV4CuzrkD/a/7AmuBW5xz18WUewHo45z7mv91J+DfwNPOuQtiyt2DF+y7OefqU712RUWFq6yszPUliYhIETOzpc65inTliqHFm8zHeC3hwEi81vGcuHJzgAPNbJD/9eFAnwTlZgO9gCNzX1URERFP0QSveTqZWS8zuwQvaP8vpshQoA5YEffU5f7n/WPKAbyZppyIiEjOFU3wAj8C6oFNeGO4E5xz98c83hP4xDXvO6+OeTz2c02ack2Y2SVmVmlmlRs3bsym/iIiIuEHr5kdb2Yug4+/xj31EWAY3qSp3wG3m9n3w6q3c+5u51yFc66iT58+Yb2siIi0MYVYQONlYL8Mym2J/cI5txEImprPmFln4FYzu8efDFUDdDczi2v1Bi3YoEUbtHR7AB+mKCciIpJzoQevc24L8E4OTlUJXAD0Az7AG6MtBQbTdJw3GLN9y/8cjOUOpWnwxpcTERHJuWIa4413NLAZ2OB//QzeGPCYuHJjgTedc1X+16/gjRMnKlcNvJSX2oqIiFAEazX747iHAX/Ba9n2As4Bvg1c7Zz7EsA5t8HMpgPXmNnnwGvAucAIvPtz8cvVm9nP8BbMWOefdwTwPeCy4HwiIiL5EPngBf4JjAZuxRuH3YS3ctUpzrmn4sr+D14reALQH3gXOMc596fYQs65u8zMAZOAq4D3gf9yzt2ZzwsREREp2pWrCkkrV4mISLz2sHKViIhI0VGLNwtmthFYU+h6tFBvvG76tk7X2bboOtuWtn6dA51zaRd6UPC2E2ZWmUkXSLHTdbYtus62pb1cZzrqahYREQmRgldERCRECt724+5CVyAkus62RdfZtrSX60xJY7wiIiIhUotXREQkRApeERGRECl42xAzu9LM5pvZh/6exlNSlD3dzF43s61mtsbMJptZxwTljjSzl82s1szWm9l0MyvL64VkwcxWJ9nX+fQEZS82s3fMrM7M3jWzHxSizi1lZnua2R/M7FMz+8zMHjOzvQpdr2yZ2TFJfmafxJXrYWa/M7NNZvaFmf3FzA4sVL1TMbM9zOx2M3vFzLb411OeoNzOZjbN/79a65c/KkG5DmZ2jf/7vdXM3jCzs8K4llRacJ3J9ls/OK5cJK8zXxS8bcvFQF9gXqpCZjYS+CPwKnASMAOYDPw8rtzXgOfxdoA6xS/zXWBWjuudK88Ch8d9LIwtYGYXAzPxrv9E4FG8DTMuDbeqLePvP/0isC/edpjjgK8AC8ysSyHrlgOX0/RndnzwgJkZMB/vZ3UZcBZQgnfde4Rf1bSG4G3iUgMsTlHu93j/X6/F+7/1IfBsfCABNwBTgDvw/q8uAR41s5NzW+0Wy/Q6wft7Ef//8l9xZaJ6nfnhnNNHG/kAOvifOwEOmJKk3OvAwrhj1wJfAv1jjj0OvAeUxBw73z/3IYW+3rj6rwbmpCnTCe9NxH1xx+/BW02nJF/1y8H1TQC2AUNijg0CGoArC12/LK/pGP936fgUZUb7ZY6NObYr3haevyr0NSSob4eYf1/k1708rsxB/vHvxhzrhLepy5Mxx/oCdcD1cc9/Afh/Ub9O/zEH3JjmXJG9znx9qMXbhjjntqcrY2Z7AgcDc+Iemo3XkjjJL1eC18qY65yrjyk3Fy+gR+eiziE7HOhD4mvvBRwZeo0ydxqwxDm3IjjgvD2mX6I4fxaZOg34t3NuQXDAOfcpXis4ctedyf9BvGuqBx6JeV4D8DAw0sxK/cMjgZ1o/vs6BzjQzAa1vsbZyfA6MxXZ68wXBW/7M9T//GbsQf+P+BZgf//QYGDnBOW2AitjykXJqf54U52ZLUkwvpvw2oHl/ucoXlNgKM3rDV7do1zvTDxgZtvM7GMzezBu3DrVde9lZl3DqWJODQWqnHNb4o4vxwugITHl6oAVCcpB8fzcL/X/T24xsxfNbHjc423lOjOm4G1/evqfaxI8VhPzeKpy1TGPR8V8vDHAkcAYYCvwuJmNjSmT7Jqq4x6Pop4k/1n0CLkuufIpcBteV+UIvHG+44FXzKyvXybVdUNxXnu6a4r9P/iJ8/tdU5SLsjnAD/F+rpfg9Sy9aGbHxJRpC9fZIp0KXQFJzMyOx5vYlM5C59wxea5OqLK5dufcZXHneBxvgsbNNO/Ckghwzr2ON98gsNDMFgH/wJtwNbkgFZOccc6Ni/lysZk9gdeDcSPRHtrJKwVvdL0M7JdBufjuqnSCd9qJWgo92PEuM1W5nuzoBsqHVl+7c26bmT0K/MLMdnPOfUjTa/owpnjwjrqa6Koh+c8iUeupKDnnXjOzfwHD/EOprjt4vNjUAAMTHI//PawBupuZxbUGi+H3NSHn3Odm9hRwYczhNned6Sh4I8of/3knD6cOAnMo8Epw0L8HrzPwln9oJd64y9CY52JmOwN7492Gkxd5uPbgP3PstccGbzCG9BbRtZy4n4Vvf6Jd72zF/sxOSPD4/sD7zrnN4VUpZ5YDZ5hZ57hx3v3xJi6uiClXijffYkVcOSjun3tswLbl60xIY7ztjHPufeANvHHQWGPxZlo+7Zf7EngGOMfMYt+gfRvvP8mT+a9t9vw6n4v3x3m9f/gVvNuGEl17Nd4M4ah6EjjMzPYODvhvlr5FxH8WLWFmFcA+eN3N4F3bADM7OqbMLsCpFO91z8e7g+Ds4EDM7+tzzrk6//AzeP8nE/2+vulPiCwq/s/uFHb8fKENXmc6avG2If4frXJ2vKHa38y+7f/7/7d3fyFSlWEcx7+/tqAMDYL+XaRbUDddRFtGgf0BL6KLIsMKqSiijExByC42KoqgiwKxItCM8kKJ/kB/8UrRLqQiKxKRhGDbbCWTYmm9ULb26eI5E4dhtx2rec/M8vvAYc/svHt4n2VmnnfOe87zbq+Nrp8APpG0CXgLuIKcT3uplqQgb2j/HHhH0qvVsV8E3ouIr7oYykmRtIK8tWQ7cAg4D3gUGAJWtNpFxKSkp8iCGWPADvKingeANdVgo1dtBlYDH0p6kvzG8BwZ76YmO/ZvSdoGjABfA+Pk63AYGANerpp9RA6Ytkp6nDwtOQwIeKF0nztRe89dWf28WdJR4GhEfBoR30h6G9hQ3bY3AjxC3pf9d/KJiF8krQeGJU2Q/6e7yNfsrYXCmdFscUpaRw6idgGHydPr64Dz6aM4u6LpG4m9/X8bWSEmZtgG29reTn7zPQH8SBbQGJjmmNeTH3zHgSPABmBe07G29fEasqrTEXLkPE4m1ZtmaP8wWTnnBFkgZFXTMXQY50Ky4tbvwARZoWyw6X79h3iGgX3k1c2T5CDiNeCCtnZnk0VOfiPn9XcClzfd/3+Ia6b34O5amzOA9cDP1XvrC+DGaY41QA6KR6vX6z5gedMxdhIneVZiD3mWaRL4lRxIXd1PcXZj87KAZmZmBXmO18zMrCAnXjMzs4KceM3MzApy4jUzMyvIidfMzKwgJ14zM7OCnHjNDABJ90uK2jYh6VtJq+vVyySdJmmVpD2Sxqsl30YkvSFpqNZuiaQtkvZL+kPSD40EZtZjXLnKzNrdAfwELKj2XwHOBZ6WdCZZVnQxsBF4HjhGriF7D1ncorWowVLgOmAvWVhhfrkQzHqXC2iYGZDfeIE3gUsi4vva73cBQxFxlqTXgXvJKkufTXOMZRHxfrV/SkRMVftbgSURMdj9SMx6m081m9lsvgQWSFoE3Adsni7pArSSbrU/Vah/Zn3FidfMZnMR8Ce5EtKp9O+qQGY9wXO8ZtZuoLqYaj5wJ7mgxsfAhdXzo011zGwucOI1s3bf1fangG3AWuChZrpjNrc48ZpZu2XkVc0TwGhEHAeQdKh6fhFwsKG+mfU9z/GaWbv9EbE3Ig62km5lNznXe0sz3TKbG5x4zawjEXEY2AKslHTtdG0k3Va0U2Z9yKeazexkrAUuBXZK2gjsIAtoXAzcDVwFfAAg6RzghurvFgLzJC2vHh+IiAMlO27WK5x4zaxjEXFM0lJgJZloHwROB8bIqlWP1ZpfBrzbdojW42eBZ7raWbMe5cpVZmZmBXmO18zMrCAnXjMzs4KceM3MzApy4jUzMyvIidfMzKwgJ14zM7OCnHjNzMwKcuI1MzMr6C8+n2NuxSv7HwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 504x504 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Based on the PCs we found to maximize Q2Y \n",
    "# we use that to continue our analysis. Here we construct \n",
    "# a scores plot\n",
    "plsda = PLSRegression(n_components = min_num_pcs)\n",
    "plsda.fit(x_zscore,y_zscore)\n",
    "print(\"Verifying X_scores shape: {}\".format(plsda.x_scores_.shape))\n",
    "\n",
    "# Extract PC1 and PC2\n",
    "scores_PC1 = plsda.x_scores_[:, 0]\n",
    "scores_PC2 = plsda.x_scores_[:, 1]\n",
    "\n",
    "# Since we know the number of observations are the same for mesenchymal \n",
    "# and pronerual the number of observations for each is half of the total \n",
    "num_observations_for_class = X.shape[0] // 2\n",
    "\n",
    "# Plot scores, separating into Mes and Pro classes based on label. \n",
    "# Because we know the order in which we structured the CSV, we know \n",
    "# the first half is mesenchymal and the second is proneural\n",
    "plt.figure(figsize =(7,7))\n",
    "plt.title(\"Scores Plot {} x {}\".format(IMAGE_DIM, IMAGE_DIM))\n",
    "plt.xlabel(\"PC1\")\n",
    "plt.ylabel(\"PC2\")\n",
    "plt.scatter(scores_PC1[0:num_observations_for_class - 1], \n",
    "            scores_PC2[0:num_observations_for_class - 1], \n",
    "            c='tab:orange',marker='D', label=\"Mesenchymal\")\n",
    "plt.scatter(scores_PC1[num_observations_for_class:], \n",
    "            scores_PC2[num_observations_for_class:], \n",
    "            c='tab:cyan',marker='D', label=\"Proneural\")\n",
    "plt.legend()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Assessing PLSDA Prediction Performance\n",
    "\n",
    "After doing some research online we discovered a method for implementing PLSDA via PLSR. We will build a model using PLSR, where our y_observed is an array of 0s (mesenchymal) or 1s (proneural). PLSR will make predictions which will result in continous values. However, we need our data to be discrete so that we can correctly classify each image into a cell type. Thus, following PLSR, we will need to threshold our y_predicted."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Leave One Out Cross Validation\n",
    "Perform LOO cross validation on our model to evaluate the prediction error of our model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Separate data and train"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "num_of_comp = opt_num_components   # The number of PCs we will use for PLS\n",
    "loo = LeaveOneOut()\n",
    "\n",
    "pls = PLSRegression(n_components=(num_of_comp))\n",
    "\n",
    "Y_pred = [] # Stores predictions on the test set based on the \n",
    "            # trained model \n",
    "Y_test = [] # Stores the true labels. This list is in the same \n",
    "            # order as Y_pred so we know the order in which we\n",
    "            # should compare\n",
    "\n",
    "# Partition into training and testing sets, taking care to normalize \n",
    "# only after separating the data to avoid leakage\n",
    "for train_index, test_index in loo.split(X, Y):\n",
    "    # Normalize X data\n",
    "    pls.fit_transform(X[train_index], Y[train_index])\n",
    "    pls.transform(X[test_index])\n",
    "    \n",
    "    # Predict\n",
    "    new_y_pred = pls.predict(X[test_index])\n",
    "    \n",
    "    # Store predictions and test order\n",
    "    Y_pred.append(np.squeeze(new_y_pred))\n",
    "    Y_test.append((Y[test_index]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Apply thresholds and build ROC curve"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def classify_and_eval_preds(thresh, num_inputs, Y_raw_preds):\n",
    "    '''Classify predictions based on the threshold'''\n",
    "    Y_class_preds = []\n",
    "    for i in range(num_inputs):\n",
    "        if Y_raw_preds[i] >= thresh:\n",
    "            Y_class_preds.append(1)\n",
    "        else:\n",
    "            Y_class_preds.append(0)\n",
    "    return np.squeeze(np.array(Y_class_preds))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Threshold our results to classify our data into two categories: \n",
    "# (0) mes and (1) pro. Then evaluate the prediction error for our \n",
    "# model. We will do this for a range of threshold values from -1 to 1\n",
    "# and plot the resulting ROC curve\n",
    "num_observations = X.shape[0] \n",
    "num_mes_labels = num_pro_labels = num_observations // 2 # We know that \n",
    "                                                        # there are equal \n",
    "                                                        # inputs of \n",
    "                                                        # mesenchymal and \n",
    "                                                        # proneural\n",
    "thresh_ranges = np.linspace(-2, 2, 100) # Apply a threshold range from -1 to \n",
    "                                        # 1 because that is the range for the \n",
    "                                        # zscored y labels\n",
    "mes_tpr_thresh_range = []   # Mesenchymal true positive rates  \n",
    "mes_fpr_thresh_range = []   # Mesenchymal false positive rates\n",
    "pro_tpr_thresh_range = []   # Proneural true positive rates\n",
    "pro_fpr_thresh_range = []   # Proneural false positive rates\n",
    "                            \n",
    "for t in thresh_ranges:\n",
    "    Y_class_pred = []   # Y_pred converted to classes 0 and 1\n",
    "    Y_class_test = []   # Stores the true labels. This list is in the same \n",
    "                        # order as Y_class_pred so we know the order in which \n",
    "                        # we should compare\n",
    "            \n",
    "    # For each observation, apply threshold and classify\n",
    "    Y_class_pred = classify_and_eval_preds(t, num_observations, Y_pred)\n",
    "\n",
    "    # Quantify correct predictions made by our model and compute true/false \n",
    "    # positive rates for the ROC curve. \n",
    "    mes_true_pos = 0\n",
    "    mes_false_pos = 0\n",
    "    pro_false_pos = 0\n",
    "    pro_true_pos = 0\n",
    "    \n",
    "    for i in range(num_observations):\n",
    "        if (Y_class_pred[i] == Y_test[i]):\n",
    "            if Y_class_pred[i] == 0:\n",
    "                mes_true_pos += 1\n",
    "            else:\n",
    "                pro_true_pos += 1\n",
    "        else: \n",
    "            if Y_class_pred[i] == 0:\n",
    "                mes_false_pos += 1\n",
    "            else:\n",
    "                pro_false_pos += 1\n",
    "        \n",
    "    mes_tpr = mes_true_pos / (num_mes_labels)\n",
    "    mes_fpr = mes_false_pos / (num_pro_labels)\n",
    "    pro_tpr = pro_true_pos / (num_pro_labels)\n",
    "    pro_fpr = pro_false_pos / (num_mes_labels)\n",
    "\n",
    "    mes_tpr_thresh_range.append(mes_tpr)\n",
    "    mes_fpr_thresh_range.append(mes_fpr)\n",
    "    pro_tpr_thresh_range.append(pro_tpr)\n",
    "    pro_fpr_thresh_range.append(pro_fpr)\n",
    "\n",
    "    if (LOG_FLAG):\n",
    "        print(\"LOO Cross Validation Results [128 CASE]:\")\n",
    "        print(\"mes cells correctly predicted: \", mes_correct)\n",
    "        print(\"mes cells incorrecly predicted: \", mes_false_pos)\n",
    "        print(\"pro cells correctly predicted: \", pro_correct)\n",
    "        print(\"pro cells incorrectly predicted: \", pro_false_pos)\n",
    "\n",
    "# Plot data\n",
    "plt.figure(figsize =(7, 7), dpi=100)\n",
    "plt.title(\"Proneural ROC {}x{}\"\n",
    "          .format(IMAGE_DIM, IMAGE_DIM))\n",
    "plt.plot([0,1], [0,1], label=\"Random\")\n",
    "plt.plot(pro_fpr_thresh_range, pro_tpr_thresh_range, label=\"Model\")\n",
    "plt.xlabel(\"False Positive Rate\")\n",
    "plt.ylabel(\"True Positive Rate\")\n",
    "plt.legend()\n",
    "plt.show()\n",
    "\n",
    "# Compute AUC\n",
    "from sklearn.metrics import auc\n",
    "model_auc = auc(pro_fpr_thresh_range, pro_tpr_thresh_range)\n",
    "print(\"AUC: {}\".format(model_auc))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Compute ROC and AUC for 8-fold crossvalidation\n",
    "\n",
    "We use the data from what was computed for the optimal run in the above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.metrics import auc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "num_mes_labels = num_pro_labels = num_observations // 2 # We know that \n",
    "                                                        # there are equal \n",
    "                                                        # inputs of \n",
    "                                                        # mesenchymal and \n",
    "                                                        # proneural\n",
    "        \n",
    "thresh_ranges = np.linspace(-2, 2, 100) # Apply a threshold range from -1 to \n",
    "                                        # 1 because that is the range for the \n",
    "                                        # zscored y labels\n",
    "mes_tpr_thresh_range = []   # Mesenchymal true positive rates  \n",
    "mes_fpr_thresh_range = []   # Mesenchymal false positive rates\n",
    "pro_tpr_thresh_range = []   # Proneural true positive rates\n",
    "pro_fpr_thresh_range = []   # Proneural false positive rates\n",
    "                            \n",
    "for t in thresh_ranges:\n",
    "    Y_class_pred = []   # Y_pred converted to classes 0 and 1\n",
    "    \n",
    "    # For each observation, apply threshold and classify\n",
    "    Y_class_pred = classify_and_eval_preds(t, num_observations, opt_Y_pred)\n",
    "\n",
    "    # Quantify correct predictions made by our model and compute true/false \n",
    "    # positive rates for the ROC curve. \n",
    "    mes_true_pos = 0\n",
    "    mes_false_pos = 0\n",
    "    pro_false_pos = 0\n",
    "    pro_true_pos = 0\n",
    "    \n",
    "    for i in range(num_observations):\n",
    "        if (Y_class_pred[i] == opt_Y_test[i]):\n",
    "            if Y_class_pred[i] == 0:\n",
    "                mes_true_pos += 1\n",
    "            else:\n",
    "                pro_true_pos += 1\n",
    "        else: \n",
    "            if Y_class_pred[i] == 0:\n",
    "                mes_false_pos += 1\n",
    "            else:\n",
    "                pro_false_pos += 1\n",
    "        \n",
    "    mes_tpr = mes_true_pos / (num_mes_labels)\n",
    "    mes_fpr = mes_false_pos / (num_pro_labels)\n",
    "    pro_tpr = pro_true_pos / (num_pro_labels)\n",
    "    pro_fpr = pro_false_pos / (num_mes_labels)\n",
    "\n",
    "    mes_tpr_thresh_range.append(mes_tpr)\n",
    "    mes_fpr_thresh_range.append(mes_fpr)\n",
    "    pro_tpr_thresh_range.append(pro_tpr)\n",
    "    pro_fpr_thresh_range.append(pro_fpr)\n",
    "    \n",
    "    if (LOG_FLAG):\n",
    "        print(\"{}-Fold Crossvalidation Results [Thresh: {}]:\".format(group_size, t))\n",
    "        print(\"mes cells correctly predicted: \", mes_correct)\n",
    "        print(\"mes cells incorrecly predicted: \", mes_false_pos)\n",
    "        print(\"pro cells correctly predicted: \", pro_correct)\n",
    "        print(\"pro cells incorrectly predicted: \", pro_false_pos)\n",
    "        print(\"\\n\\n\")\n",
    "\n",
    "\n",
    "# Plot data\n",
    "plt.figure(figsize =(7, 7), dpi=100)\n",
    "plt.title(\"Proneural ROC {}x{}\"\n",
    "          .format(IMAGE_DIM, IMAGE_DIM))\n",
    "plt.plot([0,1], [0,1], linestyle='--', lw=2, label=\"Chance\")\n",
    "plt.plot(pro_fpr_thresh_range, pro_tpr_thresh_range, label=\"Model\")\n",
    "plt.xlabel(\"False Positive Rate\")\n",
    "plt.ylabel(\"True Positive Rate\")\n",
    "plt.legend()\n",
    "plt.show()\n",
    "\n",
    "# Compute AUC\n",
    "model_auc = auc(pro_fpr_thresh_range, pro_tpr_thresh_range)\n",
    "print(\"AUC: {}\".format(model_auc))\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}