principal_component_analysis.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "### Import Libraries"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "import matplotlib\n",
    "import matplotlib.pyplot as plt\n",
    "matplotlib.style.use('fivethirtyeight')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "from sklearn.datasets import fetch_olivetti_faces\n",
    "\n",
    "import time\n",
    "import timeit\n",
    "%matplotlib inline\n",
    "\n",
    "import scipy\n",
    "import scipy.stats\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "from ipywidgets import interact\n",
    "\n",
    "plt.style.use('fivethirtyeight')\n",
    "%matplotlib inline\n",
    "\n",
    "from numpy.testing import assert_allclose"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "**Compute the sample mean vector of a matrix**\n",
    "\n",
    "$\\overline{X}$ = $\\displaystyle\\frac{1}{N}$ (${X}^T$ ${j}$), where \n",
    "- ${j}$ $\\in$ $\\mathbb{R}^{n \\times 1}$\n",
    "- $\\overline{X}$: Sample mean vector\n",
    "- ${X}^T$: Transpose of the input matrix of shape (N, D)\n",
    "- ${N}$: Size of the dataset\n",
    "- ${D}$: Dimensionality of the dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "def mean_x(X):\n",
    "    N, D = X.shape\n",
    "    j = np.ones(N)\n",
    "    XTJ = X.T @ j\n",
    "    mean_vector = XTJ / N\n",
    "\n",
    "    return mean_vector"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "**Compute the covariance matrix of a dataset**\n",
    "\n",
    "${S}$ = $\\displaystyle\\frac{1}{N}$ ${X}^T$ (${I}$ - $\\displaystyle\\frac{1}{N}$ ${J}$) ${X}$, where\n",
    "- ${X}$: Input matrix of shape (N, D)\n",
    "- ${N}$: Size of the dataset\n",
    "- ${D}$: Dimensionality of the dataset\n",
    "- ${X}^T$: Transpose of the input matrix\n",
    "- ${I}$ $\\in$ $\\mathbb{R}^{n \\times n}$: Identitiy matrix\n",
    "- ${J}$ $\\in$ $\\mathbb{R}^{n \\times n}$: Unit matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "def cov(X):\n",
    "    N, D = X.shape\n",
    "    I = np.identity(N)\n",
    "    J = np.ones((N,N))\n",
    "    S = (1/N) * X.T @ (I - (1/N) * J) @ X\n",
    "\n",
    "    return S"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "**Data normalization**\n",
    "\n",
    "$\\bar{X}$ = $\\frac{X - μ}{σ}$, where\n",
    "- $\\bar{X}$: Normalized dataset with mean zero and unit variance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "def normalize(X):\n",
    "    μ = mean_x(X)\n",
    "    X_bar = X - μ\n",
    "\n",
    "    return X_bar, μ"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "**Compute eigenvalues and eigenvectors**\n",
    "\n",
    "Where\n",
    "- ${S}$: Covariance matrix\n",
    "- ${\\lambda}$: Eigenvalues sorted in descending order of the eigen values\n",
    "- ${V}$: Eigenvectors sorted in descending order of the eigen values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "def eig(S):\n",
    "    λ, V = np.linalg.eig(S)\n",
    "    sort_indices = np.argsort(λ)[::-1]\n",
    "\n",
    "    return λ[sort_indices], V[:, sort_indices]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "**Compute Projection Matrix onto the space spanned by `B`**\n",
    "\n",
    "${P}$ = ${B}$ (${B}^T$ ${B})^{-1}$ ${B}^T$, where\n",
    "- ${B}$: Basis of the subspace with dimension (D, M)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "def projection_matrix(B):\n",
    "    P = B @ np.linalg.inv(B.T @ B) @ B.T\n",
    "\n",
    "    return P"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "**Principal Component Analysis**\n",
    "\n",
    "\n",
    "The function below takes advantage of the functions implemented above to perform PCA."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "def PCA(X, num_components):\n",
    "    # Normalization of the dataset so that they have zero mean\n",
    "    X_normalized, μ = normalize(X)\n",
    "\n",
    "    # Compute covariance matrix S\n",
    "    S = cov(X)\n",
    "\n",
    "    # Compute the eigenvalues (λ) and eigenvectors (V) of S in descending order of the eigenvalues\n",
    "    λ, V = eig(S)\n",
    "\n",
    "    # Take the top \"num_components\" of λ and V\n",
    "    principal_values = λ[:num_components]\n",
    "    principal_components = V[:, :num_components]\n",
    "\n",
    "    # Compute the projection matrix onto the space spanned by the principal components\n",
    "    P = projection_matrix(principal_components)\n",
    "\n",
    "    # Reconstruct data using the basis spanned by the principal components\n",
    "    X_reconsted = (P @ X_normalized.T).T + μ\n",
    "\n",
    "    return X_reconsted, μ, principal_values, principal_components"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "Generate a dataset X from a 2D Gaussian distribution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "def draw_vector(v0, v1, ax=None, label=None):\n",
    "    \"\"\"Draw a vector from v0 to v1.\"\"\"\n",
    "    ax = ax or plt.gca()\n",
    "    arrowprops=dict(arrowstyle='->',\n",
    "                    linewidth=2,\n",
    "                    shrinkA=0, shrinkB=0, \n",
    "                    color='k')\n",
    "    ax.annotate('', v1, v0, arrowprops=arrowprops, label=label)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "D = 2\n",
    "N = 10\n",
    "\n",
    "mvn = scipy.stats.multivariate_normal(\n",
    "    mean=np.ones(D, dtype=np.float64), \n",
    "    cov=np.array([[1, 0.8], [0.8, 1]], dtype=np.float64)\n",
    ")\n",
    "\n",
    "X = mvn.rvs((100,), random_state=np.random.RandomState(0))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "**Compare the results of our PCA implementation with the implementation of scikit-learn**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "from sklearn.decomposition import PCA as SKPCA"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "difference in reconstruction for num_components = 1: 3.404427844094429e-29\n",
      "difference in reconstruction for num_components = 2: 1.1327549560907966e-29\n"
     ]
    }
   ],
   "source": [
    "for num_component in range(1, 3):\n",
    "    # We can compute a standard solution given by scikit-learn's implementation of PCA\n",
    "    pca = SKPCA(n_components=num_component, svd_solver=\"full\")\n",
    "    sklearn_reconst = pca.inverse_transform(pca.fit_transform(X))\n",
    "    reconst, _, _, _ = PCA(X, num_component)\n",
    "    # The difference in the result should be very small (<10^-20)\n",
    "    print(\n",
    "        \"difference in reconstruction for num_components = {}: {}\".format(\n",
    "            num_component, np.square(reconst - sklearn_reconst).sum()\n",
    "        )\n",
    "    )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "**As we can see above, our almost achieves identical results as the scikit-learn implementation.**\n",
    "\n",
    "**We can also visualize PCA:**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1AAAALjCAYAAAAcMuReAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAABYlAAAWJQFJUiTwAACgVElEQVR4nOzdfZzUZb3/8dc197uywMKCd+yKICa5Bi7eAQbY3bKY2OmknpOVdsJzqKDyV2ZZlnUK7dZKO1Fomud4TMtTgbGIlYFyY8IqioKiiOyqJQsrIjv3c/3+GHbZ2dvvzM7s3Oz7+XjwQGe/O99rZnaZ73uuz/W5jLUWERERERERGZgr3wMQEREREREpFgpQIiIiIiIiDilAiYiIiIiIOKQAJSIiIiIi4pAClIiIiIiIiEMKUCIiIiIiIg4pQImIiIiIiDikACUiIiIiIuKQApSIiIiIiIhDClAiIiIiIiIOefI9gGJ38OBBm+8xiIiIiIjI4IwaNco4OU4zUCIiIiIiIg4pQImIiIiIiDikACUiIiIiIuKQApSIiIiIiIhDClAiIiIiIiIOKUAJALt27WLXrl35HoYMgl7D4qfXsPjpNSx+eg2Ln17D0lDIr6MClIiIiIiIiEMKUCIiIiIiIg4pQImIiIiIiDikACUiIiIiIuKQApSIiIiIiIhDClAiIiIiIiIOefI9gOEqGo1y6NAhgsEgiUQi38PB7XYD0NzcnOeRSKaG02vocrkoKyujoqICr9eb7+GIiIjIMKIAlQfRaJTXX3+diooKjj32WNxuN8aYvI4pFAoBEAgE8joOydxweQ2ttcTjcdrb23n99dcZP368QpSIiIgMGZXw5cGhQ4eoqKhg5MiReDyevIcnkWJijMHj8TBy5EgqKio4dOhQvockIiIiw4gCVB4Eg0HKy8vzPQyRoldeXk4wGMz3MERERGQYUYDKg0Qi0bleRUQy53a7C2INoYiIiAwfClB5orI9kcHT75GIiIgMNQUoERERERERhxSgREREREREHFKAkoJ0xhlncMYZZ+R7GCIiIiIiKRSgpGSNHj2aCy+8MN/DEBEREZESogAlIiIiIiLikAKUiIiIiIiIQwpQJchay5Z9Ea54eD/H3/UqlXe8wvF3vcqVD+9n674I1tp8DxFIjvMXv/gF5513HsceeyxTp07lmmuu4eDBgz2OPXjwID/5yU+46KKLePvb3864ceOYPHky//Iv/8Lf/va3lGPvvvtuRo8eDcCGDRsYPXp0558bb7wx5biPfvSjTJs2jeOOO47q6mrq6+u59957c/q4RURERKR4efI9AMmuaMKyeH0bjc0hQjFLxxajwbhl5Z4Qa1vCNFQHWD6nEq8rv3vofOlLX+LnP/85xx13HFdccQVer5fVq1ezZcsWotEoXq+389jnn3+e//zP/2TWrFm8733vY/To0bS0tNDY2Mif/vQnfv3rX/Oe97wHSDaguPbaa/nOd75DdXU1H/7whzvv5/zzz+/8789//vOcdtppzJo1i+OOO44DBw7w0EMP8R//8R/s2rWLr371q0P3ZIiIiIhIUTCFMhtRrA4ePJj2E9jc3Ex1dXXWx2KtZdG6NlbvDRKM931cmRsW1JRx29zKzo1IQ6EQAIFAIOvj6s1jjz1GfX09J598Mn/5y1+orKzsHMdFF13E448/TnV1NU8//TSQnIGKxWKMHTs25X5eeeUV3v3udzNy5MgeM1GjR49m9uzZ/PGPf+x1DC+99BInn3xyym2RSIQPfehDbNy4kaeeeooTTjghWw8554b6NSwUufp9yoddu3YBMGXKlDyPRDKl17D46TUsfrl6Da21bG2Ncsv2Q6xtDhOKWwJuQ321n6W1FdRVebXBexbl43dx1KhRjl5AlfCVkK2tURqbQ/2GJ4BgHBqbQzS1RodmYL24++67geQsUEd4guTF/9e+9rUex48aNapHeAI48cQTWbhwIc8//zzNzc1pjaF7eALw+XwsWrSIWCzGunXr0ro/ERERKU3RRPJD6oVrWlm1J0QwbrEcrfC5aE0ri9a1EU1oYmI4UAlfCbl1+yFCMWe/uKGY5dbtb3HHBWNyPKrebdu2DUgtqeswc+ZM3G53j9s3b97M8uXLefzxx9m3bx+RSCTl66+99lpaMxHNzc38+Mc/Zt26dbS0tBAMBnvcn4iIiAxv1iaXR/RV4ZMA2mOW1XuDLF5PSoWPlCYFqBLyYHO4c83TQBLAgy2hXA6nX2+++SYA48aN6/E1j8fTY7Zp1apVXHHFFQQCAebNm8fJJ59MeXk5LpeLRx99lA0bNhAOhx2ff8+ePbzrXe/ijTfeYObMmVxwwQWMHDkSt9vN3r17ueeee9K6PxERESlNmVT4zBjnG5rBSV4oQJWQUDy9aeOgw9mqXBg5ciQA+/bt45hjjkn5WiwWY//+/Snrj5YtW4bP5+Phhx/mbW97W8rxn/vc59iwYUNa5//pT3/KgQMH+OlPf8rll1+e8rXf/va33HPPPWndn4iIiJSmYqrwkaGhNVAlJOBOb7q4zJO/6eVp06YB8Oijj/b42qZNm4jHUz/m2b17N29729t6hKdEIsHmzZt7PYfL5SKR6H1Obvfu3QAsXLiwx9fSDWMiIiJSuoqpwkeGhgJUCamv9jt+QV1A/YT8dWvraC3+gx/8gLa2ts7bQ6EQ3/zmN3scX1NTw+7du1PWJVlrufHGG9m5c2ev5xgzZgwtLS29fq2mpgboGeD+/Oc/c9ddd6X3YERERKRkFVOFjwwNlfCVkCW1FaxtCdPu4Bc34DYsqR0xBKPq3Xnnnce///u/84tf/IKZM2eycOHCzn2gRo8ezXHHHZdy/Kc+9Smuvvpq5syZw8KFC/F4PDz22GM899xzzJ8/nzVr1vQ4x9y5c7n//vu57LLLmDZtGl6vl1mzZjF79mw+8YlPcPfdd3PllVdy8cUXc9xxx7Fjxw7+9Kc/8U//9E/83//931A9FSIiIlLAAm5DMI0Qlc8KHxkamoEqITOqvDRUByjr2cAuRZkbGmoC1FV5+z8wx77zne/w3e9+l5EjR3LnnXdy//338+53v5vf//73KZvoAnz84x/npz/9Kcceeyz33HMPv/nNbzjxxBP505/+1FkO2N1NN93Ehz70IZqamvj+97/Pt7/9bdavXw9AbW0tq1at4pxzzuHBBx/kl7/8JYcOHeK///u/+fjHP57zxy4iIiLFoZgqfGRoaCPdQSqkjXQhuU/B4vVtNDaHCMVsSs2ui+SnKA01AZbPqcTrOvoJyXDdhLWUDNfXUBvpSiHRa1j89BoWv2y/hlv2RVi4ptVRhU+527CqoUpd+LJAG+nKkPG6DLfNrWTV/CoWTgxQ7jEYoNxjuHhiGQ80VHH7vDEp4UlEREREeldsFT6Se1oDVYKMMcwY5+POC8YOfLCIiIiI9MkYw/I5lSxej6MKH22iW/oUoERERERE+tFR4dPUGuWW7YdY2xImGLOUeQz1EwIsrR1Bncr2hg0FKBERERGRAajCRzpoDZSIiIiIiIhDClAiIiIiIiIOqYRPRERERCTXrMW1ewfexnvxbNsMkQj4fMSmzSTacBmJSaeBGlAUBQUoEREREZFcisXwr1iGp2kjRMOYjn1YI2E8W9bh2baZWN0swlddBx5dnhc6lfCJiIiIiOSKtUfC0wZMJHQ0PB1hrMVEQniaNuBfsQzswBv2Sn4pQImIiIiI5Ihr9w48TRsxkXC/x5lIGE/TRly7dw7RyCRTClAiIiIiIjnibbwPov2Hp07RMN419+V2QDJoClAiIiIiIjni2bapR9leX4y1eJ7clOMRyWApQImIiIiI5Eokkt7xTmerJG8UoESkoIwePZoLL7ww38MQERHJDp8vveO9/tyMQ7JGAUokz2688UZGjx7NI488ku+h9KrQxyciIlLIYtNmYh3u72SNITZ9Zo5HJIOlACUiIiIikiPRhkudzyp5/UTnX5rbAcmgKUCJiIiIiORIYtJUYnWzsL7+Q5T1+YnVzSIx6bQhGplkSgGqFFmL68Vn8d/6dY65qp5jrriAY66qx3/rDbhe3FEwG7S9/PLLjB49mk9+8pO88MILfPzjH+eUU06hsrKys1zsz3/+M5dccgmTJk1i/PjxTJ8+neuvv5433nij1/t85ZVX+OIXv0hdXR3HHXccEydO5F3vehff/e53exz75JNP8tGPfpRTTjmF8ePHU1tby+c//3n+/ve/9zj2k5/8JKNHj+bll1/mjjvuYNasWRx77LFMmTKFz372sxw8eLDH92zfvp1PfOITnHHGGYwfP57JkyczZ84cvvSlLxGNRgE444wz+M53vgPARRddxOjRozv/dD/3nj17+PnPf86sWbM47rjjOtcJ3X333YwePZpf//rXvT4nfa0pisfj/PKXv6S+vp6amhqOO+44zjzzTJYuXcqLL77oeHwA7e3t/PCHP+T888/nhBNO4MQTT+S9730vv/3tb3sdUyQS4bvf/S7Tp09n/PjxvOMd7+Bb3/oW4bAWzoqISIkxhvBV1xGrm431BXqU81ljsL4AsbrZhK+6DhyW+0n+ePI9AMmyWOzIbtcbIRo+2jYzEsazZR2ebZuJ1c1K/oJ6CuPl37NnD+9+97s55ZRTuOSSSwiFQlRUVHDTTTdx0003UVlZSX19PePGjeOZZ57hlltu4aGHHmLt2rWMHDmy836eeOIJPvjBD9LW1sasWbO46KKLaG9v57nnnuOmm27ii1/8Yuexa9as4WMf+xjWWi6++GKqq6t58sknuf3221m9ejWNjY1MnDixx1i//vWv85e//IX58+dzwQUX8Mgjj/CrX/2K3bt3s2rVqs7jtm/fznve8x6MMTQ0NHDSSSdx6NAhdu/eze23385Xv/pVvF4vn/zkJ/njH//Ihg0b+Nd//Vdqamr6fJ6uvfZaNm3aRH19Pe973/twu90ZP+eRSITLLruMhx9+mAkTJvChD32IiooK9u7dywMPPMB5553H5MmTHY3vjTfeYOHChTz11FNMmzaNyy+/HGstf/7zn1m0aBE7d+7kq1/9aufx1lquvPJKVq9ezcknn8xVV11FJBLhf/7nf3jmmWcyfkwiIiIFy+MhvPh6ort34m28F8+2zclue14/sekziTZcppmnIlIYV9CSHdYeCU8bet3t2lgLkRCepg2wYhnhxdcXxKccmzZt4v/9v//H1772tc7b1q9fz0033cQ555zDfffdlzLjcffdd/PpT3+aG2+8kRtvvBFIBoIrrriCtrY2VqxYwSWXXJJyjldeeaXzv9966y0++clPEovFeOCBB5g1a1bn1370ox9xww03cPXVV/O73/2ux1i3bNnChg0bqK6uBiAWi3HRRRfxyCOPsHXrVmbMmAHAPffcQygU4u677+4x+/PGG29QXl4OwKc+9SkOHjzIhg0b+PCHP8w73/nOPp+np556ivXr1/ca7NJ100038fDDDzN//nx+9atf4fcfLSsIh8McOnTI8fi+/OUv89RTT/GNb3yDz372s523h0IhLr/8cn7wgx+wcOFC3vGOdwDw29/+ltWrV3P22WezatUqAoEAANdddx0XXHDBoB+biIhIQTKGxOSphJfcgOotiptK+EqIa/cOPE0bew1PXZlIGE/TRly7dw7RyPo3fvx4rr322pTbfv7znwPw4x//uEe52OWXX84ZZ5zBb37zm87bGhsb2bt3Lw0NDT3CE8CJJ57Y+d+rV6+mra2ND37wgynhCWDJkiXU1NTw8MMP09zc3ON+vvjFL3aGJwCPx8Pll18OwNatW3scX1ZW1uO20aNH43Kl/6v3mc98JivhKR6Pc9ttt1FWVsbNN9+cEp4A/H4/VVVVju7rwIED3HfffZx55pkp4QkgEAhwww03YK1NKeW7++67Abj++us7wxNAZWUl11xzTaYPS0RERGRIaAaqhHgb73O++Vo0jHfNfYQ//fXcDsqB2traHhfxjz/+OF6vl9///vf8/ve/7/E90WiU1tZWDhw4wJgxY9iyZQsA733vewc837Zt2wCYM2dOj695PB5mzZrF3r17eeqpp1LCEsD06dN7fM+ECRMAUtZlffCDH2T58uVcfvnlLFy4kHnz5nHeeedx8sknDzi+vnTMbg3W888/z5tvvslZZ53F8ccfP6j7ampqIh6PY4zpnA3sKhaLAfDcc8913rZt2zZcLhczZ/Zs03r++ecPajwiIiIiuaYAVUI82zYdXfM0AGMtnic3FcQU8vjx43vcduDAAWKxWGcDg7689dZbjBkzprOJg5NA8OabbwJw7LHH9vr14447DqDXxhCjRo3qcVvHWqR4PN5524wZM2hsbOQHP/gBK1eu5N577wVgypQpXHvttXzoQx8acJzd9fY8ZSKd52ogBw4cAJJBqqmpqc/jDh8+3Pnfb775JpWVlXi93h7H9fWaiIiIiBQKBahSEomkd7zT2aocM72swxo5ciSJRII9e/Y4uo+OYPPaa68NeGxH44l//OMfvX69owtf1wYVmTjnnHO49957CYfDPPnkk/zpT3/iF7/4BYsWLaKqqop58+aldX+9PU9AZzlg1wDXobduhek8VwPpeI4+9alPsWzZMsff09bWRjQa7RGi+npNRERERAqF1kCVEp8vveOdbuqWB2effTZvvPEGO3bscHT8WWedBcBDDz004LEdzQweffTRHl+LxWJs2rQJgGnTpjkdbr/8fj/nnnsuX/nKVzpn1P74xz92fr1jBiuRSGR0/x1rxLo2yujw5JNP9rjt1FNPZdSoUTzzzDOOQlR/45sxYwYul6vzOXNi2rRpJBKJXr+nt9dEREREpJAoQJWQ2LSZPfYW6Is1htj0nmtQCsWnPvUpAD772c/2epF/+PBhHn/88c7/b2hooKamhsbGxl73HuoaLi688EIqKyu5//77U+4D4Gc/+xkvv/wy8+bN67H+KR2PPfYYwWCwx+379u0D6OzCBzBmzBiAXptWOHHmmWficrn43e9+R3t7e+ftbW1tKZ0NO7jdbhYtWkQwGOTqq6/usfdSJBKhtbXV0fjGjRvHJZdcwhNPPMF3v/vdXmfBXnrppZSZxI6mG9/61rcIhUIp4/3+97/v8FGLiIiI5IdK+EpItOHS5L4CkdDAB3v9ROdfmvtBZWju3LnccMMNfOMb32DGjBm8973v5aSTTuLw4cM0NzezYcMGzjvvPO6//34AfD4fv/rVr/inf/onFi1axB133MHZZ59NKBTi+eefZ926dezfvx+AESNGcOutt3LllVdy4YUX8oEPfIAJEybw5JNP8pe//IVjjz2WH/3oR4Ma/49//GMeeeQRZs6cyUknncQxxxzDjh07+NOf/sTo0aO58sorO4995zvficvl4pvf/CY7duzonFFy2pHuuOOO45JLLuHee+/lPe95D/X19Rw6dIiHHnqIWbNm8dRTT/X4nmuvvZYtW7awZs0aZsyYwfz58xkxYgSvvPIKf/nLX/jmN7/ZGXQGGt/3vvc9du/ezbJly7j33ns577zzGD9+PK+99hrPP/88TU1N3H777Z0dBD/0oQ/xf//3fzQ2NjJz5kwWLFhANBpl5cqVnHnmmbz00kuZPekiIiIiQ0ABqoQkJk0lVjerz32gOlifn1jdrILfsO1zn/sc5557Lj//+c/ZvHkzq1evZuTIkRx//PFcccUVPdqVn3nmmTzyyCP86Ec/4qGHHuJvf/sbI0aMYNKkSVx33XUpx1544YU8+OCD/OAHP+DPf/4zb775Jsceeyz/9m//xjXXXDPoBguLFi1i9OjRbN26lc2bNxOLxTjhhBP4xCc+0dkqvcPb3vY2fvazn3HLLbdw++23d87KpNPS+yc/+Qljxozh97//PbfddhsTJkzgP/7jP/jMZz7T635WPp+P+++/n1/+8pf8+te/5p577sFay3HHHcf73//+lA55A41v5MiR/PGPf+TOO+/kt7/9LatWrSIUCjF+/HgmTZrEsmXLUvZ3Msbwq1/9iptvvpn//d//ZcWKFRx77LFcfvnlfPGLX1QjCRERESloxjrs2ia9O3jwYNpPYHNz86DKw/oVix3ZTHcjRMMpXfmsMckdr+tmEb7qOvAczc8dF8Vd9+WR4jJcX8Oc/j4NsV27dgHJbo1SnPQaFj+9hsVPr2FpyMfrOGrUKEdrYTQDVWo8HsKLrye6eyfexnuTJX3RcDI4TZ9JtOGygp95EhEREREpVApQpcgYEpOnEl5yQ0Hs8yQiIiIiUirUhU9ERERERMQhBSgRERERERGHFKBEREREREQc0hooERERGfastWxtjXLjDh8b2tyEH32FgNtQX+1naW0FdVVejMPN6kWktClAiYiIyLAWTVgWr2+jsTlEMObGkgxKwbhl5Z4Qa1vCNFQHWD6nEq9LIUpkuFMJX55o/y2RwdPvkYgMlrXJ8LR6b5D2mO0MTx0SQHvMsnpvkMXr2/TvjogoQOWDy+UiHo/nexgiRS8ej+Ny6Z8xEcnc1tZocuZpgLflYBwam0M0tUaHZmAiUrB05ZEHZWVltLe353sYIkWvvb2dsrKyfA9DRIrYrdsPEYo5m1UKxSy3bn8rxyMSkUKnAJUHFRUVHDp0iDfffJNYLKZyAJE0WGuJxWK8+eabHDp0iIqKinwPSUSK2IPNYRIOj00AD7aEcjkcESkCaiKRB16vl/Hjx3Po0CH+8Y9/kEg4/ac7d6LRZEmC1+vN80gkU8PpNXS5XJSVlTF+/Phh8XhFJHdC8fQ+xAw6nK0SkdKlAJUnXq+XMWPG5HsYnXbt2gXApEmT8jwSyZReQxGR9AXchmAaIarMoy58IsOdSvhERERk2Kqv9ju+GHIB9RMCuRyOiBQBBSgREREZtpbUVhBwOKsUcBuW1I7I8YhEpNApQImIiMiwNaPKS0N1gDJ3/8eVuaGhJkBdldZdigx3ClAiIiIybBljWD6nkgU1ZZR7DIbU9VAuoNxtWFBTxvI5lRijNVAiw52aSIiIiMiw5nUZbptbSVNrlBs3vcajbW7CCUOZx1A/IcDS2hHUjfPle5giUiAUoERERGTYM8YwY5yPG6dGAJgyZUqeRyQihUolfCIiIiIiIg4pQImIiIiIiDikACUiIiIiIuKQApSIiIiIiIhDClAiIiIiIiIOKUCJiIiIiIg4pAAlIiIiIiLikAKUiIiIiIiIQwpQIiIiIiIiDilAiYiIiIiIOKQAJSIiIiIi4pAClIiIiIiIiEMKUCIiIiIiIg4pQB1hjPmOMebPxphmY0zQGHPAGPOEMebrxpix+R6fiIiIiIjknwLUUVcDxwAPAT8G7gZiwA3AU8aY6vwNTURERERECoEn3wMoICOttaHuNxpjvg1cB3wZ+NSQj0pERERERAqGZqCO6C08HXHfkb+nDNVYRERERESkMClADeyiI38/lddRiIiIiIhI3hlrbb7HUFCMMV8ARgCjgLOA80mGp/dYa/d1P/7gwYO9PoG7du3K5TBFRERERCQDU6b0Xlg2atQo4+T7tQaqpy8Ax3b5/zXAlb2FJxERERERGV40A9UHY8yxwCzgJqACeL+1tqn7cX3NQBWbjhmzvhK5FD69hsVPr2Hx02tY/PQaFj+9hqUhH6+j0xkorYHqg7X2H9ba3wHvA8YCd+V5SCIiIiIikmcKUAOw1r4MPAucboypyvd4REREREQkfxSgnDnhyN/xvI5CRERERETySgEKMMacaowZ1cvtriMb6Y4HNlpr24Z+dCIiIiIiUijUhS9pAXCjMeZR4CVgP8lOfHOBScDfgavyNzwRERERESkEClBJfwJOIbnn05nAaOAw8Dzw38BPrLUH8jY6EREREREpCApQgLV2O7Ak3+MQEREREZHCpjVQIiIiIiIiDilAiYiIiIiIOKQAJSIiIiIi4pAClIiIiIiIiEMKUCIiIiIiIg4pQImIiIiIiDikACUiIiIiIuKQApSIiIiIiIhDClAiIiIiIiIOKUCJiIiIiIg4pAAlIiIiIiLikAKUiIiIiIiIQwpQIiIiIiIiDilAiYiIiIiIOKQAJSIiIiIi4pAClIiIiIiIiEOefA9ARERECpu1lq2tUW7Zfoi1zWFCcUvAbaiv9rO0toK6Ki/GmHwPU0RkSChAiYiISJ+iCcvi9W00NocIxSyJI7cH45aVe0KsbQnTUB1g+ZxKvC6FKBEpfSrhExERkV5ZmwxPq/cGae8SnjokgPaYZfXeIIvXt2GtzccwRUSGlAKUiIiI9Gpra5TG5hDBeP/HBePQ2ByiqTU6NAMTEckjlfCJiIhIr27dfohQzNmsUihmuXX7W9xxwZgcj0o6aG2aSH4oQImIiEivHmwO9yjb60sCeLAllMvhSBdZXZtmLa7dO/A23otn22aIRMDnIzZtJtGGy0hMOg0UxEQ6KUCJiIhIr0Lx9NY0BR3OVsngdF2b1lt5ZeraNLhtbmXfM1GxGP4Vy/A0bYRoGNOxji0SxrNlHZ5tm4nVzSJ81XXg0WWjCGgNlIiIiPQh4E5v1qHMo1mKoZC1tWnWHglPGzCR0NHwdISxFhMJ4WnagH/FMlCTEBFAAUpERET6UF/td3yh4ALqJwRyORw5IpO1ab1x7d6Bp2kjJhLu9z5MJIynaSOu3TvTHqtIKVKAEhERkV4tqa0g4HBWKeA2LKkdkeMRCWRvbZq38T6I9h+eOkXDeNfc5/CsIqVNAUpERER6NaPKS0N1gDJ3/8eVuaGhJkBdlXdoBjbMZWttmmfbph5le30x1uJ5clNa5xUpVQpQIiIi0itjDMvnVLKgpoxyj+lx0eACyt2GBTVlLJ/TT6MCyaqsrU2LRNI7sdPZKpESpwAlIiIiffK6DLfNrWTV/CoWTgxQ7jEYoNxjuHhiGQ80VHH7vDEDt8qWrMna2jSfL70Te/3pHS9SotSPUkRERPpljGHGOB93XjA230MRkmvT1raEaXfQSKK/tWmxaTPxbFnnqIzPGkNs+sy0xypSihSgRERERIpIx9q0vvaB6tC5Nm2sB9eLz/bYKDd+yung9kLMQSmf1090/qXZexAiRUwBSkRERKSIdKxNW7w+uc9TKGZTuvK5SM48NdQEWD6rgsDPv9XrRrnuHU+AMVjjwti++/pZn59Y3SwSk07L6eMSKRYKUCIiIiJFpmNtWlNrlFu2H2JtS5hgzFLmMdRPCLC0dgR1Yz0Evv9F3M829RqQjLVgLdblwuICm6DrSjZrDHiT4Sl81XWgJiEigAKUiIiISFHqd21aLEbg+9fgfmYrA8Uek0hgvT7ip9TifvHZZLc9r5/Y9JlEGy7TzJNINwpQIiIiIqXEWvwrluF+9okBw1OnWBRbMYrDK9bkcmQiJUFtzEVERERKiGv3DjxNG/td19SdNsoVcU4zUCIiIiJ5ZK1la8dapuYwobgl4DbUV/tZWltBXZU3rU2KvY33ZbbprTbKFXFEAUpEREQkT6IJy+L1bT266QXjlpV7QqxtCdNQHWD5nErHmxV7tm1ytLdTD9ooV8QRlfCJiIiI5IG1yfC0em+Q9m6tyAESQHvMsnpvkMXr27BOQ1HEwb5O3ceijXJFHFOAEhEREcmDra1RGptD/W6GCxCMJ/d7amqNOrtjny/9wWijXBHHFKBERESk5Fhr2bIvwhUP7+f4u16l8o5XOP6uV7ny4f1s3RdxPpuTQ7duP0Qo5mwcoZjl1u1vOTo2Nm1mcg8nh6xxaaNckTRoDZSIiIiUlFysK8qFB5vDPcr2+pIAHmwJOTo22nApnm2bITLw8RaIv71OG+WKpEEzUCIiIlIycrauKAdC8fTOHXQ4W5WYNJVY3Sysr/+mENa4iNeeReia74FHn6mLOKUAJSIiIiUjZ+uKciDgTm/Gp8zj8HhjCF91HbG62VhfoEc5nzUG6wsQO3ceoatv0syTSJoUoERERKRk5GpdUS7UV/sdX4i5gPoJAed37vEQXnw9wS/dTOysuZ1ByvoCxM6eR/DLPyL8ya9p5kkkA/qtERERkZKRq3VFubCktoK1LWHaHQS+gAu+XNGM/9YfH1nfFAGfj9i0mUQbLks2gOg+k2QMiclTCS+5AW2RK5I9ClAiIiJSMnK1rigXZlR5aagOsHpvsN+SwwoTY9VLKzjz0S0QDR/dJDcSxrNlHZ5tm4nVzUo2gsjCjJK1lq2tUW7Zfoi1zWFCcUvAbaiv9rO0toK6Ki9GZX8yjKmET0REREpGztYV5YAxhuVzKllQU0a5x/S4KHMB5S5YtWcFs195HBMJHQ1PHfdhLSYSwtO0Af+KZTDIphjRhGXRujYWrmll1Z4QwbjFcrSD4UVrWlm0ro1oIv9t4EXyRQFKRERESkZO1xXlgNdluG1uJavmV7FwYoByj8EA5R7DxRPLWHd6K7Nf2YKJ9F+EZyJhPE0bce3emfFYiqmDoUg+KUCJiIhIyVhSW0HA4axSAli9N5j3zXWNMcwY5+POC8by6kdPoO3jJ/LqR0/gjgvGULv5dxB1uIIpGsa75r6Mx1FMHQxF8kkBSkREREpGx7qigMMrnHCClNK0mNMOFEPEs21Tj7K9vhhr8Ty5KeNzFVMHQ5F8UhMJERGRAqYF/ekxxnDL+aN57PUILYcHmEo5omtp2qFDPr71tkhuB5mOSJpjcTpb1Yti6mAokk8KUCIiIgUqmkiuSWlsDhHqsialY0H/2pYwDdUBls+pxOtSiIJk4Fz66BvsDzkLT10F47D+gJtn33Jxag7GlhGfDwZY/5TC68/4VMXUwVAkn1TCJyIiUoC0oD8zTtfx9CWcgLtfKZzPl2PTZmIdzjBaY4hNn5nxuYqpg6FIPhXOvxAiIiLSKZMF/TPG+YZmcHnWX1ljWzjheB1Pr/eN4ZED7iyOdnCiDZce2TjXQbmc1090/qUZn6u+2s/KPSFHZXyF0MFQJF80AyUiIlKAtKC/dwPtU7TutYjjdTx9CRdQI4nEpKnE6mZhff2X5lmfn1jdLBKTTsv4XOl0MAy4DUtqR2R8LpFipgAlIiJSgLSgvycnZY3Z4C+kqyNjCF91HbG62VhfoEc5nzUG6wsQq5tN+KrrYBANRTo6GJYNMAFX5oaGmgB1Vd6MzyVSzFTCJyIiUoC0oL+nwa5vcsJgeeeYHJ4gEx4P4cXXE929E2/jvcmSvmgYvH5i02cSbbhsUDNPHYwxLJ9TyeL19GhcAslP3QNuQ0NNsnGJuj/KcKUAJSIiUoACbkMwjRA1HBb0p1PWmCm/Cy4/MZbTc2TEGBKTpxJecgOZNyofmNdluG1uJU0da8xawgRjljKPoX5CgKW1I6gbJmvtRPqiACUiIlKASnlBf6Z7W6VT1piJMje8szLO20cU0CKoPDDGMGOcjzsvGJvvoYgUJAUoERGRArSktoK1LWHaHcy4FNOC/sHsbZVuWaNTXUvTvnD8/sEsIxKRYaCQlkmKiIjIEaW4oH+we1ulu08R9LzQMYDbJP8AlHsMF08s44GGKm6fNwaProxEZACagRIRESlApbigf7B7W6VT1miAOcf7qfQbreMRkaxSgBIRESlQpbagP5O9re64YEznbemUNZa5DV+bMXJINhfOdE2XiBQnBSgREZECVkoL+ge7t1VHWePqvcF+Z7GGsqyx3zVdLwU58MwzfHtfI+e92gSRCPh8xKZ1aT2uYCVSdBSgREREZEgMdm+rQitr7Lqmq3ug8yRi3LFzOQtbmwgkIhiOPJZIGM+WdXi2bSZWNyu5+a1Hl2MixURLJUVERGRIpNsEore9rTrKGlfNr2LhxADlHoOhZzOI7h38cqHPNV2JBCuf/h6XvL6ZYxJh3HQLgtZiIiE8TRvwr1gGtvQ3QRYpJfrIQ0RERIZEtva2KpSyxt7WdHkSMVY+/T3e27adgSKciYTxNG0kunsniclTczdQEckqzUCJiIjIkFhSW0Ggl1ml3hTD3lY91nQdmXl6j4Pw1CkaxrvmvuwPTkRyRgFKREREhkSp7W3VdU2XJxFj9dPf4b1t29O6uDLW4nlyU/YHJyI5owAlIiIiQ6KjCcSCmjLKPabHRYgLKHcbFtSUFcXeVp1rujKZeeoqGs7msEQkxxSgREREZMgUUhOIwaqv9uPLcOYphdefzWGJSI6piYSIiIgMqUJpAjFYS04fwWWN3+ddbc9kNvMEWGOITZ+Z1XGJSG4pQImIiMiArLVsbY1yy/ZDrG0OE4pbAm5DfbWfpbUV1FV5C77kLqus5dwdD+F7/bEebcrT4vUTnX9p9sYlIjmnACUiIiL9iiaSG8Z237w2GLes3BNibUuYhurk5rXFUHo3aLEY/hXL8Pztr2CdNGXvnTUuYnWzSEw6LXtjE5Gc0xooERER6ZO1yfC0em+Q9i7hqUMCaI9ZVu8Nsnh9G7bUN4W1NhmemjZgEonMS/eA+Ol1hK+6DobTzJ1ICdAMlIiIiPRpa2uUxuYQwXj/xwXj0Ngcoqk1yoxxvqycuxDLBl27d+Bp2oiJZN45zwLx2rMIfeF7Ck8iRUgBSkREpIDlO0Tcuv0QoZizWaVQzHLr9re444Ixgz5voZYNehvvG1Tb8eTM01mErr5J4UmkSKmET0REpEBFE5ZF69pYuKaVVXtCBOMWy9EQcdGaVhatayOayF3Z3IPN4R5le31JAA+2hAZ9zkIuG/Rs24TJ8HydM0/XfA88+gxbpFgpQImIiBSgQgkRoXh69xt0OFvVn0zKBnPGWlwvPov/1q9zzFX1kGHpnmaeREqHPv4QEREpQPlce9RVwG0IphGiyjyDDwf5KhvsIRolcPOXce9ogsE0jDAu4qfXac2TSInQDJSIiEgByiRE5EJ9td/xxYILqJ8QGPQ581E22EM0Svm1l+N+ZkvG3fYsR1qVnztPM08iJUQzUCIiIgWoIEIEsKS2grUtYdodhLmA27CkdsSgz5mPssEUiQSBby/F7H8941knAFwuwv92DbF3NmRrZCJSADQDJSIiUoDyHiKOmFHlpaE6QJm7/+PK3NBQE6Cuyjvocwbc6cWWbJQNdorHmHzPj3G/tHNQ4cn6/MTOmUfs/PlZG5qIFAYFKBERkQKU1xDRhTGG5XMqWVBTRrnH9LhwcAHlbsOCmjKWz6nMSkv1fJQNAmAtJ638JRV7dgxivZPB+gLE6mZrk1yREqUSPhERkQJUX+1n5Z6QozK+rIaIXnhdhtvmVtLUsR9VS5hgzFLmMdRPCLC0dgR1WWxgkY+yQUhukjvq+W0ZtSm3AL4AsekziTZcRmLSaVkZk4gUHgUoERGRApSvENEXYwwzxvm484KxOT0PHC0bXL032G8XwmyWDUJyk1xXLMOW6L4Ah1esyco4RKSwqYRPRESkAOVj7VGhGPKywSN7PXm2rM949ik2febgxpDuOa1ly74IVzy8n+PvepXKO17h+Lte5cqH97N1X2RINxcWGW40AyUiIlKAOkLE4vXJfZ5C3TbTdZGceWqoCWRt7VEhGbKywVgM/4pleJo2gnXa97Abl4vo/EsHPxaHoonkJsvdfy6CccvKPSHWtoRpqE7+XHhdpfVzIVIIFKBERGTYs9ayteNCvTlMKG4JuA311X6W1lZQV+XNS0AZ6rVHhSbnZYPWHglPGzCRcGZ3AcSn1g3Zmidrk+Gpr/LGBNAes6zeG2TxerhtbumFa5F8U4ASEZFhrdA/zR/KtUfDhrW4du/A95sVuHc8kXnHPcCOHU/o6huHrNve1tYojc2hfteGAQTjyZnLptYoM0o4ZIvkg9ZAiYjIsNX10/z2biVy0P3T/DatKykFsRj+5f9J2U3/b9DhKT7pNNq/+7/gHbr1Z7duP0TI4Z5foZjl1u1v5XhEIsOPApSIiAxbmXyaL0UspWQvlHF4ihsXsdPPIvS1n4FnaIt5HmwOO2ptD8kPAB5sCeVyOCLDkgKUiIgMW/o0fxixFs/61Xj+9tdBrXc67PJx//hzueKsL5KP+chQPL2zBh3+fIuIc1oDJSIiw5Y+zR8molECN38J9zNbB3U3cVy8e/pX2TJyMuWvxPKyvijgNgTTCFFlHjWQEMk2zUCJiMiwpU/zh4FolPJrL8f9zFYMZF62h+H+cWezZeRkIH8zkvXVfscXby6gfkIgl8MRGZYUoEREZNgKuNO7nNan+UXGWgI3fxmz//WMg1OHkMvLj6ov7Pz/fM1ILqmtIODw5zDgNiypHZHjEYkMPwpQIiIybOnT/NLm2r0D946mQYenwy4fK6tm8HjFpJTb8zEjOaPKS0N1gDJ3/8eVuaGhJkBd1dB1CBQZLhSgAGPMWGPMImPM74wxLxhjgsaYg8aYR40xnzDG6HkSESkh1lq27IuwP5xwvAZKn+YXEWtxvfgsgVu+Bgmnr3BPCY6Gp4+ftrjHXk/5mJE0xrB8TiULasoo95geF3IuoNxtWFBTxvI52kRXJBfURCLpEuBnwGvAw8Be4Fjgg8BtQIMx5hKrDUBERIpe141znc4g6NP8IhKLHWlVvhEG0arcAn+urOX6ky/tXPfUVT5nJL0uw21zK2lqjXLL9kOsbQkTjFnKPIb6CQGW1o6gTpvniuSMAlTS88BC4I/W2s6Pqowx1wF/A/6ZZJi6Pz/DExGRbOi6ce5Aez9B8iI54DY01AT0aX4xSNnnKbNW5XB01umjUz/dY9apQ75nJI0xzBjn484LxuZtDCLDlQIUYK39Sx+3/90Ysxz4NjAPBSgRkaLmdOPcDnOO9/G1GaP0af4QsNaytWNGpTlMKG4JuA311X6W1lZQV+UdMMC6du/A07RxUPs8JYyLP1bV9Vqy10EzkiLDmwLUwDq2nY/ldRQiIjJo6Wyc6wIq/W6FpyHQtawyFLOd69KCccvKPSHWtoRpqE7OAnpdfYcob+N9EM185gng5Qs/xp2BC/C94SXRZSygGUkRSTJa1tM3Y4wHeAKoBeZbax/sfszBgwd7fQJ37dqV49GJiEi63rmxjFDC+UVvwGV5ZFYwhyMSa+Grz/lYf8Dd72sTcFnmjInzrbdF+poY4h3f+TTuaCSzcQBvnjyV3R++Govh2bdc/E+Lh0fb3IQT4HfBO8fEufzEGKdXZN6YQkTyb8qUKb3ePmrUKEdvEJqB6t9NJMPT6t7Ck4iIFJdwmte96R4v6XvmLdeA4QkglDCsP+Dm2bdcfQYY1yDCU2TkGHZfthSMwQCnVyS4cWpm99fneWzy8f5Pi4cNXYLZ+ZVxPjIhxttHJPoMhyJSOBSg+mCM+QzweWAn8NF0v7+vZFuoOmbMim3ccpRew+Kn1zD3ApteJRh3XnlR5nGl9XroNUzftx/eTyThbEPaSMLwh4OVfKBuTO8H+PyQxvonC+ByEX97HdGrb2KKx5Oz17CvMsVQAv6y38PGg15HZYoyMP0eloZCfh21v1EvjDFLgB8DzwIXWGsP5HlIIiKSBdo4t/A82Bx2vBdXAniwpe+wFZs2E+twCscCdsw4gtf/F6Frvg+e3H2m3LX7Y3u3dVWQfFztMcvqvUEWr29DyytECptmoLoxxnwOuBnYDrzbWvt6fkckIiLZsqS2grUtYdodNJLId5vq4SKUxowgkNy7y1pcu3fgbbwXz7bNEImAz0d88ung8YKTUj5fgNCSb5KYdFqGI3fOaffHYBwam0M0tUaZoeYlIgVLM1BdGGOuJRmeniQ586TwJCJSQmZUJcukytz9H6c21UMn4E6vXK3CxAh87xrKvrUEz+PrMJEwBouJhHHvaIJYDOvq//LG+vzE6mYNSXiC9Lo/hmKWW7e/leMRichgKEAdYYy5nmTTiK0kZ55a8zwkERHJMmMMy+dUsqCmjHKP6fEm6ALK3YYFNWVqU+2QtZYt+yJc8fB+jr/rVSrveIXj73qVKx/ez9Z9kQHL0RyXVVrLrDd28ur6q3A/swWTSND91TGAsQmwYF2uHuV81hisL0Csbjbhq67rc5+nbMtmmaKI5J9K+ABjzBXAN4E48AjwmV7eNPdYa+8c4qGJiEiWeV2G2+ZW0tSxaWtLmGDMUuYx1E8IsLR2hPZ+cigb+zc5Kav0JGLcseNnXLLvMdzYHsGpO2MTWLeP+Km1uF94Nrk3lNdPbPpMog2XDdnMU4eMyhRFpGApQCWdfORvN/C5Po5ZB9w5FIMREZHcMsYwY5yPOy8Ym++hFK2ujRF6W9uT2hgBbpvb+4xeR1llX/eDtdyx42d8cN/jeEgjWMSj2BGjOLxijfPvyZGA26TZ/VEznyKFTCV8gLX2BmutGeDPvHyPU0REpFBk0hihN/2WVVrLolf+wmX7NuNngBN1v19rcT2xsSA62qn7o0hpUYASERGRtGWzMUJHWeWq+VUsnBig3GPwJmL8eset/OyFXw5YstcXE42waF0b0UR+Q9SS2goCDmeV1P1RpPApQImIiEjast0YoWtZ5asfOZ43w7/in/c9lmwMkeEYgy5fQeytpO6PIqVFAUpERETSlsvGCK7dO/A0PYoZROixwANjpw9YQjgU1P1RpLQoQImIiEja0t2/KZ3GCN7G+5xthtuPOC5+VH0hUBh7K/VWpmiAco/h4ollPNBQxe3zxvTZrVBECoe68ImIiEja6qv9rNwTclTGl25jBM+2TRmX7UFy9ukvlafzeMUkoHD2VlL3R5HSoBkoERERSVtOGyNEMp99ssDL/rEsrP18yka52ltJRLJFM1AiIiLDiLWWrR2bCDeHCcUtAbehvtrP0toK6qq8jtbgDLh/0xEDNkawFteLz+L77W24n9sGCaetKbrdzZG/146u5eJ3XEPMlXqJo72VRCRbFKBERESGiWgiufltY3OIUMx2lt8F45aVe0KsbQnTUB1g+ZzKAdfidDRGWLyeHvcHyRKXgNvQUBPouzFCLIb/59/Gs2UdJBKDKtsDuOrURdx5wgU9btfeSiKSTSrhExERGQasTYan1XuDtHcLO5BcJ9Qes2m1/R5UYwRr8f9iGZ4t6zCDDE8WWFtZy53Hz+v169pbSUSySTNQIiIiw8DW1iiNzaF+y+2AlLbfM8b5BrzfTBsjdLYqz7Bkr4MFXvaN4eLaL6SseeqgvZVEJNs0AyUiIjIM3Lr9ECGHjRRy3vbbWny/vW1QrcrtkT9PHf8O6s7/EQl3akDS3koikiuagRIRERmErk0ZHtxbRjgBgU2vpt2UIdcebA47ajkOOW77HYvhX7EM97NNaZftdcY/4yI+9Uwil1zFpJPfxu86mmK0hAnGLGUeQ/2EAEtrR1DnYBZNRCQdClAiIiIZ6tmUIRkJMmnKkGuheHptvHPS9tta/CuW4WnakNmaJ2M4fOfDqTeB9lYSkSGlEj4REZEM5KIpQy4F3OlFlqy3/U4k8P72Njyb/4KJhDO7D68/u2MSEcmAApSIiEgGMmnKkE/11X7Hb/pZb/sdClH+hX/B98DdGd+FBWLTZ2ZvTCIiGVKAEhERyUBBNWVwYEltBQGHs0pZbfudSFB+3ccw+1/HQObtyj0+ovMvzc6YREQGQQFKREQkAwXTlMGhGVVeGqoDlLn7Py7bbb89jzR2hqdMWeMiNmM2iUmnZWVMIiKDoQAlIiKSgYJoypAGYwzL51SyoKaMco/pcQGQq7bfvt/dmfH3Wo6Ep7PnEP73r/S6z5OIyFBTFz4REZEMBNyGYBohKutNGTLgdRlum1tJU67bfluLa/cOvI33Ytr2ZTT7ZAE7cjShq2/SzJOIFBQFKBERkQzUV/tZuSfkqIwv600ZBsEYk9u230f2efI0bYRoOPPSPV+A0OduVHgSkYKjEj4REZEM5K0pQyFLJAj88Fo8j/0VEwlhMmzdbn1+YnWzFJ5EpCApQImIiGQgX00ZClY0SuCbn8L9zFaMddpeI5U98idWN5vwVddpzZOIFCQFKBERkQzkqylDQYpGKb/2ctwv7RxUtz2AyMKPEv7k18CjVQYiUpj0r5OIiEiGujdlWLM3SDgBZR5XdpsyFDJrCdz8pcG3Kgfs2PFE/+nj2RqZiEhOKECJiIgMQtemDLt27QJgypQpeR7VEEkk8N5/e7JsL8O76FglZceOp33ZXeBScYyIFDYFKBEREUlfMEj5Fz+MebMt47uwgB0zjsgHP0HsnfOzNzYRkRxSgBIREZH0RCIc85kPQCTzNuVxDI9NnMnUr38bryv/68OstWxtjXLjDh8b2tyEH32FgNtQX+1naW0FdVXe0l7HJiKOaZ5cREREnLOWwLLPDCo8AYRcXq6rms/i9W3YDNudZ0s0YVm0ro2Fa1r5y343oYTBAsG4ZeWeEBetaWXRujaiifyOU0QKgwKUiIiIOObavWPQ3fZiGFZWzeDRYybR2ByiqTWatfGly1rL4vVtrN4bpD1msd0eWQJoj1lW7w0WRNgTkfxTCZ+IiIg45m28b1Dfb4E/V57Ox09bDMYQilm+ufUgo/0u1jaHCcXtkJbObW2N0tgcIhjv/7hgnM6wN6PUOyuKSL8UoEREZEh0rDG5ZfuhvFwoS3Z4tm0aVMe9xyomc+E7vtS5SW4CWPdaBNeR/4ajpXNrW8I0VAdYPqcyZ+ukbt1+iFDM2axSKGa5dftb3HHBmJyMRUSKg0r4REQk57quMVm1J0QwbrXGpFhFIhl9mwVe9o9l3vTrO8NTV4le/n8oSucebA73OHdfEsCDLaGcjENEiocClIiI5FT3NSb5ulCWLPGlV75mj/xZW3kGp537Q2Jub1rf37V0LhdC8fR+3oIOZ6tEpHQpQImISE5lssZEClds2kzSjRDfqvkAF077EjFXZisHOkrnciHgTq80sMyjMlOR4U4BSkREciqTNSZSuKINl4LH2SyUBf7uHcU3Jv7zoM6Zy9K5+mq/44shF1A/IZCTcYhI8VCAEhGRnNIak9KSmDSV2FnnYz39l+JZ4LDLx5RzfgCuwV9u5Kp0bkltBQGHs0oBt2FJ7YicjENEiocClIiIZMxay5Z9Ea54eD/H3/UqlXe8wvF3vcqVD+9n674I1lqtMSk1xhC+6jpiZ70T6/H1KOfrWPMUmzSVf//o/+AqK+9xsZHJxUeuSudmVHlpqA5Q5h7g/G5oqAlQV5XeGi4RKT1qYy4iIhmJJpLNIRqbQ4S6NIfo3oI64DYE0whRWmNSBDwewouvJ7p7J97Ge/Fs2wzRMHj9xKbPJNpwGYlJp/Fza2nqaF3fEiYYs5R5DPUTAhwIx3nktYij2clclc51tNaPJiyRzoFY6NKo3UVy5qmhJtlOvdBb7Wu7AJHcU4ASEZG0de2s11tziK6d9cb4XbzW3rP7Xm+0xiRPEgk8jzTi+90dmLbWzptt5TgiH/w3YufX9yzDM4bE5KmEl9xAuI+7NcYwY5yPOy8Y2+NrW/ZFWLimlXYHM465KJ3r6wOAZHhKhii3gffXBPjsGRXUFcHmuU4/1Mjlvloiw4FK+EREJG3pdNZrDSfwDVAe1UFrTPIgFKL885fh/+X3MG2tGDj6p20f/tu/Q/nnL4NQdtemOS2dA2iPWy5cvS+lNHQwBmqt3zED5XOBx2U4swjK9rRdgMjQUYASEZG0pdNZLxKHqoBba0wKUSJB+Zc/ijmwrzM0ddUZpA7so/zLH4WE03YgAzPGsHxOJQtqyij3mAEvSEIJsrbpcim21i/FxyRSqBSgREQkbel01rPA/lC8zwtlF1DuNiyoKSuKNSbddW2k8c6NZZzzaFmPRhqFyrN+dWd46k9HiPI8siar5/e6DLfNrWTV/CoWTgxQPsD6t2zNopRia/1SfEwihUoBSkRE0pZuZ71QnB4XygYo9xgunljGAw1V3D5vTNGty4gmLIvWtbFwTSur9oQIJQwW07nmJBuzJbnk++1taR7/i6yPoes6qZXzqwYMUTD4WZRSbK1fio9JpFCpiYSIiKQtk856/TUUKEbpNNJYvD4ZIAttds0cemPA2afOYwHefCMn4+joHPexv+x31FQCjs6i3HHBmLTPV4qt9UvxMYkUKs1AiYhI2uqr/Y7fQEq1s17RrjmxFteLz+K/9ev5HgmQOov3arvzNVaDmUUJuNMLssXQWr8UH5NIoVKAEhGRtC2prSDg8AKsVDvrFeWak1gM//L/pOym/4dny3rHs0+50r1zXLoynUUpxQ8ASvExiRQqBSgREUmb0xbUpdxZr+jWnCQSBH74JTyP/RUTCWEKoLmF01m8vmQ6i1KKHwAsqa3Ar+0CRIaEApSIiKRtoBbUxd5Zz4miWnMSixH4/jW4n9mCsZm1IrdA/OTTsjqsdGbxuhvMLEopfgDwjjEexgYGTlABV/E8JpFCpQAlIiIZ6a0FdSl01nOqaNacWIt/xTLczz4x6JK9yEeWZmVIHdKZxetuMLMoxhh+9s7RnHesn95eRoOl3E3RfABgreWTj7zB/tDAU3lVZW5+9s7RBf+YRAqZuvCJiEjGSq2zXjrqq/2s3BNyFADysubkSLMI33//GPee5wcVniwQO3seiclvz9bogPRn8ToMdmYomkgGjs3/CNOzw7zFAOce62f5nMqi+AAgnVLIA+EETx+IMWOcL/cDEylRmoESERHJQEGvo4nF8P/XNyj71tJBhScLRF0eIue+i/Dir0KWZy3SncWDwZeGdm8/3zPCGRIYNv8jPKjNeodSUTY0ESliClAiIiIZKNh1NIkEgR9ci+dvf8XYxKDC033jzuU9M77Gx96+BOt22KEgDel0jgM4sdw96NLQom0/34+ia2giUuQUoERERDJQkI00OppFPLt1UCV7cQz3jTuXy0//DBuOmZyzIJHOLF6523DXu8ZQN8jSs1KcrSmqhiYiJUABSkREJEPdG2kEXDbZgCAfjTSy2Cwi5PLyo+oLj/5/joJEPmbxSnG2pmgamoiUCDWREBERGYSujTR27doFwJQpU4Z8HK7dO/A0bcy4TXmHGIaVVTN4vGJS5225ChIds3iL1yfL5UIxmxJuXCTDQUNNIGuzeKU4W1PwDU1ESoxmoERERIqYTSR4bss22r//VYgMLuRY4C+VtXz8tMU9GkbkKkgMdTv8UpytKeiGJiIlSDNQIiIiRSoaifLMsm8w4+XHCCSigyrdi2H4c2UtF77j2l677eUySAxlO/xSnK3pKIXs6CzYl2LaGFikkGkGSkREpAjZRIJnv30D5720kbJBhCcLxHDxm/HncfEZX+g1PBVLkHCiFGdrCrKhiUgJ0wyUiIhIEXq+6WnO3bMJj+OWCL2LY/j3t13FXcfP6fOYYgkSTpTqbE1HKWRTa5Rbth9ibUuYYMxS5jHUTwiwtHbEoDsYikiSApSIiEgR8v72NtyDDU/Gxf+NO4e7jntnn8cMJkhYa9nacUHfHCYUtwTchvpqP0trK6ir8g75bMhAjSsMFr8LFtSUF91szVCWQooMZwpQIiIixcRaXLt3UPva04Mq28O4iJ01h9+//VOUvxLLege8aMKyeH1bj5ASjFtW7gmxtiVMQ3XyvoekzXsX/c3WzB4d4/ITY3ygbsKQjklEiocClIiISLGIxfCvWIanaWPGd2EBO7KS0OeWkZg8lZ9bm/WyL2uT4amvMrkE0B6zrN4bZPF6uG3u0M/09DVb09GKXkSkLwpQIiKStkIszSp5RzbK9TRtwETC6X87gC9ArG4W4auuA0/yEiAXZV9bW6M0Nof6XWMEEIwny+iaWqPM0PocESkSClAiIpKWQi7NKmWdG+VmEJ4gGaBCX/4RiUmnZXdgvbh1+yFCDveNCsUst25/izsuGJPjUYmIZIfamIuIiGNdS7Pau62Zge6lWW1Ym5vNV4cVa3G9+CyBW76e8Ua5FnhjyplDEp4AHmwOO25vkQAebBncBsAiIkNJM1AiIuKYSrOGWNc1T5FQxk0j4hh8l30iq0PrTyieXnAOOpytEhEpBJqBEhERxzIpzZIMpax5yjw8RTH8fvw5bB01JavD60/And5oyxxubCsiUggUoERExDGVZg2dbKx5ChoPvx1/Hh877VPc+szh7A6wH/XVfscXGC6gfkIgl8MREckqlfCJiIhjKs0aOt7G+yCaeXhq9o3h0trPsWXkZCC9MDvYLotLaitY2xKm3cHrH3AbltSOcDw2EZF80wyUiIg4ptKsoePZtgmTYROOdpePy7qEJ3AeZqMJy6J1bSxc08qqPSGCcZuczTrSZfGiNa0sWtdGNNH3/c2o8tJQHaDM3f+5ytzQUBOgrsrraGwiIoVAAUpERBxTadYQikQy+rbDLh8rq2bweMWklNudhNlsdVk0xrB8TiULasoo95gePzMuoNxtWFBTxvI5Q7+JrojIYChAiYiIY0tqKwg4nFVSadYg+dLrXmg5Gp4+ftpi6BJKnIbZTLos9sXrMtw2t5JV86tYODFAucdggHKP4eKJZTzQUMXt88ZorzARKTpaAyUiIo51lGat3hvs9yJbpVmDF5s2E8+WdY7K+Hpb89SV0zCb7Q1wjTHMGOfjzgvGOrpPEZFioAAlIiKOdZRmLV6fnIEIdSvzcpG8WG+oCZR0aVZfTRZmj/bxkQkxTrF20I892nApnm2bHW2e29uapw7phFl1WSwMg23iISK5pQAlIiJp6SjNauq4wGsJE4xZyjyG+gkBltaOoK6EN8+NJpLrhLoHyGDc8uf9bh5tc7PqzTaWz6kcVHlaYtJUYnWzjuwD1Xc3vrDbx+qqGWzttuYpkzCrLov519/P18o9Ida2hGmoDgz486UQJpI7ClAiIpK24Vqa1bXJQm8ljBZDKMGRJgtw29xBzMIZQ/iq62DFMjxNGyEaTinns8aA14/rzFkc+8EvsHBHcNBhNuA2BNMIUeqymF0D/XylNvHo++crWyFMRHqnACUiIuJQJk0WZgxmNs7jIbz4eqK7d+JtvDdZ0hcNg9dPbPpMog2XkZh0GnXAnceVZ36eI+qr/azcE3JUxqcui9mXjZ+vbIUwEembApSIiIhD2W6y4IgxJCZPJbzkBjLbVtc5bYCbX9n4+RrykC8yDKmNuYiIiEOl3mRBG+DmVzZ+vjIJYSKSHgUoERERh0q9yYI2wM2vbPx8lXrIFykEKuETERHpor/uZT4XhJ1enVKcTRaGe5fFfMpGE49SD/kihUABSkRE5IiBupelM+FSzE0WstVlUa2005ONJh7qpCiSewpQIiIiOOteRhof1gdc8OWKZvy3/vjIhrgR8PmITTvaPS+tRFZk1Eo7fdlo4qFOiiK5pwAlIiKaKcB59zInKojy8I6bmfbw02ATdD5zkTCex/+K58lNxGbMTu7z5Cm9t2K10s5MRxOPvp63Dv018VAnRZHcUxMJEZFhLpqwLFrXxsI1razaEyIYt1iOzhRctKaVRevaiCZKe61EOt3LANym55uowTLCRnju8S8w7bVtmK7hqfMYMNEwnq2P4l+xDGzpPa+ZtNKW7DTxUCdFkdxTgBIRGca6zhS0dymz6pA6U9CGLcGL/Q7pdC8D8Lpg4cQA5R6DAco9hveMibJr5zcZd7i1R3DqzkQjeLZuwLV75yBGXZjUSjtzHU08Vs2v6vHzdfHEMh5oqOL2eWP6LHtUJ0WR3Cu9ugEREXEs15tuFlNpYLrdy8JxUpssxGLEv7WUkf94acDw1CkaxrvmPsKf/npa5y50Q91Ku5h+zpwYbBMPdVIUyS0FKBGRYSyTmYI7Lhjj6PhiayIwqO5l1uJfsQz3np3OwxPJcj7Pk5sIp/E9xWAoW2kX28/ZUMlWJ0UR6UklfCIiw1iuZgqKsTSwvtrv+E2xe/cy1+4deJo24srkcURLLT4lw2g6Mm2lXYw/ZyJS/BSgRESGsVzNFBRjE4EltRUEHF7Id3YvsxbXi88SuOXrEMmwDM3rz+z7Cthgwmg6ivHnTESKnwKUiMgwlquZgmJsIpB297JREPjeNZR9awmmbV9apXsdLBCbPjOT4Ra0jMJoBorx50xEip8ClIjIMJarmYKhbiKQDWl1L5s5gmO+9BHcz2zBJHq2Knd+UhfR+ZcObuAFaKhaaRfjz5mIFL+8ByhjTIUxpibf4xARGY5yNVMwlE0EsslRC+m5lVT8+MuY/a9nHpxIzj7FT68jMem0LI2+cAxVK+1i/TkTkeJWCF34PgN8ExjgcyoREcm2jpmC1XuD/a4jSXemYFAd7boZ6hbVA3Uvc734LO5nnxh0eLJjxxP63I1QRO210zEUrbSz+XMmIuJUIQQoERHJk46ZgsXr6dEGGpIzBQG3oaEmkNZMQX21n5V7Qo7Kq/orDSy4FtXW4vvNbWDT2XK3210A8UlTCX3lFvCU9ttwrltpZ+vnTEQkHTn5l9sY88U0Dn9nLsYgIiLO5GKmYEltBWtbwrQ7KJnqqzSwa4vq3mbHUltUw21z0ysFS3tmKxZL7vW0oynj2SdrXMRPryP0he+V7MzTUMrGz5mISLpy9dHXTSQ/ZHP67qCiZBGRPMr2TEE2SgMzaVE9w2HQS3tm68hGuZ6mDRl328MXIFY3i/BV1yk8ZUmuSlBFRPqTqyYSrcBDwMkO/vwwR2MQEZE8yUYTgVy1qM5k89WOjXJNJLnp7d72CP/y+B5mr3+e37/6Rv/nA+yYcQS//CPCn/xayZftDaWhalYhItJVrv4Vfxyotda+PNCBxpi2HI3BMWPMh4C5wHRgGlAB3G2t/Ug+xyUiUswGWxqYaYvqgUrzEtY6m9mKWfZt305w60OMe3pD57qnl9sjvHvDC+xpjwBw2ZY97Gs4g5HePnohuVyElnyzJLvtFYKhaFYhItJVrgLUFqDBGHOctfbvAxz7BrA3R+Nw6qskg9NbQAugdzkRkSwYTGlgJi2qnZTmVfrNgO2sPYkYd+xczsLWJgKJCOZIpXlLMMJ7u4QngISFN2PxXgOUBeJTS7NVeSHJdbMKEZGuclXC912S5Xn7BzrQWvtTa+3JORqHU1cDpwIjgU/meSwiIkJy0X86yjzGUWneK4cT/S+8TSRY+fT3uOT1zRyTCOM+cvSrwSjv2fACu7uEJ4CrJo5lQlnPGY7OVuVXl26rchGR4SgnAcpae9ha+7K1NpqL+882a+3D1tpd1lo1sxARKRD11X7Hb1Iu4OxxXv44QDOBgXgSMVY//R3e27YdT5eY9fdQMjy9cDg1PF1ywmh+csaEHvdjgXjtWbR/93/Bq8YFIiKlZNAByhjzHWOM43cHY8zJxpj1gz2viIiUtiW1FQQcbnzqd8OzB6KEBhGesJY7di7nXW3PpHTaez0c5b0bXuD5w+GUwz94/CjumnESnm77T1kgfvoMQtd8Xw0jRERKkBnspIsxJgFsAz5irX1mgGM/QbLr3ghrbR+rbfPLGDMPeBiHTSQOHjzY6xO4a9eu7A5MRGSYsRa++pyP9QfchBJ9B6mAyzLaY/l7xOB894yeznnzBdY+uYwRiaNBqTUc4z0bXmD7oVDKsQuPG8mvz56Iz9Xzc8i418cLH/kC7SfmuzpdRER6M2XKlF5vHzVqlKM3kWyV8L0D2GKM+X+9fdEYM84Y8wfgFyQ73ImIiPTLGPjGqRHmjIkTcNnORg6dX8cScFmmjYzTFh1ceAK4unk1ZYmjJXoHIjHqN/YMTwuOHck9Z/URnjw+Dp46jfYTJg5qLCIiUriyWVvgB75njHk/cIW1thnAGPMB4OdAFYN9dysifSXbQtUxY1Zs45aj9BoWP72Gvbv3VNtHi+oyltaO4CfbD/H4G6GB72gAC/Y/0dkwoi0SY/7GF9n2Zur91o+v4L6zJ+J3p4Ynawx4/RyccgYvL/w3ppx66qDHI/mh38Pip9ewNBTy65iNAPUCcArJsm8DzAOeMsZ8ATgf+NiR222XY+7LwnlFRGQYGKhFdTr7RXVwAceXu2mLJDrbnXfMPh2MxmnYtJumg8GU73n3uBH89pyTCXQJTxbA5SJ21lyiDZfxcrwgq9NFRCSLshGg3gH8J/A5ku9JBhhFslwPjoYnA+wDPmWtvT8L5xUREUl7vyhItkj/1bzRjGx+jtDKXzO9eSsGeDMaZ8GmF9nyRnvK8fOqRvC7cyZR1m3mCV+A4JduJjF5avL/tf5VRKTkDTpAWWtDwDXGmN8AtwOnczQwdfUb4NPW2tbBnlNERKRDwG0IphOirOWz/j3MufGzmLZ9QPIN661YnPdv3s1jbanh6Z1jj+EP555Muad72Z6LWN0sbZIrIjLMZG0NlLX2b8aY60mW57mhc7WvAfYC1yk8iYhIttVX+1m5J+SojM+TiHHv88tZ+I/HwCY6P+k7HItz0ebdbDxwOOX4WWOOYeW5kzjGk1qal2xVXkf4quu0Sa5knbWWrR3r/prDhOKWgNtQX+1naW0FdVVejH7uRPImKwHKGDMW+AnwL7182QI1wLYjAetHhbZh7ZFGFx848r/HHfl7pjHmziP/3Wqt/cIQD0tERBxYUlvB2pYw7bEB3lqs5a4dP+P9+x/H2KNxqz2W4OLHXuKR/anh6dzKch44bxIV3m7hybiIn15H6AvfSwlP1lq2H3LxPy0eNm56VRe9WTAcg0Q0YVm8vo3G5lDn+jyAYNyyck+ItS1hGqoDLJ9TiddVWo9dpFgMOkAZYz5Ccm+nsaSud3qZZBjxH7mtHPg+cIkx5uPW2ucGe+4smg5c0e22SUf+QPKxKECJiPQhnxe6M6q8NFQHWL03SLCfjXTnvfkcH9z3WGenPYBgPME//W03f219K+XYs0aXsXrmZEZ2CU8dnfZidbN6zDx1XPT+8WU/4QTYI+cY6oveUgocwzFIWJt8zH39LCeA9phl9d4gi9fDbXMri+b1FCkl2dgH6i6S4QmOhqefk1wLdRbJTXa7BqvzgCeycN6ssdbeYK01/fyZmO8xiogUqmjCsmhdGwvXtLJqT4hgPBkfOi50L1rTyqJ1bUQTuSk+MMawfE4lC2rKKPeYHm9sLmAkMR546jsp4SkUT/Chv73En/elhqczR5XROHMyo46EJwtYX4DY2fMIfvlHhD/5NfAc/fyx60VvKGGw3ZYAp170tpGrIox8vw7Z1PU5be8SnjoM1XM61La2RmlsDvX7QQBAMA6NzSGaWqNDMzARSZGtjXQhGY7+AVxorf2ktbbdWvsMcDZwE8l/7zr+hfNn8bwiIpInhXKh63UZbptbyar5VSycGKDcYzBAucdwcY2fPS8uwx+PdEabcDzBpY/v4cHXD6Xcz7SRAdbMmkylr0uBhjEcXrGG8Ke/3mvDiEK46C2U1yFbCuE5zYdbtx8iNFAp6hGhmOXW7W8NfKCIZF22ApQB7gdqrbWNXb9grY1Za68D5gK7s3Q+EREpAIV0odt1v6hXP3oCbR8/kVf/dTz3bL2RiubnOsNTNGH51y17WP2PN1O+v7YiwJpZpzDW16263dv/Z36FcNFbSK9DNhTCc5oP6explgAebBn8BtIikr5sBKg3gY9Zay+x1h7o6yBr7UaSa41uy8I5RUSkABTsha61uF54hvJr/hX3M1tTwtPlW/aw8u+p4WnqCD9rZ09mnD81PFljiE2f2e+pCuGit2BfhwwVwnOaD+nuaRZ0+JqLSHZlowvfGdbaZicHWmsPA/9ujPl9Fs4rIiJ5VpAXurEY/hXL8Gx9FKJHy/ZiCcvHtr7M/712MOXwt43w89DsUxjv9/a8L6+f6PxL+z1dIVz0FuTrMAiF8JzmQ7p7mpV51EBCJB8GPQPlNDx1+57Vgz2viIjkX8Fd6FqbDE9NGzBdwlPcWj7+xF5+8+obKYefcoyPh2adwnGBnuHJeryONsoNuNO7iM3FRW/BvQ6DVAjPaT7UV/sdX5i5gPoJgVwOR0T6kM0mEiIiMswU2oWua/cOPE0bMZFw521xa1n0xF7uaWlLOXZSuY8/zT6FE8pSw5M98ic2452ONsothIveQnsdBqsQntN8WFJbQcDhaxNwG5bUjsjxiESkNwpQIiKSsUK70PU23gfRo+EpYS2Ln2zmv5tTw9NJZV4emn0KE8p8vd5P+OPXEP5UarvyvhTCRW+hvQ6DVQjPaT507GlW5u7/uDI3NNQEqKvqpexURHJOAUpERDJWaBe6nm2bMEdadCes5dPbWrhjb2p/o+oyL3+afQonlfcMTxaI155FbO4Cx+cshIveQnsdBqsQntN8cLKnWbnbsKCmjOVztImuSL4oQImISMYK7kI3EgGS+yJ99qlXWPHy/pQvnxhIhqeTj+nZmtwCdux4Qp+7ccCyva66XvQGXBZD6vqiobjoLbjXYZCGc5Dod0+ziWU80FDF7fPG4HWVzmMWKTbZ6MInIiLDVMeF7uL1yf2FQt02cXWRnPFoqAkMzYWuz4cNh/j89lf42Z7WlC8d7/fw0OzJTO4jPMUnTSX0lVscle1113HR+4cn2vifFg8bD3oJxixlHkP9hABLa0dQN673csFsKLjXIQs6ntOm1ii3bD/E2pbwkD6n+dR1TzMRKTwKUCKSU9ZatnZcADWHCcUtAbehvtrP0toK6qq8RXExJ30rpAvd6DvO4yu/uoef7E4NT8f6PTw0+xROHdFz7Y81LuKn1xH6wvfSmnnqzhjD6RUJbpwaYcqUkzK+n0wV0uuQLQoSIlKIFKBEJGeiCcvi9W09PhEPxi0r94RY2xKmoTr5ibjKUYpbIVzoWmv5ygv7ufmF11NuH+fzsHbWZE6rSA1PFsDjI3bW+Y667RWDQngdRERKnQKUiOSEtcnwtHpvkGC859cTQHvMsnpvkMXr4ba5xVFWJAXAWly7d+BtvBfPts3JdU8+H9fvs9z817+lHDrW52bt7MmcPrIs9S4AO2Y8oSXfIDF56hAOXkREip0ClIjkxNbWKI3NoV7DU1fBeHLNRlNrlBlFVl4keRCLHdkodyNEw50d9/7z6Ze5aeffUw6t9LpZM+sUzugSnqwx4PUTq5uVnHXKYL2TiIgMb3rnEJGcuHX7IUIxO/CBQChmuXX7W9xxwZgcj0qKWiJB4Idfwv1sE8YebZFw4/P/4BvdwtOogJ/fr/gZtS8/hd22Obk3lNdPbPpMog2XkZh02lCPXkRESoQClIjkxIPN4ZQuYP1JAA+2hHI5HCl2sRiBH16L+5mtdC30/N6uf3D9jtdSDh3pcbF69hTOePtUwhd9kDAiIiLZo32gRCQnQnFns08dgg5nq2QYshb/imW4n30iJTz96MXX+fKzqeFphNvFH2dO5pwKL9419w3tOEVEZFhQgBKRnAi402sIUeZRAwnpnWv3DjxNG1PK9m7dvY8vbH815bhj3C4emDmJmWOOwViL58lNQz1UEREZBlTCJyI5UV/tZ+WekKMyPhdQP6Hn/jzSt+G0v5a38b7kGqYjlr/UyueefiXlmDK3YeV5kzh/7IijN0ZVvCciItmnACUiObGktoK1LWHaHZTmBdyGJbUjBjxOkobF/lpdW5U/vq6zdO+2PftZ8lRLyqEBl+EP505iblW3nyGvP0dD6zu8fmCUi7ePcLr6T0REipEClIjkxIwqLw3VgT73gepQ5oaGmgB1Vd6hG1wRGxb7a3VvVX7k5jtf3s/ibc0ph/pdhv8792TeNa4i5XZrDLHpM7M+tIHC6xqXnzlj4tx9ii3e8CoiIv3SGigRyQljDMvnVLKgpoxyj+nxj40LKHcbFtSUsXxOEV7k50km+2vlg7WWLfsiXPHwfo6/61Uq73iF4+96lSsf3s/WfRGs7WNmMpEg8MNr8Tz2V0wk1LnP0/80H+CqJ1PDk89l+O05J/O+8SN73o/XT3T+pVl/TB3htb1LeOocOhBKGNYfcLN4fVvfj1FERIqaApSI5IzXZbhtbiWr5lexcGKAco/BAOUew8UTy3igoYrb543RJ/VpyGR/raEWTVgWrWtj4ZpWVu0JEYxbLEdnaS5a08qidW1EE90eRyxG4PtfTLYq79Iw4p6WNv6taS9dj/Yaw71nT6Th2J7hyRoXsbpZWd/ryWl4DSVMXsOriIjklkr4RCSnjDHMGOfjzgvG5nsoJaHQ99fKuMTwyMxT932efvNKG1dsfTnlMXsM/Prsk7jouFE9zw/E315H+KrrIMuzmqW2OfRwakQiIpJNClAiIkWk0PfXyqTEcMZoQ+DbS3G/tDMlPP3u1Tf4SLfw5DZw91kTufj40T3u0xoX8dPrCH3he1kPT1D44TUdw6IRiYhIjqiET0SkiBT6/lrpztL8bFsb5dde3iM8rXrtIP+6ZQ9d86IL+O8ZJ/HPJ4xOuR9rDNYXIHbuPEJX35ST8ASFH16dcrKW6+gsodZyiYh0pxkoEZEiUuj7a6U1S2MtV/3xPzEHXk8JT3/8+0EufXwPXfOHAe6sq+HSEytT7sMCsbPnEW24LOtrnroLuA3BNEJUoW4OndEs4Tjf0AxORKQIaAZKRKSILKmtIODwwjwf+2ulM0tzzqEXmXPg2ZTw9ODrb3LJ43uIdpn1MMDtZ9bw4erU9UTWGGLnXED401/PeXiCZHh1+qZZyJtDF0MjEhGRQqYAJSJSRDr21ypz939cvvbXclRiaC3nvPkC927/EZ4u81V/ev0Q//zYS0S6def7+fRqPlbTSzOGHLQq70+hh1enSmktl4hIPihAiYgUkULfX2ugWRpPIsZ/7/gpDz25jAmRts7Zp7+2HuKf/rabULfw9F/TJvBvJ/Xs4JirVuX9cRpeAy5b0JtDl8paLhGRfFGAEhEpMoW8v1a/szTWcsfO5Sxs3coxiXBneHpk/1ss3PxSj/VFPznjRP59YlXPuwHip+emVXl/nITXgMsyZ0y8oDeHLvRGJCIihU5NJEREilCh7q/VMUvTYx8oa7ny7+v40OuP4e1SQLbxwGHev2k37fHUorIf1p7IpyaN63H/FojXnpWzVuUD6QivTR37J7WECcYsZR5D/YQAF49q4/SKREG3/i70RiQiIoVOAUpERLKmY5Zm8Xo69xhyJWLcsXM5l7z+GO4ul+2PHTjMhZte5HC38PSd00/gM5P7CU85bFXuRH/hddeu/XkYUXqW1FawtiVMu4PSvEJeyyUiki8q4RMRkaxKKTE8yc9/P5cs2/OQ6Czb29LWTsOmFzkUSw1P3556PJ8/ZXyP+7RAfNLU5MyTR5/9DUahNyIRESl0ClAiIpJ1HbM0d530Oh9qe4JjEpHOrzW90c78TS/yZrfwdMNpx3Htqcf2uC8L2LHjCV33k7zOPJWKQm9EIiJS6PQxnoiI5Iy38T6Ihjv/f9vBIPM3vsgb0dRdXL/6tmP56tuOS7nNAhgX8dPrkmV7mnnKmoHWci2tHUGdNs8VEemV3o1ERCRnPNs2YY5sirv9zSD1G1/gQLfwdO2U8Xy9l/Bkx4wjtPQ/h7RV+XBSqI1IREQKnQKUiIhkj7W4du/A23gvnm2bIZKcfdpxKMT7Nr5IayQ1PH3+lPF8a+rxPcvEfAFCS76p8CQiIgVHAUpERLIjFsO/Yhmepo0QDXfOPD13KMR7N7zA6+FYyuGfmTSOm97eMzxZn3/IN8kVERFxSgFKRCTPrLVsbY1y4w4fG9rchB99hYDbUF/tZ2ltBXVV3sJfyG/tkfC0ARM5uubphbfCvHfji/y9W3j69MlV/KD2hJTH1bHmKVY3e8g3yRUREXFKAUpEJI+iCcvi9W00NocIxtzYI42+g3HLyj0h1raEaagOsHxOZUFvzuravQNP08aU8LT7cJj3bHiBV0PRlGP/Y+JYfnTGiT1DoctF+N+uIfbOhqEYsoiISEbUxlxEJE+sTYan1XuDtMdsZ3jqkADaY5bVe4MsXt+GtQNvfJov3bvt7WlPhqeWbuHpEyeN4ZZ3TOi9bO+cecTOnz8k4xUREcmUApSISJ5sbY0mZ57i/R8XjENjc4im1mj/B+ZR1257zcEI793wInuDqeO9onoMP5tWjatr2Z4xWF9AZXsiIlI0VMInIpInt24/RCjmbFYpFLPcuv0t7rhgTEbn6lhndcv2Q6xtDhOK28Gvs+rSca+j294rwQjv2fACL7VHUg69fEIlvzjzaHiyAL4AsekziTZcpoYRIiJSNBSgRETy5MHmMAmHxyaAB1tCGZ2n6zqrUMx2nnNQ66y6d9wDXgtFec+GF3nxcGp4uuzE0fyyrgZ314DmC3B4xZqMHo+IiEg+qYRPRCRPQvH01jQFHc5WddV9nVX3wJbROquUjnshjLX8IxTlvRteYNfhcMqhHzphNL+qOyklPFljiE2fmfZjERERKQSagRIRyZOA2xBMI0SVedJfH5TJOqsZ43y9H3SkZM/3mxW4dzzR2fJiXzjGeze+yM63UsPTB44fxX/POAlP91ktr5/o/EvTfiwiIiKFQDNQIiJ5Ul/td/yPsAuonxBI+xyZrLPqVSyGf/l/UnbT/0sJT63hGO/b+ALPHkotL3z/cSP537NO6lESqE1yRUSk2GkGSkQkT5bUVrC2JUy7g4ATcBuW1I5I+xxZWWeVSBD44bW4n30CY4/e24FIjPmbXuTpN1O/p+HYkdx71kR8rqPx0BoD3mR46ui2l5PGFiIiIjmmACUikiczqrw0VAdYvTfYb4ldmRsaagLUVXnTPseg11nFYgR++CXcz2xN2aXqjWiMho0v8uTBYMrh7x1XwW/Onojf3SU8AbGz56V028tJY4sh1Ff4mz3ax0cmxDjFWoU/EZESpRI+EZE8McawfE4lC2rKKPcYDKnhxQWUuw0LaspYPqcyowvygDu970lZZ3WkWYT72aaU8HQwGqdh4262dgtPF1SN4P/OPZmAu9tbiy9A+NNf7wxPOWlsMYSiCcuidW0sXNPKqj0hgnGLJRn+/rzfzeKn/Sxa10Y0UVjjFhGR7FCAEhHJI6/LcNvcSlbNr+LdY+MEXBYDlHsMF08s44GGKm6fNybjWZjBrLNy7d6Bp2ljStneoWicCze9yONvtKd875yxx/CHcydR1i089dZxr5g3EB4o/FkMoYQp2PAnIiKDpxI+EZE8M8YwY5yPG6cm90+aMmVK1u57MOusvI33QfRoZ723YnHev3k3m9tSw9PsMcew8rxJlHt6iWq9dNwbyg2Esy2rXQ1FRKQoaQZKRKSEdayzKnP3f1zKOitrcb34LJ6t6zFHZlAOx+Is3LybDQcOp3zfeZXlPHDeJEZ4ep6gr457Q7WB8ECstWzZF+GKh/dz/F2vUnnHKxx/16tc+fB+tu6L9Dp7lLWuhiIiUrQ0AyUiUsI61lktXk+Phg2Q/BQt4DY01CQbNph4/MgmuRshkTwyGE/wgcdeYv3+1PB09uhy/jhzMhXe1PBkAXyBlI57XQ3FBsIDybSJRaGEPxERyR8FKBGREtexzqqpo2tcS5hgzFLmMdRPCLC0dgR143ydTSM8TRswkWTpXiie4IOPvcTDrakzKTNGldE4axKjeglP8bfPIHLJVT1mnjo617kMpJOhMtlAuD9d1zH1VoqX2sQCbpt7tIFHIYQ/ERHJLwUoEZFhoGOd1Z0XjO3zmM6mEUfCUzie4EN/e4mH9h1KOW76qDIaZ01mtDf1LcQaF7Fz5xFefH2PWaeuMz7pZJBMNxDuz2DWMQXchmAaDyDb4U9ERPJPa6BERARIbRoRSSS4bMse1ryeGp7OGBngwVmTGePrFp6A+Ol1vZbsde9cl45MNxDuz2DWMQ2mq6GIiJQGBSgREQHAs20TxlqiCcuHt7zMA39/M+Xrp1cEWDvrFMZ2D0/GRbz2LEJf+B54ehY2OJ3x6W4wGwj3ZzDrmJbUVhBwOKuUi/AnIiL5pxI+EZFhyFrL1n0R1vy1iXMe/z31rU9AIkIsYfno1pf5/WsHU44/bYSftbMmM85/9G1joGYRHdKZ8elQ3rWxRQYbCPdnMOuYOroa9rV+qkOuwp+IiOSfApSIyDATTVg+/dd9/NOam7l+XxOBRAQ3lljCckXTy/z21TdSjj/1GD8PzT6FYwPdwoBxEfzyj3o0i+gunRkfALeBBxqqko0tcmAw65gG6mposPhdsKCmPCfhT0RE8k8lfCIiw4i1lsXrDvCBNTezYN9WjkmEcWOJW8snntjLva+8kXL85GN8PDR7Msd3C0/WGGJnzx0wPEH6Mz4JS87CEwx+HVNHV8NV86tYODFAucdggHKP4T1VcZafEeb2eWNS2p+LiEjp0AyUiMgwsrU1yr5nnmHBviaOSUQASFjLvz/RzN0tbSnHTiz38dCsUzixrJcw4/UTnX+po3MWWue6JbUVrG0JO2po0dc6pr66Gu7atStr4xQRkcKkGSgRkWHCJhLc9+AWfrXtZo5JJLvtJazlk9ua+VXzgZRja8q8/Gn2ZGrKe4Yn6/MTq5vlaPYJCq9zXcc6pjJ3/8dpHZOIiPRGAUpEZBiIhiO8eN3V/OShrzIh0oYhWc639KkWbn85NTxNCHj50+xTmFjuT7ndGoP1BYjVze63aUR3hda5rmMd04KaMso9pscboYtkE4sFNWVaxyQiIj2ohE9EpMTZSITo1R9m2uFWOqKAtZbPPf0KP9+zP+XYEwJeHpp9CpOOORqeOrvtTZ9JtOEyxzNPHQqxc13HOqam1ii3bD/E2pYwwZilzGOonxBgae0IR+uwrLVs7biP5jCheBl+F8xv2c/S2grqqrwKYCIiJUYBSkSklFlL9LtfYly38PT57a/y05daUw49zu/hoVmTmTIideYJX4DDK9ZkPISBOte5SM485apteX/j6m0dk1PRRHKD4NTHZAglYOWeEGtbwjRUJx+TGkqIiJQOlfCJiJQw1+4djHrhyZTw9KVnX+Mnu/elHDfe72HtrFN4W0Xq+qM4htj0mYMeR3+d6y6eWMYDDVVF1bkukUhw6dr9/O6lIO3dAiEkN+Btj1lW7w2yeH0b1qbXiVBERAqXZqBEREqYt/E+XDZ5eW+t5as7XuMHL7yeckyVz83aWZN5+8hemjd4fY677Q1ksDM+hSKasFz60AEefi084LHBeHLWrak1yowctmYXEZGhoxkoEZES5tm2qXP26RvP/Z3v7EoNT2O8bh6cdQq1I8t6fG/Y7SNRNzvtNU+lzNpk2d56B+GpQyhmuXX7WzkclYiIDCXNQImIlLJIcq+nbz33d7713D9SvjTa62bNrMlMG5UaniwQdvsxZ80m/O/Ou+0NB1tbozQ2h0hnb+AE8GBLKGdjEhGRoaUAJSJSynw+vrN9Lzfs/HvKzaM8LtbMnEzd6PKU2y3wevlYKq75tmaeenHr9kOEHGzA210wg+8REZHCpBI+EZES9v1DPr6y47WU2yo8LlbPnMxZleU9jo9jeGXR1xWe+vBgc7hHwwgnyhzugyUiIoVPM1AiIgXCWnjmLRffenj/kT2FLAG3ob7a3/+eQtbi2r0Db+O9eLZtTpbt+Xz86HCAL695JOXQEW4Xf5w5mXPHHNPzboDtJ7yDKWfW5ugRFr9QOrV7R7iA+gm9NOgQEZGipAAlIgWh54akDsNDiYgmLF99zsf6A24iiVDnLEcwbvvfUygaJXDzl3HvaIJEorNhxH/tbOGap19JOUe528WqmZOY1Ud42ndMFSfd8D2MS8UJfQm4DcE0Q1TAbVhSOyJHIxIRkaGmd0kRybtowrJoXRsL17Syak+IYNxiORoeLlrTyqJ1bUQTpbmOpLOz2wE3oYRxvqdQNEr5tZfjfmYLpkt4+sWeVj7TLTyVuV2sPP9tnF9VkXpuIG5ctJ06g/Kf/BqvX622+1Nf7U/rjdMNNNQEqKvy5mpIIiIyxBSgRCSvOsLD6r3Dd0PSjs5uoUT/M2xd9xTCWgI3fwmz/3W6ftcvX97Pp7a1pHxfwGX43axTOe/zXyN21lysL4A1BusLEDvnAsJf+y98X/kBeFSUMJAltRUE0ljPNOcEP8vnVJb07KmIyHCjd0uRItWz5K0Mvwvmt+wvqpK3jvAQjPd/XClvSJpOZ7eOPYV+VfN33M8+kRKe7tp7gP94sjnleJ/LcP+5J/PusWXEnv4b4SU34HwHI+luRpWXhuoAq/cG+/2ZdRuYe7yf+983tih+D0VExDnNQIkUod5L3gyhhCm6krdMwkOpSaezW8eeQt7G+8Ae/a7/bT7AJ57YS9dn0msMvzl7IvXjR2KsxfPkpmwOe1gyxrB8TiULasoo95geb6IGS8Bl+cDEMu59r8KTiEgpUoASKTKlVvKWSXgoNel2dgvGLJ5tmzpnn+57pY0rm1LDk8fAfWdP5MLjRh29Maq5p2zwugy3za1k1fwqFk4MUO4xGKDcY3hPVZzlZ4S5fd6Y1GYfIiJSMlTCJ1JkSq3kLZPwUGocdXazlnMOvcjVzX/kwv1PQiICwP2vvsFHt76cEkLdBu45ayIXHT8q9T68/qyOezgzxjBjnI87LxibcvuuXbvyNCIRERkqmoESKTKlVvIWcKf3KX0pbkg6UGc3TyLGf+/4KQ89uYwP7nuc8kQEA/zhtTe4fMseumYvF/A/M07in04YnXIfFohNn5n9wYuIiAwzClAiRabUSt7SaQtdqhuS9tvZzVru2Lmcha1bOSYRxn2kUO+Bvx/kXx5/mVi38PSrGSdxyYmVPe/H5SI6/9LsD15ERGSYUYASKTKlVvKWTlvoUt2QtKOzW8DV7bUKh5j980v55v/+kh/u2Nt5c+M/3uTSx/cQ7bK+zQC/rKvhXyf0DE8WiE+tIzHptBw9AhERkeFDAUqkyJRayVtHeChz939cmbt0NyTt6Ow2Z0ycgMviAkZE3uL6uy5mw3PPsutwmG/s/DuP7H+Lta+/yYf+9hKRbh0WV0yv5iPVY3rctwXs2PGErr4R1BFORERk0BSgRIpMqZW8DdQW2gWUuw0LaspKekNSr8vwrbdFWH5GmIXVXl7ZvJTf7N2XcsyjrW/xwcdeItwtPC2fVs2VJ6U2M7CANS7itWfR/t3/BW/pBU8REZF8UBc+kSKzpLaCtS1h2h2U5hVLyVtHW+imjo2BW8IEY5Yyj6F+QoCltSOoK+BOgtliDJxekeCSF//AC2+8yc63jrYddxtY9vw/CHULT7e+YwKLJh4NTxbAFyA2fSbRhstUticiIpJlClAiRaaj5G313mC/rcyLreStr7bQw048hu+P9/C7195Ivd1CsNueXj8640QWn1yVepwvwOEVa3I7RhERkWFMJXwiRUYlbyXMWk5a+UuwCX7/2sGUL3XPyldPHse5lcekbJRsjVGrchERkRxTgBIpQh0lb6vmV7FwYoByj8FgCbgsF08s44GGKm6fNwavS+GpmJS/+hKjnt/GS4fDPHkw2Odxozwubn5xHzPXP88nt7Uc/YLXr1blIiIiOaYSPpEi1b3kbdeuXQBMmTIhn8OSQRi/eS2uWJT/e/Vgv8cdjB3dCeyOvfv5Ye2JWF8AX90srXkqENZatnas6WsOE4pbAm5DfbWfpbUV1FV5NTssIlKkFKBERArEyBeexljL/zQfcPw9bxsRoMzj4rfjZtBw1XVqVV4AYglYtK6NxuYQoZjt3Pg6GLes3BNibUuYhuoAy+dUapZYRKQIqYRPRKRAuKIRXngrxPZDoQGPLXe7+PCEStbMmswvj53LfQ3XgEefieWbtfD1532s3hukvUt46pAA2mOW1XuDLF7flrKGTUREioPebUVECkTC6+P+V//R59f9LkPDsSO57MTRLDh2JOUeN4ddPj439RN8aqSbqfe+xmvtRy/ZTyh38ZUzR/Kvp5ThcunzsqHwzFsu1h9wE+qenLoJxqGxOURTa5QZw6BFv4hIKdE7qohIgXjzlDNYcNyolNsMMH98BXecWcOr82v57Tknc8mJlZR73MQwTJl9C7jd/OCpt1LCE8Cr7Qk+veEN3vHb12mPDXBFL1lxd4uHsMOnOhSz3Lr9rdwOSEREsk4zUCIiBeL1897H2194mj+cczI3v/g67xhVxpenHMu4QM+9vCK4qD/rBt7yjeh3PzCAlsNxzvm/13nqQ+M1E5Vjj7a5sThb15QAHmwZuFxTREQKiwKUiEiBaD/hZA6eOo0F5ikuPH5U38e5fDSddA7nzXwHjz592NF9txyOc88LQS4/9ZhsDVd64XT2qUMwpjVQIiLFRh9FipQ4ay1b9kW44uH9HH/Xq1Te8QrH3/UqVz68n637IlrEXkiM4eWF/0asbjbWF8B266hnjcH6AnjPOZ/pX72B/32x772ierPsiUPZHK30wp/mu2qZR134RESKjWagjjDGTAC+CcwHxgKvAb8HvmGtbcvj0EQyFk1YFq9XO+Wi4vYQXnw90d078Tbei2fbZoiGwesnNn0m0YbLOvd6erU9vemOV9oHqPWTQTu/Ms6f9zsr43MB9RMCuR+UiIhklQIUYIyZDGwExgN/AHYC5wCfBeYbY2Zba/fncYgiabM2GZ5W7w32ukYmtZ0y3Da3Uht7FgpjSEyeSnjJDYTzPRZJy+UTYjzaNnAXPoCA27CkdkTuByUiIlmlEr6k/yIZnj5jrf2AtfZL1tp3ATcDbwO+ndfRiWRga2uUxubQgA0GurZTFpHBOX1Egjlj4pS5+z+uzA0NNQHqqno2CBERkcI27APUkdmn9wF7gJ92+/LXgcPAR40xWnktReXW7YcIOVygrnbKPRXD2rETytP7J/zE8gGu6h0qhucmX4yBb5waYUFNGeUe0+NN1gWUuw0LaspYPkezviIixUglfHDBkb/XWmtTii6stYeMMRtIBqzzgD8P9eBEMvVgcxinK2TUTjlVsawd+/KZI1m64Q3Hx193ZsWgz1ksz00+eVzJktim1ii3bD/E2pYwwZilzGOonxBgae0I6rR5rohI0TLD+ZNCAGPM94AvAF+w1v6gl6/fCnwa+JS19mfdv37w4MFen8Bdu3Zle6giaTnn0TLH+9EAGCx/Oz+9rm6lyFr46nM+1h9wE0r0/fwFXJY5Y+J8620RUiYRrKX81ZcYv3ktI194Glc0QsLr481TzuD18+ppP2EiZGnWIZGAi7cE+HvEQL+vteU4n+UPZ4UYzDZQg35uRERECsCUKVN6vX3UqFGO3rU0AwUdm60c7OPrHbePzv1QRLLH78LRQvauxws885ZrwIAAEEoY1h9w8+xbLk6vSCSDU8sLTPzdCnxvJht3dtyDOxph9M4mRr7wNAdPncbLC/8N3IP/59flgt/UhbikKcDfIx23pqQ5AI7zWX5TN7jwBIN4bkREREqIAlSO9JVsC1XHjFmxjVuO6v4azm/Zz8o9IUdlfC6goaacKVMm5G6AReLbD+8nknBWzhhJGP5wsJIPvGMk/l98G8/f1oFN9DoXZKzFHY1QuespKh7+DeHF1/eYicr09/DZtyW454Ug337izZTW5ieWe/hKXQUfnpKdJZzpPDfhhOEPB0fzgbqxWTl3sdC/pcVPr2Hx02tYGgr5dVSAOjrDNKqPr3fc/kbuhyKSPUtqK1jbEqbdQSMJtVM+Ku21Y81B/CtuxbP1UYwd+DtNJIynaSPR3TtJTJ46qLF2cLlcXH7qMVx+am573aTz3Fhg5cshogk7bNdCiYhIaVLRDjx35O9T+/h6R+x9fgjGIpI1M6q8NFQH1E45TaF4eutCz2h7AU/TRkwsjTbw0TDeNfelObL8S/e5iVtYvL5tWHflExGR0qMABQ8f+ft9xpiU58MYUwHMBtqBzUM9MJHBMMawfE6l2imnKeBO43mwlm+9dB9E0utgaKzF8+SmNEeWf2k9N0dojzERESk1wz5AWWtfBNYCE0l22+vqG8AxwH9baw8P8dBEBs3rMtw2t5JV86tYODFAucdggHKP4eKJZTzQUMXt88aoxKqL+mq/o38YPfEoq5+6iXltz6TR67CLaDiT78orp89NV9pjTERESo3WQCV9CtgI/MQY825gB3AuyT2inge+ksexiQyKMYYZ43zcecHwWsyfKSdrxzzxKDv/9nlOCu/PLDwBeP2ZfmfepLOurkMx7zFmrWVrx15OzWFCcUvAbaiv9rO0toK6Kq9mbkVEhqFhPwMFnbNQZwF3kgxOnwcmAz8GzrPW7s/f6ERkKA24dsxaVm7//qDCkzWG2PSZmQ4xbzqem3QF0whchSKasCxa18bCNa2s2hMiGLdYjm4YfNGaVhatayOaKL7HJiIig6MAdYS1ttla+3Fr7fHWWp+19iRr7eestW35HpuIDJ2B1o6de+hF3tX2bOYzTwBeP9H5lw7mHvKi47kZoC9JD2We4pqlsdayeH0bq/cGaY/ZHp0HE0B7zLJ6b1BNMkREhiEFKBGRbnpdO2Ytcw6/yB92/gS342bePVmPl1jdLBKTTsviiIeO12V4/0l+xwHSBdRPSH/WKp+2tkZpbA4RjPd/XDCuJhkiIsORApSISC+6rh179V/H0x7+JX/etoyx7ZmV7lnAulzEZryT8FXX9dhEt5gsPWOk41mlYtxj7Nbthwg5LDtUkwwRkeFHTSRERPpjLf4Vy/A0bcBEMuucZwE7ZjyhT32NxJTa7I4vDzrWQq3eG+x3lqZY9xhLezPlIm2SISIimdEMlIhIP1y7dyQ3yh1EeIqfPoP2m+8rifAEpb/HWLobBhdjkwwREcmcZqBERPrhbbxvcHs2GReRf16UvQEViI51Yk0dbb5bwgRjljKPoX5CgKW1I6gb58v3MDMScBuCaYSoYmuSISIig6MAJSLSD8+2TZgMu6wlZ5/qiJ/8Nrbui5TcfkKlusdYfbWflXtCjsr4irFJhoiIDI4ClIhIfyKRjL7NAnbseA59ZhmL179BY3OIUJeW2B37Ca1tCdNQHWD5nMqsDVkGJ50Ng4uxSYaIiAyO1kCJiPTHl14ZWke3vXjtWRz+zt0s3vRWGvsJZWvQMhgDbqZ8RLE2yRARkcFRgBIR6Uds2kysw/K6ZLe9cQSv/y9C13yfrW/YtPYTevYt/ZNcCEq9SYaIiAyO3q1FRPoRbbgUvH5nB/sChJZ8s3OT3HT3E7r7FVVVF4peN1MGyj2GiyeW8UBDFbfPG4PXpfAkIjLc6N1aRKQfiUlTidXNGnAfKOvzE6ub1RmeIP39hB45MEDNmAypUm2SISIig6MZKBGR/hhD+KrriNXNxvoCPcr5rDFYX4BY3WzCV10HXb6e7n5CYadpS0RERPJGM1AiIgPxeAgvvp7o7p14G+/Fs21zcm8or5/Y9JlEGy5LmXnqkO5+Qn59pCUiIlLwFKBERJwwhsTkqYSX3IDTbXXT3U/onWMG6DYhIiIieafPO0VEcmRJbQUBj7MmAwG34fITYzkekYiIiAyWApSISI6ku5/Q20doEZSIiEihU4ASEcmR9PcTyscoRUREJB1aAyUikkMd+wk1tUa5Zfsh1raECcYsZR5D/YQAS2tHUDfOl+9hioiIiEMKUCIiOab9hEREREqHSvhEREREREQcUoASERERERFxSAFKRERERETEIQUoERERERERhxSgREREREREHFKAEhERERERcUgBSkRERERExCHtAyWSIWstWzs2R20OE4pbAm5DfbWfpbUV1FV5Mcbke5giIiIikkUKUCIZiCYsi9e30dgcIhSzJI7cHoxbVu4JsbYlTEN1gOVzKvG6FKJERERESoUClEiarE2Gp9V7gwTjPb+eANpjltV7gyxeD7fNrdRMVC5Yi2v3DryN9+LZthkiEfD5iE2bSbThMhKTTgM97yIiIpJlClAiadraGqWxOdRreOoqGIfG5hBNrVFmjPMNzeCGi2iUwM1fxr2jCRIJOmNSJIxnyzo82zYTq5tF+KrrwKN/5kRERCR71ERCJE23bj9EKGYdHRuKWW7d/laORzTMRKOUX3s57me2YLqGpyOMtZhICE/TBvwrloF19lqJiIiIOKEAJZKmB5vDnWueBpIAHmwJ5XI4w4u1BG7+Emb/6z2CU3cmEsbTtBHX7p1DMjQREREZHhSgRNIUiqc3oxF0OFslA3Pt3oH72ScGDE+domG8a+7L5ZCKmrWWLfsiXPHwfo6/61Uq73iF4+96lSsf3s/WfRGsZu9ERER60OIAkTQF3IZgGiGqzKNGBllhLb7f3AbW6fxfspzP8+QmwjkcVrFSJ0kREZHMaAZKJE311X7HvzguoH5CIJfDGR5iMfzL/xP3jibns08doopP3XXtJNneJTx1SO0k2aaZKBERkS4UoETStKS2goDDWaWA27CkdkSOR1TirMW/Yhmepg3phycArz/bIyp6mXSSFBERkSQFKJE0zajy0lAdoMzd/3FlbmioCVBX5R2agZUo1+4deJo2YiLpzyRZIDZ9ZvYHVeTUSVJERCRzClAiaTLGsHxOJQtqyij3mB6/RC6g3G1YUFPG8jnaRHewvI33ZV6G53IRnX9pdgdUAtRJUkREJHNqIiGSAa/LcNvcSppao9yy/RBrW8IEY5Yyj6F+QoCltSOo0+a5WeHZtgmTwRocC8Sn1pGYdFr2B1Xk1ElSREQkcwpQIhkyxjBjnI87Lxib76GUtkgk7W+xgB07ntDVN4JmAHtIt5Oky8DWfRHqqryaURURkWFPAUqkwFlr2dox09UcJhS3BNyG+mo/S2srSv+i1ueDNNY/WSBeexahq28Cj0fPXy/qq/2s3BNyXMYXt3DRmla1NRcREUFroEQKWjRhWbSujYVrWlm1J0QwbrEc3avnojWtLFrXRjRRuiVWsWkzsQ4DjgXip88gdM33wePR89eHdDpJdlBbcxERkSQFKJEClf5ePfkYZe5FGy513orcFyDyz4sA7XXUH6edJLtTW3MREREFKJGCle5ePc++VZq/zolJU4nVzcL6+g9R1ucnVjers2mE9jrqW/dOkulQW3MRERnuSvOKS6QEpLtXz92vlOiSRmMIX3UdsbrZWF+gRzmfNQbrCxCrm034qus6m0Zor6P+dXSSXDW/Kq03ArU1FxGR4a5Er7hEil+6e/U8ciDNeqxi4vEQXnw90d078Tbei2fb5uTeUF4/sekziTZc1qNdufY6GlhHJ8l0ixfV1lxERIYzBSiRApXuXj1hp2mhWBlDYvJUwktuwElPPu115Fy6bc3L0iz7ExERKSUKUCIFKt2LWn+WC3KLvf23QoFz6bQ1dwH1EwK5HpKIiEjB0hookQJVX+13/AvqAt45ZoBuCWkohfbf6T5/wzkUpNPWPOA2LKkdkeMRiYiIFC4FKJECle5F7eUnxrJy3lJp/61Q4JzTtuZlbmioCVBX5R2agYmIiBQgBSiRApXuRe3bR2RnEVSptP9WKHCue1vz7m8MLqDcbVhQU8byOZUFXbopIiKSawpQIgUq/Yva7Jy3VNp/KxSkp2tb84UTA5R7DAYo9xgunljGAw1V3D5vDF7X8H6eRERE1ERCpIB1XNQ2dTRzaAkTjFnKPIb6CQGW1o6gbpwvq+cspfbf+Xj+illHW/M7Lxib76GIiIgULAUokQI31Be1pdb+W6FAREREskklfCKSIuBOr0RrOLf/FhERkeFHAUpEUqj9t4iIiEjfFKBEJIXaf4uIiIj0TWugZFiw1rK1o5FAc5hQ3BJwG+qr/SytraCuyjvsu7B16Gj/vXpvsN9W5mr/LSIiIsORApSUvGgiuTFsY3OIUJeNYYNxy8o9Ida2hGmoDrB8TqVaNHO0/ffi9fR4ziA5bR1wGxpqAkfbf1uLa/cOvI334tm2GSIR8PmITZtJtOEyEpNOI2t91kVERETySAFKSpq1yfDU12xKAmiPWVbvDbJ4Pdw2V/sBQZrtv2Mx/CuW4WnaCNEwxh7pyhcJ49myDs+2zcTqZhG+6jrw6J8cERERKW66mpGStrU1SmNzqN9SNIBgPDnb0tQaZYb2BQIctv+29kh42oCJhHveh7UQCeFp2gArlhFefL1mokRERKSoqYmElLRbtx8i5HCfolDMcuv2t3I8otLi2r0DT9PGXsNTVyYSxtO0EdfunUM0MhEREZHcUICSkvZgczhl/U5/EsCDLaFcDqfkeBvvg2j/4alTNIx3zX25HZCIiIhIjilASUkLxZ3NPnUIOpytkiTPtk1H1zwNwFiL58lNOR6RiIiISG4pQElJC7jTW29T5nD/IzkiEknveKezVSIiIiIFSgFKSlp9td/xD7kLqJ8QyOVwSo8vzYYbXn9uxiEiIiIyRBSgpKQtqa0g4HBWKeA2LKkdkeMRlZbYtJlYh131rDHEps/M8YhEREREcksBSkrajCovDdUBytz9H1fmhoaaAHVV3qEZWImINlzqfFbJ6yc6/9LcDkhEREQkxxSgpKQZY1g+p5IFNWWUe0yPH3gXUO42LKgpY/kcbaKbrsSkqcTqZmF9/Yco6/MTq5tFYtJpQzQyERERkdzQRrpS8rwuw21zK2lqjXLL9kOsbQkTjFnKPIb6CQGW1o6gTpvn9mQtrt078Dbei2fb5mTDCJ+P2LSZRBsuS4YhYwhfdR2sWIanaSNEwyld+awx4E2Gp/BV12kTXRERESl6ClAyLBhjmDHOx50XjM33UIpDNErg5i/j3tEEiQSdsScSxrNlHZ5tm4+GIo+H8OLrie7eeTRsRcPJ4DS9S9gSERERKQEKUCKSKhql/NrLMftfp7f5ImMtREJ4mjbAimWEF18PxpCYPJXwkhtQo3IREREpZVoDJSJHJRIEvr20z/DUlYmE8TRtxLV755AMTURERKQQKECJSFIsRuD7X8T90s4Bw1OnaBjvmvtyOSoRERGRgqIAJSJgLf4Vy3A/2+Q8PJEs5/M8uSlnwxIREREpNApQIsOdtXjWr8bz2F8xNpH+90e16klERESGDzWREBnOYjH8K5bh+dtfIZPwBM430hUREREpAQpQIsPVkbI9T9MGTCKz8GSB2PSZ2R2XiIiISAFTCZ/IMOXavQNP00ZMZBAleC4X0fmXZm9QIiIiIgVOAUpkmPI23jeo9UsWiE+t0ya58v/bu//guOv7zuPPt7QrrYx9roNJ0xz2JaYQPOcZHJGS2mRsKNcYudc01zZwN0NJmeKeLzX5MeXanC8QyLWUXK9pG7ieU0jio7RXSNpOycUC37WNXTBuBit26ykkBNexHAiNg881tXa1K33uj10VoR/2V5Z2V6t9PmY0q/3sV5oXXrTalz6f7+crSVJbsUBJbSp36OnqRXHPQwLShW+k+JFfg5jJvn2SJEmtzQIltavh4fP6sgSMrFrNmf/6B5DPz20mSZKkec4CJbWrrq4Zf0mKDkbWvIPinb8DOfegkSRJ7ccCJbWpyhXrSBmX3yWq5anyzmsofuRel+1JkqS2ZYGS2lS574bs13Dq6KD0c/+R0n+405knSZLU1ixQUpsaXbWaSu96UtfZS1Tq6qZy1TVU3nV9g5JJkiTNXxYoqV1FUNqynUrv1aSuwqTlfCmC1FWg0ns1pS3bXbYnSZIEuBZHame5HKWtd1A+8hz5/kfIHdpfvTZUvpvK2nWU+270Ok+SJEnjWKCkdhfB6CWrKW27i/O/rK4kSVJ7cAmfJEmSJGVkgZIkSZKkjCxQkiRJkpSRBUqSJEmSMrJASZIkSVJG7sKnlpdS4sCJMvcdPs3uwRLFkUShM9i0opvb1iyhd3me8BpGkiRJmgMWKLW08mhi696T9A8WKVYSo7XxoZHEY0eL7D5eom9FgR0blpHvsERJkiRpdlzCp5aVUrU87To2xJlx5WnMKHCmkth1bIite0+SUmpGTEmSJC0gFii1rAMnyvQPFhkaOftxQyPQP1hk4ES5McEkSZK0YLV9gYqIfER8KCI+HxEHI2I4IlJE3NrsbDq7+w+fpljJNqtUrCTuP/xqnRNJkiRpofMcKLgA+K3a5y8D3wFWNC2NMntisDRp2d50RoEnjhfrGUeSJEltoO1noIAzwGbgzSmlNwGfa3IeZVQcmdk5TUMZZ6skSZKk6bT9DFRKaRjob3YOzVyhMxiaQYnqybXYLnwp0XHkWfL9j5A7tB+Gh6Gri8oV6yj33cjoqsvB7dklSZIaqu0LlFrXphXdPHa0mGkZXwew6eJCvSPNnUqF7gfuITewD8olYmwHweESuWf2kDu0n0rvekpbtkPOH2NJkqRGcQmfWta2NUsoZJxVKnQG29YsrnOiOZJSrTw9RQwXXytPNZESMVwkN/AU3Q/cA27PLkmS1DDhtXFeLyLuAj4ObEkpPXiu40+dOjXlP+Dzzz8/x8k0UUrwsa93sfeVToqj0xepQkdiwxtG+JW3DbfEirdF3z7CDz78G3SWh8957Ei+i2/edDtn/vlbG5BMkiSp9V166aVTji9dujTTO8UFMQMVEUdrW49n/Xi42Zk1exFw92XDbHjDCIWORDBhpob0T+Xp7staozwBvHH/bjoq2a5Z1VEpc9H+3XVOJEmSpDEL5eSJF4CZ7FH9Yr2CjJmu2c5XYzNmrZYb4JHLEgMnytx3+DS7j5cYqiR6csGmi3u4bc1iei/qanbEGbnghcOTlu1NJ1Ji2QuH6br00pZ+DlXlc9j6fA5bn89h6/M5XBjm8/O4IApUSum6ZmdQ80QEV17Uxc5rL2x2lLkxfO6le69TLtUnhyRJkiZZEAVKWlC6umB4BqUo312/LPNQSokDYzOOgyWKI4lCZ7BpRTe3rVlC7/I80SrrNSVJUsuxQEnzTOWKdeSe2ZNpGV+KoLJ2XQNSzQ/l0cTWvSfpHyxSrKR/2sJ+aCTx2NEiu4+X6FtRYMeGZeQ7LFGSJGnuLYhNJGYrIj4aETsjYifw3trwLWNjEXFr89Kp3ZT7bsg+q5Tvpnz9DfUNNE+kVC1Pu44NcWZceRozCpypJHYdG2Lr3pO4w6gkSaoHZ6Cqrgc2ThhbX/sYc84tzaW5MLpqNZXe9bXrQE2/lC91dVPpXc/oqssbmK55Dpwo0z9YZGjk7McNjUD/YJGBE2WubLENRCRJ0vznDBSQUrompRRn+fjZZmdUG4mgtGU7ld6rSV0F0oTzeVIEqatApfdqSlu20zL7s8/S/YdPU6xkm1UqVhL3H361zokkSVI7cgZKapSU6DjyLPn+R8gd2l/dba+ri8oV6yj33VidSRorQ7kcpa13UD7y3GvHl0uQ76aydtzxbeSJwdKkZXvTGQWeOD6TKxtIkiRlY4GSGqFSofuBe8gN7INy6bUNIoZL5J7ZQ+7Qfiq966szSrnaj2UEo5esprTtLtyoHIojMzunaSjjbJUkSdJMuIRPqreUauXpKWK4OGl3vUiJGC6SG3iK7gfuATc/mFKhc2ZLFXty7bG0UZIkNZYFSqqzjiPPkhvYd9YNIQBiuERuYB8dR55rULLWsmlFd+YXrA5g08WFesaRJEltygIl1Vm+/9Hq+UtZlEvkH3+0voFa1LY1SyhknFUqdAbb1iyucyJJktSOLFBSneUOPZ3porhQXc6XO/h0nRO1piuX5+lbUaCn8+zH9XRC38oCvcvzjQkmSZLaigVKqrfh4Zkdn3W2qs1EBDs2LGPzyh4W5WLSi1cHsKgz2Lyyhx0blhFtsr27JElqLHfhk+qtqwvOcf7T6+S765elxeU7ggc3LmPgRJn7Dp9m9/ESQ5VETy7YdHGB29YspteL50qSpDqyQEl1VrliHbln9mRaxpciqKxd14BUrSsiuPKiLnZee2Gzo0iSpDbkEj6pzsp9N2SfVcp3U77+hvoGkiRJ0nlzBkqqs9FVq6n0rq9dB2r6pXypq5tK73pGV13ewHQLW0qJA2PL/QZLFEcShc5g04publuzhN7lec+VkiRJM+IMlFRvEZS2bKfSezWpq0Ca8IY9RZC6ClR6r6a0ZTv4hn5OlEcTt+45yXseP8GXjhYZGkkkYGgk8djRIj/++Alu3XOS8qgXLpYkSdk5AyXNRkp0HHmWfP8j5A7tr+6419VF5Yp1lPturM4mRUAuR2nrHZSPPPfaseUS5LuprB13rOZESomte0+y69gQQyOTHx8FzlQSu44NsXUvPLjRXfskSVI2FijpfFUqdD9wD7mBfVAuvbZJxHCJ3DN7yB3aT6V3fXVWKZeDCEYvWU1p2124UXl9HThRpn+wOGV5Gm9oBPoHiwycKHOlu/dJkqQMXMInnY+UauXpKWK4OGmHvUiJGC6SG3iK7gfugYwX0tXcuP/waYqVbP/mxUri/sOv1jmRJElaKCxQ0nnoOPIsuYF9Z90UAiCGS+QG9tFx5LkGJRPAE4MlRjMeOwo8cbxYzziSJGkBsUBJ5yHf/2j1HKYsyiXyjz9a30B6neLIzGb8hjLOVkmSJFmgpPOQO/R0pgvjQnU5X+7g03VOpPEKnTPbEKIn5wYSkiQpGwuUdD6Gh2d2fNbZKs2JTSu6M7+4dQCbLi7UM44kSVpALFDS+eia4Y5t+e765NCUtq1ZQiHjrFKhM9i2ZnGdE0mSpIXCAiWdh8oV6yZdEHc6KYLK2nV1TqTxrlyep29FgZ7Osx/X0wl9Kwv0Ls83JpgkSWp5FijpPJT7bsg+q5Tvpnz9DfUNpNeJCHZsWMbmlT0sysWkF7oOYFFnsHllDzs2eBFdSZKUnQVKOg+jq1ZT6V1P6jp7iUpd3VR61zO66vIGJdOYfEfw4MZlfOn65bznLQUW5YIAFuWCn3hLD/+7bzmfveYN5DssT5IkKbtcswNILSmC0pbt8MA95Ab2Qbn0ul35UgTkq+WptGU7OMPRFBHBlRd1sfPaC5sdRZIkLRAWKOl85XKUtt5B+chz5PsfIXdof3W3vXw3lbXrKPfd6MyTJEnSAmOBkmYjgtFLVlPadhduVC5JkrTweQ6UJEmSJGVkgZIkSZKkjCxQkiRJkpSRBUqSJEmSMrJASZIkSVJGFihJkiRJysgCJUmSJEkZWaAkSZIkKSMLlCRJkiRlZIGSJEmSpIwsUJIkSZKUkQVKkiRJkjKyQEmSJElSRhYoSZIkScrIAiVJkiRJGVmgJEmSJCkjC5QkSZIkZWSBkiRJkqSMLFCSJEmSlJEFSpIkSZIyskBJkiRJUkYWKEmSJEnKyAIlSZIkSRlZoCRJkiQpIwuUJEmSJGVkgZIkSZKkjCxQkiRJkpSRBUqSJEmSMrJASZIkSVJGuWYHkLQwpJQ4cKLMfYdPs3uwRHEkUegMNq3o5rY1S+hdnicimh1TkiRpVixQkmatPJrYuvck/YNFipXEaG18aCTx2NEiu4+X6FtRYMeGZeQ7LFGSJKl1uYRP0qykVC1Pu44NcWZceRozCpypJHYdG2Lr3pOklJoRU5IkaU5YoCTNyoETZfoHiwyNnP24oRHoHywycKLcmGCSJEl1YIGSNCv3Hz5NsZJtVqlYSdx/+NU6J5IkSaofC5SkWXlisDRp2d50RoEnjhfrGUeSJKmuLFCSZqU4MrNzmoYyzlZJkiTNRxYoSbNS6JzZrno9OXfhkyRJrcsCJWlWNq3ozvxC0gFsurhQzziSJEl1ZYGSNCvb1iyhkHFWqdAZbFuzuM6JJEmS6scCJWlWrlyep29FgZ7Osx/X0wl9Kwv0Ls83JpgkSVIdWKAkzUpEsGPDMjav7GFRLia9qHQAizqDzSt72LFhGRGeAyVJklpXrtkBJLW+fEfw4MZlDJwoc9/h0+w+XmKokujJBZsuLnDbmsX0XtTV7JiSJEmzZoGSNCcigisv6mLntRc2O4okSVLduIRPkiRJkjKyQEmSJElSRhYoSZIkScrIAiVJkiRJGVmgJEmSJCkjC5QkSZIkZWSBkiRJkqSMLFCSJEmSlJEFSpIkSZIyskBJkiRJUkYWKEmSJEnKyAIlSZIkSRlZoCRJkiQpIwuUJEmSJGVkgZIkSZKkjCxQkiRJkpSRBUqSJEmSMrJASZIkSVJGFihJkiRJyihSSs3O0NJOnTrlP6AkSZLU4pYuXRpZjnMGSpIkSZIyskBJkiRJUkYWKEmSJEnKyAIlSZIkSRlZoCRJkiQpI3fhkyRJkqSMnIGSJEmSpIwsUJIkSZKUkQVKkiRJkjKyQEmSJElSRhYoTRIRl0bEL0fEn0fEYEQMR8TLEfGnEXFts/Pp3CIiHxEfiojPR8TB2nOYIuLWZmfTZBFxcUR8LiJejIhSRByNiN+KiGXNzqZzi4ifjoj7IuIvI+Ifaj9rDzc7l7KJiAsj4taI+JOI+GZEDEXEqYh4MiJ+LiJ8r9QiIuKTEfFntfcuQxHxSkR8LSI+HhEXNjufZi4ibqq9ps6r9zDuwqdJIuIPgRuBvwWeBF4B3ga8B+gEPpRS+nTzEupcIuL7gJO1uy8Dw8AKYEtK6cFm5dJkEXEJsA94I/CnwHPAVcC1wNeBq1NK32teQp1LRBwErgBeBY4DlwO/n1K6qZm5lE1EbAX+B/AS8BfAMeD7gZ8ElgJ/BLwv+YZp3ouIYWCA6vuXvwcuAH4YeAfwIvDDKaXB5iXUTETECuBvqL73XMw8eg+Ta3YAzUuPA59MKX1t/GBEbAT+D/DrEfGFlNJLTUmnLM4Am4GDKaWXIuIu4OPNjaRp/A7V8vTBlNJ9Y4MR8SngI8CvAlublE3ZfIRqcfomsJHqm3C1jm9Q/QPhl1NKo2ODEbEd+CrwU1TL1B81J55m4J+llIoTByPiV4HtwH8CPtDwVJqxiAjg88D3gD8Gbm9uotdzWlqTpJR2TixPtfE9wFeALmB9o3Mpu5TScEqp35I7v9Vmn94NHAX++4SHPw78I/AzEXFBg6NpBlJKf5FSet4ZitaUUvrzlNKXxpen2vh3gB21u9c0PJhmbKryVPNo7fbSRmXRrH0Q+BHgFqq/C+cVC5Rmqly7rTQ1hbQwjJ1TuHuKN2+ngaeARVSXoEhqPH/nLQw/Xrv966amUCYRsRq4F/jtlNLeZueZikv4lFlE/AvgOqrLw+bl/9BSi3lb7fYb0zz+PNUZqsuAP2tIIkkAREQOuLl29/FmZtHMRMTtVM+ZWUr1/Kd3US1P9zYzl86t9nP3e1TPRdze5DjTskApk4joBn4f6AZ+KaV08hxfIuncltZuT03z+Nj499U/iqQJ7gXWALtSSk80O4xm5HaqG4GMeRz42ZTSd5uUR9ndCbwdeFdKaajZYabjEr4FqrYNcprBx7Rb7kZEJ9W/BlwNPAL8t0b9d7SzuXwOJUnZRcQHgV+kuivmzzQ5jmYopfSmlFIAb6K6Acgq4GsR0dvcZDqbiHgn1Vmn30gpPd3sPGfjDNTC9QIw3cmUU3lxqsFaeXoYeB/VkzBv8kTphpmT51Dz2tgM09JpHh8b/3/1jyIJICK2Ab9NdSvs61JKrzQ5ks5TSull4E8iYoDqUumHqM4qap6pLd17iOrzdEeT45yTBWqBSildN9vvERF5qsv23gf8AXBzSmlktt9X2czFc6h57+u128umeXxsx6jpzpGSNIci4sPAbwKHqZanv29uIs2FlNK3IuJvgbURsTyldKLZmTTJYl77XVis7mI+yQMR8QDVzSU+3KhgU7FAaUoR0UV1xuknqP5F4JaJu4RJmrWx6wW9OyI6JlyDZgnVZbNngP3NCCe1k4j4ZarnPR0EftQ32QvOm2u3/iF4fioBn53msV6q50U9SfUPj01f3meB0iS1DSP+mOqFWD8L/LzlSZp7KaUXImI31Z32fgG4b9zDdwMXAJ9JKc27a2BIC0lE3AF8AjgAvNtle60nIi4DXk4pnZow3gH8F6oXLN/nJljzU23DiFuneiwi7qJaoP5nSunBRuaajgVKU9lBtTydAL4N3DnFVOpXUkpfaXAuzUBEfBS4vHZ3be32loh4V+3zJ+fLC1Gb+wCwD/h0RFwHPAu8k+o1or4B/OcmZlMGEfFe4L21u2+q3a6LiJ21z0+klG5vcCxlFBHvp1qeRoC/BD44xe+8oymlnQ2OppnZDPxaRDwJ/B3wPao78W2kuonEd4AtzYunhcQCpam8tXa7nOp2ktP5Sv2jaBaup/qLY7z1tY8xFqgmq81CvYPqG7jrqb4JeInqSex3+9fSlrAWeP+EsVW1D4BvUd1WWfPT2O+8TuDD0xyzB9jZiDA6b/8X+EGq13x6O9XLP/wj1T9E/R7waWcWNVfCDdUkSZIkKRuvAyVJkiRJGVmgJEmSJCkjC5QkSZIkZWSBkiRJkqSMLFCSJEmSlJEFSpIkSZIyskBJkiRJUkYWKEmSJEnKyAIlSZIkSRlZoCRJkiQpIwuUJEmSJGVkgZIkSZKkjCxQkiRJkpSRBUqSJEmSMrJASZLaSkSsjIhTEZHGfdw84ZieiPjGhGN+vVmZJUnzhwVKktRWUkrHgNsmDH8qIi4ad/8TwKXj7v8N8LF6Z5MkzX+RUmp2BkmSGi4ivgj81LihP0wp/buIuArYB3TWxoeBH0op/XWjM0qS5h9noCRJ7erfA98Zd//fRsS/AT7Ha+UJ4GOWJ0nSGGegJEltKyI2A18eN1QBcuPu7wWuTSmNNjSYJGnecgZKktS2Ukq7gM+MGxpfnv4BuHmq8hQR3RFxb0R8OyKGIuKrEbGp3nklSc1ngZIktbtfBL41xfjtKaWpxgF21r7ufwEfAsrAlyNiY10SSpLmDZfwSZLaWkT8S+AA0D3hoYdSSu+f4virgL8CPppS+mRtrAAcBl5JKV1V58iSpCZyBkqS1LYiIg88zOTyBHBzRPzkFOM/DYwCvzs2kFIqAp8Ffigi3lKHqJKkecICJUlqZ3cDa8fdPzbh8c9ExPdPGHs78EJK6eSE8a+Oe1yStEBZoCRJbSki1gO/NG5oCPhXwBfGjS0HHpzwpT8AvDTFtxwbe/NcZZQkzT8WKElS24mIxcBDvP56T3emlJ4HfgH47rjxfx0RW8bd7wFKU3zb4rjHJUkLlAVKktSOPgVcMu7+V4HfBEgpfRf4wMTjI2JV7fMhpj5nqjDucUnSAmWBkiS1lYj4MWD8jNIwcEtKaWRsIKX0ReDRcccsBh6KiA6qS/V+YIpvPTb24twmliTNJ7lzHyJJ0sKRUvoyEBmOuxG4ceJ4RBwEfiQilk3YSOKdtduDcxBTkjRPOQMlSdLMfJHq78+fHxuIiG7gFuBASunvmhVMklR/zkBJkjQDKaW/iogvAL8SEcuB54GbgbcCP9rUcJKkuouUUrMzSJLUUiKiAHwCuAl4A3AYuCOl1N/UYJKkurNASZIkSVJGngMlSZIkSRlZoCRJkiQpIwuUJEmSJGVkgZIkSZKkjCxQkiRJkpSRBUqSJEmSMrJASZIkSVJGFihJkiRJysgCJUmSJEkZWaAkSZIkKSMLlCRJkiRlZIGSJEmSpIwsUJIkSZKUkQVKkiRJkjKyQEmSJElSRhYoSZIkScro/wOv7+k2uva3MQAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 432x432 with 1 Axes>"
      ]
     },
     "execution_count": 13,
     "metadata": {
      "image/png": {
       "height": 369,
       "width": 424
      },
      "needs_background": "light"
     },
     "output_type": "execute_result"
    }
   ],
   "source": [
    "num_components = 1\n",
    "X_reconst, mean, principal_values, principal_components = PCA(X, num_components)\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(6, 6))\n",
    "# eig_vals, eig_vecs = eig_vals[:num_components], eig_vecs[:, :num_components]\n",
    "ax.scatter(X[:, 0], X[:, 1], label='data')\n",
    "for (princial_variance, principal_component) in (zip(principal_values, principal_components.T)):\n",
    "    draw_vector(\n",
    "        mean, mean + np.sqrt(princial_variance) * principal_component, \n",
    "        ax=ax)\n",
    "ax.scatter(X_reconst[:, 0], X_reconst[:, 1], label='reconstructed')\n",
    "plt.axis('equal');\n",
    "plt.legend();\n",
    "ax.set(xlabel='$\\mathbf{x}_0$', ylabel='$\\mathbf{x}_1$');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "**Principal Component Analysis for high-dimensional datasets**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
    "def PCA_high_dim(X, num_components):\n",
    "    N, D = X.shape\n",
    "\n",
    "    # Normalization of the dataset so that they have zero mean\n",
    "    X_normalized, μ = normalize(X)\n",
    "\n",
    "    # Find the covariance matrix\n",
    "    M = np.dot(X_normalized, X_normalized.T) / N\n",
    "\n",
    "    # Compute the eigenvalues (λ) and eigenvectors (V) of S in descending order of the eigenvalues\n",
    "    λ, V = eig(M)\n",
    "\n",
    "    # Compute eigenvectors for the original covariance matrix\n",
    "    XV = X_normalized.T @ V\n",
    "\n",
    "    # Choose eingenvectors with the `num_components` largest eigenvalues to be the basis of principal subspace U\n",
    "    U = XV[:, :num_components]\n",
    "\n",
    "    # Compute the orthogonal projection of the data onto the subspace spanned by columns of U\n",
    "    P = projection_matrix(U[:, :num_components])\n",
    "\n",
    "    # Reconstruct the data\n",
    "    X_reconstructed = (P @ X_normalized.T).T + μ\n",
    "\n",
    "    return X_reconstructed, μ, λ, V"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "Next, we will compare the implementation from the PCA in higher dimensions with the model from scikit-learn."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "difference in reconstruction for num_components = 1: (3.5091984330690947e-29+0j)\n",
      "difference in reconstruction for num_components = 2: (1.3971466188562764e-29+0j)\n"
     ]
    }
   ],
   "source": [
    "for num_component in range(1, 3):\n",
    "    pca = SKPCA(n_components=num_component, svd_solver=\"full\")\n",
    "    sklearn_reconst = pca.inverse_transform(pca.fit_transform(X))\n",
    "    reconst, _, _, _ = PCA_high_dim(X, num_component)\n",
    "    print(\n",
    "        \"difference in reconstruction for num_components = {}: {}\".format(\n",
    "            num_component, np.square(reconst - sklearn_reconst).sum()\n",
    "        )\n",
    "    )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 0,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
   ],
   "source": [
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (system-wide)",
   "language": "python",
   "metadata": {
    "cocalc": {
     "description": "Python 3 programming language",
     "priority": 100,
     "url": "https://www.python.org/"
    }
   },
   "name": "python3",
   "resource_dir": "/ext/jupyter/kernels/python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}