{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Multitask GP Regression\n",
    "\n",
    "## Introduction\n",
    "\n",
    "Multitask regression, introduced in [this paper](https://papers.nips.cc/paper/3189-multi-task-gaussian-process-prediction.pdf) learns similarities in the outputs simultaneously. It's useful when you are performing regression on multiple functions that share the same inputs, especially if they have similarities (such as being sinusodial). \n",
    "\n",
    "Given inputs $x$ and $x'$, and tasks $i$ and $j$, the covariance between two datapoints and two tasks is given by\n",
    "\n",
    "$$  k([x, i], [x', j]) = k_\\text{inputs}(x, x') * k_\\text{tasks}(i, j)\n",
    "$$\n",
    "\n",
    "where $k_\\text{inputs}$ is a standard kernel (e.g. RBF) that operates on the inputs.\n",
    "$k_\\text{task}$ is a lookup table containing inter-task covariance."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import math\n",
    "import torch\n",
    "import gpytorch\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "%matplotlib inline\n",
    "%load_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Set up training data\n",
    "\n",
    "In the next cell, we set up the training data for this example. We'll be using 100 regularly spaced points on [0,1] which we evaluate the function on and add Gaussian noise to get the training labels.\n",
    "\n",
    "We'll have two functions - a sine function (y1) and a cosine function (y2).\n",
    "\n",
    "For MTGPs, our `train_targets` will actually have two dimensions: with the second dimension corresponding to the different tasks."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "train_x = torch.linspace(0, 1, 100)\n",
    "\n",
    "train_y = torch.stack([\n",
    "    torch.sin(train_x * (2 * math.pi)) + torch.randn(train_x.size()) * 0.2,\n",
    "    torch.cos(train_x * (2 * math.pi)) + torch.randn(train_x.size()) * 0.2,\n",
    "], -1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Define a multitask model\n",
    "\n",
    "The model should be somewhat similar to the `ExactGP` model in the [simple regression example](../01_Exact_GPs/Simple_GP_Regression.ipynb).\n",
    "The differences:\n",
    "\n",
    "1. We're going to wrap ConstantMean with a `MultitaskMean`. This makes sure we have a mean function for each task.\n",
    "2. Rather than just using a RBFKernel, we're using that in conjunction with a `MultitaskKernel`. This gives us the covariance function described in the introduction.\n",
    "3. We're using a `MultitaskMultivariateNormal` and `MultitaskGaussianLikelihood`. This allows us to deal with the predictions/outputs in a nice way. For example, when we call MultitaskMultivariateNormal.mean, we get a `n x num_tasks` matrix back.\n",
    "\n",
    "You may also notice that we don't use a ScaleKernel, since the MultitaskKernel will do some scaling for us. (This way we're not overparameterizing the kernel.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "class MultitaskGPModel(gpytorch.models.ExactGP):\n",
    "    def __init__(self, train_x, train_y, likelihood):\n",
    "        super(MultitaskGPModel, self).__init__(train_x, train_y, likelihood)\n",
    "        self.mean_module = gpytorch.means.MultitaskMean(\n",
    "            gpytorch.means.ConstantMean(), num_tasks=2\n",
    "        )\n",
    "        self.covar_module = gpytorch.kernels.MultitaskKernel(\n",
    "            gpytorch.kernels.RBFKernel(), num_tasks=2, rank=1\n",
    "        )\n",
    "\n",
    "    def forward(self, x):\n",
    "        mean_x = self.mean_module(x)\n",
    "        covar_x = self.covar_module(x)\n",
    "        return gpytorch.distributions.MultitaskMultivariateNormal(mean_x, covar_x)\n",
    "\n",
    "    \n",
    "likelihood = gpytorch.likelihoods.MultitaskGaussianLikelihood(num_tasks=2)\n",
    "model = MultitaskGPModel(train_x, train_y, likelihood)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Train the model hyperparameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Iter 1/50 - Loss: 1.220\n",
      "Iter 2/50 - Loss: 1.180\n",
      "Iter 3/50 - Loss: 1.138\n",
      "Iter 4/50 - Loss: 1.096\n",
      "Iter 5/50 - Loss: 1.054\n",
      "Iter 6/50 - Loss: 1.013\n",
      "Iter 7/50 - Loss: 0.972\n",
      "Iter 8/50 - Loss: 0.932\n",
      "Iter 9/50 - Loss: 0.892\n",
      "Iter 10/50 - Loss: 0.854\n",
      "Iter 11/50 - Loss: 0.815\n",
      "Iter 12/50 - Loss: 0.777\n",
      "Iter 13/50 - Loss: 0.739\n",
      "Iter 14/50 - Loss: 0.700\n",
      "Iter 15/50 - Loss: 0.660\n",
      "Iter 16/50 - Loss: 0.620\n",
      "Iter 17/50 - Loss: 0.579\n",
      "Iter 18/50 - Loss: 0.538\n",
      "Iter 19/50 - Loss: 0.497\n",
      "Iter 20/50 - Loss: 0.456\n",
      "Iter 21/50 - Loss: 0.415\n",
      "Iter 22/50 - Loss: 0.376\n",
      "Iter 23/50 - Loss: 0.338\n",
      "Iter 24/50 - Loss: 0.301\n",
      "Iter 25/50 - Loss: 0.265\n",
      "Iter 26/50 - Loss: 0.231\n",
      "Iter 27/50 - Loss: 0.197\n",
      "Iter 28/50 - Loss: 0.165\n",
      "Iter 29/50 - Loss: 0.134\n",
      "Iter 30/50 - Loss: 0.104\n",
      "Iter 31/50 - Loss: 0.076\n",
      "Iter 32/50 - Loss: 0.050\n",
      "Iter 33/50 - Loss: 0.027\n",
      "Iter 34/50 - Loss: 0.006\n",
      "Iter 35/50 - Loss: -0.012\n",
      "Iter 36/50 - Loss: -0.027\n",
      "Iter 37/50 - Loss: -0.040\n",
      "Iter 38/50 - Loss: -0.051\n",
      "Iter 39/50 - Loss: -0.059\n",
      "Iter 40/50 - Loss: -0.065\n",
      "Iter 41/50 - Loss: -0.068\n",
      "Iter 42/50 - Loss: -0.070\n",
      "Iter 43/50 - Loss: -0.069\n",
      "Iter 44/50 - Loss: -0.068\n",
      "Iter 45/50 - Loss: -0.066\n",
      "Iter 46/50 - Loss: -0.063\n",
      "Iter 47/50 - Loss: -0.060\n",
      "Iter 48/50 - Loss: -0.057\n",
      "Iter 49/50 - Loss: -0.055\n",
      "Iter 50/50 - Loss: -0.053\n"
     ]
    }
   ],
   "source": [
    "# this is for running the notebook in our testing framework\n",
    "import os\n",
    "smoke_test = ('CI' in os.environ)\n",
    "training_iterations = 2 if smoke_test else 50\n",
    "\n",
    "\n",
    "# Find optimal model hyperparameters\n",
    "model.train()\n",
    "likelihood.train()\n",
    "\n",
    "# Use the adam optimizer\n",
    "optimizer = torch.optim.Adam(model.parameters(), lr=0.1)  # Includes GaussianLikelihood parameters\n",
    "\n",
    "# \"Loss\" for GPs - the marginal log likelihood\n",
    "mll = gpytorch.mlls.ExactMarginalLogLikelihood(likelihood, model)\n",
    "\n",
    "for i in range(training_iterations):\n",
    "    optimizer.zero_grad()\n",
    "    output = model(train_x)\n",
    "    loss = -mll(output, train_y)\n",
    "    loss.backward()\n",
    "    print('Iter %d/%d - Loss: %.3f' % (i + 1, training_iterations, loss.item()))\n",
    "    optimizer.step()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Make predictions with the model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAd4AAADMCAYAAADd98LnAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAABZsklEQVR4nO2dd3wUZfrAv7ObZDd9U0ijaAI2emAtWClRBAvKoSCeh3cKeoqnHmI5RH9W1DtUxHKiWE5EKaKggCJBpQpuCL2TQIB0Npu+fX5/bGE32Q2BhM0G3u/nwweYeWfmmdl55nnLUyRZlhEIBAKBQBAYFG0tgEAgEAgE5xLC8AoEAoFAEECE4RUIBAKBIIAIwysQCAQCQQARhlcgEAgEggAiDK9AIBAIBAHknDa8kiRpJEmaIElSP0mSsiRJmuCxL0uSpAUBlud1SZKe9LMvS5KkigYyPilJ0oeSJGl8tP35DMg3yvl3P0mScpzyZnjsz3A9s+Y8P085fRzbKvJ7PlPn7z2qNc4rCC6ELp+yfEKX25Bz2vACC4D5sixvlmV5JaB3/bDO/weaef52OOWZ1WDzZlmWH5Bl2eCjrYFWRJKkLGCz8/ybgTxgnizLeR7XzZNl+Q4PGZrEU04fx7aW/O5n6npOnh8YwVmD0OVmInS57TlnDa8kSf3gxA/o/PdC4BmPZhnOHtso58vq6iG6tzl7Xk+6etmuHp5z/5POv3Oc7bIkSfrcs73znE86z591ErE/BB7w+L/Gebz72j7u07Mn+qQkSa87/91Qbq/78nHt6z0V098z9dW79Xgemob33tSxznZPun4r5zaXrBNOss3nM3X+xp7PUNDOEbosdLm9cc4aXkCLo6fXCOnEdI9eluWVzh/4dee20eD+0fNwKPdKZ8+uv/PvDFmWF8qy/Iaznd75UdADxZ7tnS+Yq5feZM/SpSySYypHA+RJjh5fhizLs4CnfBzj2eNc6LHLS24f99UQjY9tDa+1mQa9W6fCPODxIfS8ZlPH9nO2W+iSTXKMYHSuXrpTSX1tO9kzbZe9ZIFfhC4LXW5XnMuGNw/fP5rGo+ds8NiudyrGNOABSZIO4niB+wHxzp7ch862mxucc4HzpY330f56/Hw0/ODqKWc5p9XyZFmeJTVYG2oGDeVoeF+nhJ+edYZTVtf5fD0rfzR8huB4Vgbnv/Oc//e3ralnqj/JtQXtC6HLQpfbFees4XX2oOI9X3LJsVj/hkczjce/Dc5eapZz/aI/jqmPn53nc62V+GI+jpdW76P9HzheYnAo88nkngW4nQqc0zI+nTh84PlxaihHw/s6Vfr52KZzntM15decZ9UUmzlxDxk4np2vbaf0TAXtG6HLQpfbGyFtLUAbMwSYIEmSq0cWL8uy5xSPa/pHw4le3aWSJIHjZVsoy3Kecw0CAOff/SRJ6ud8KZFl2SBJkt75/82e7WVZfsNj7aMfcL0kSbMaOlk0wHPaJQNHLzHDKe8op2z9JEnKcH5g/nCePwPIkiRJ43Fd1zm97svHNd3yeMg62vl8rgcMru2uKTNA63p+kiS9LsvyUw2eld6jvcbPv7Oc/9Y0OL6fLMtvOM/jb5u/Z9rUsxW0T4QuC11uN0iyqE4kaAbO6ae8kzllBDtny30IBKfL2aID7fk+ztmpZsGp4ZzO8zUF1W5wTUW2R0UVCFoLocttjxjxCgQCgUAQQMSIVyAQCASCANIi5yqtVpuFY/HcoNPp2iI7jEAgaCWEPgsEgeG0R7xarbYfoNfpdO02e4hAIHAg9FkgCBwtWuPVarUaYAKwUqfT+QqS5umnnxaLyAJBM3nttdekk7c6Mwh9FghaF3/63KKpZp1OZ9BqtQtxpDfz20t+4YUXTnqu0tJSkpKSWiJOqxNsMgWbPCBkai7Nken5558PkDS+aS19bq/PP9AImZpHe5WpKX1uyVTzk1qtVqPT6fKAO0/3PAKBoO0R+iwQBI6WjHhXAhlarVaLj4TeAoGgXSH0WSAIEKdteD3WgHyuBQkEgvaD0GeBIHCc67maBU1QVVVFRUUFFovFbxu73U5VVVUApTo57UEmSZKIjY0lMTGxDaUSnEsIfW49GsoUGhpKXFwcMTExzTpeGF6BX8rKyujYsSMqlQqPBOxeWCwWQkNDAyxZ07QHmWw2G4cOHRKGVxAwhD63Hp4yybKMyWTi2LFjwvAKWo4sy6jV6rYW46xEqVQi0rUKAonQ5zODJEmo1epT0meRMlLQKhQVFZGVlUVxcfFpnyM3N5dFixaRnZ3Nxx9/TF6eI//5okWLmDJlSmuJ6heDwcDw4cN9ynXJJZeQnZ1NdnY206dPx2AwnHF5BIK2QujzmUUYXkGrMG3aNNavX8+rr756WscbDAY+/vhjRo4cyZAhQ7j//vt59tlnARg8eHBriuoXjUZDenp6o+2ZmZmkp6czZMgQhgwZwqRJkxg7dqzPcxgMBqZPn36mRRUIzihCnx2cKX0WU82CFpGYmIjJZHL/f9asWcyaNQu1Wn1KvciFCxeSmZnptS0uLo7c3FzS09PJzc0lNzeXVatWcd9995GTk0N8fDyrVq3i9ttvZ9WqVcTHx9O3b19ycnL49ttv6d+/P+np6SxcuJC5c+fy8MMPM2nSJACv9vHx8SxcuJD+/fuTn5/fLHk1Gg0GgwG9Xs+qVauorKx0y5WTk0Nubi6xsbHufX/5y1/o0KFDs5+HQNAWCH0OjD6LEa+gRezYsYPRo0cTHh4OQHh4OGPGjGHPnj2tep3MzEwyMzMZPHgws2fPZtWqVSxcuJC+ffvy7LPPkp6eTnp6OrNnz2bgwIHExcUxadIkRo4c6T7HyJEjycjIaNR+ypQpDB482N0Tbi56vZ6MjAzuv/9++vbty+zZs90fh8zMTK99n376aas+D4HgTCD0OTD6LEa8ghaRkpJCTEwMJpMJtVqNyWQiJiaGlJSUUzrPqFGjeOihh7j//vvd2/Lz88nMzPTZ077vvvsAePbZZzGbzfTv3x+NRkNmZiZlZWVoNBp328GDBzN9+nT3MYBX+4cffpj4+HgAKioqmiWvwWAgIyOD7Oxs8vPz6d+/v9f+vLw88vPzfe4TCIIVoc+B0WdheAUtprS0lPHjx3Pfffcxe/bs03LI0Gg0TJ48mY8//pj09HTy8/N59913vdq4pqYmTZrE9OnT6du3L6NGjaJv377uqSVwhE3k5+djMBjQaDSMGjWKKVOmuJX35Zdf9mo/adIkr6mp3Nxcr2kyl9JlZ2cDeMlWWVlJfHw8+fn55OXlodfrMRgM5Ofne+1z7c/IyDjlZyMQBBKhz2den1tUnag5PP3007IoktA6BFqeAwcO0K1btybbBHuMXbDgS6aGz/f5559v0+pEzaE5+hxsegNCJhD63Jq0VJ/FGq9AIBAIBAFEGF6BQCAQCAKIMLwCgUAgEAQQYXgFAoFAIAggwvAKBAKBQBBAhOEVBAV5eXkMGDCA3NzcJrcJBILgR+hz04g4XkGzUKtVfvb42+4bo9Hkc3tGRoY77u+9994DHDF1rqwxAoGgdfCvyyD0OTCIEa8gaIiNjfW7Ly8vj48//phFixaRl5fn/r+rskh2djbDhw8nNzeXqVOnBlBqgUDgC6HP/hGGV9AsjEaTzz/V1TV+9/n6czJGjhzJokWLGmWbaZiPtWFO1SFDhrhTxlVWVp7JRyEQtGua0k+hz4FBTDULgoohQ4YwduxYJk+e3GifZz5WXzlV4+LiAimqQCA4CUKffXPahler1WqADEAL6HU63cLWEkpw7pGXl+euBtK/f39iY2PJzc1150RtmI+1YU7VH374wZ3P9dChQ4162IKmEfosaE2EPjdNS0a8dwIrdTrdLK1WexAQiio4bTIyMpg7dy6Au8YmwIYNG7zauPBUQlepsJtvvhmAJUuWBF1u13aA0GdBqyH0uWlO2/DqdLpZ4O4pb24tgQQCQeAR+iwQBI7WWON9BhjfVIPS0tKTnsRXjca2JthkCrQ8drsdi8XSZBubzRYgaZpPe5HJbrc3SzcCTIv1Odj0BoRMIPS5NWmpPrfI8Gq12lHANCAeMPhr19zSV8FWtguCT6ZAylNVVdWsKZ5gnAZqDzIpFIqger9aU5+D6b5cnOsyCX1uXVqizy1xrsrC0Tse7dx0x+meSyAQtC1CnwWCwNGSNd6VQP+TNhScNcz85WCjbTabHaWyeeHgjwzq2uT+3NxccnJySE9Pp7KykpycHF555ZVmnXv69OkMHjyYhQsX8n//93/u7QaDgbFjx7Js2bJmnedcRejzuYUvXQahz4FCxPEKggKDwcC///1vtyckQE5OTrOP1ev1ZGZmkp6e7rVPo9E02iYQCM4sQp+bRmSuEgQFCxcuZPDgwV7bJk+eTF5eHtOnTyc7O5tFixZ5pZKbMmUK4FDo/Px8srOzGTt2LOBQ3o8//tgdOwiN09T5Ohc4etuuazU8RiAQnByhz00jDK8gaNFoNDz77LPcd999DBkyhIULF3qlknN5hfbv35/09HSGDBni7g1PmTKFwYMHe/WaG6ap83Uu1xTXyJEjqaysbHSMQCA4PYQ+n0AYXkFQMGrUKFatWuW1LTs722fb5qaSi4+PB6CiosK9rX///mRmZrrXmhqey7UmBScC+RseIxAImkboc9OINV5BUKDRaNxlxFzOGIMHDyY9PZ2FCxeSnp7OqFGj3FNNBoOB/Px8cnNz0ev1Xtvy8/OZNGmSOyWdq13DNHVAo3NNnjzZfb34+PhGx5xNaesEgjOF0OemkWRZPqMXePrpp+UXXnjhpO1KS0uDLs4u2GQKtDwHDhygW7duTbaxWCxBF2PXXmRq+Hyff/55XnvtNSnQsp0KzdHnYNMbEDKB0OfWpKX6LKaaBQKBQCAIIMLwCgQCgUAQQIThFfhFkiTO9FLEuYosy0hSUM8qC84yhD6fOU5Vn4XhFfhFpVJhMBiEsp4BTCYTISHCt1EQOIQ+nxlkWcZgMKBSqZp9jNB8gV9SUlIoLi5Gr9f7VVa73Y5CEVz9t/YgkyRJdOjQoQ0lEpxrCH1uPTxlkiQJlUpFSkpKs48Xhlfgl5CQEDp16tRkG+Et2jyCUSbBuYXQ59bBbLVj0Je3SCZheAUCgUAg8IPdLlNUZeTQ8ToOldfRIVpF38SWnVMYXoFAIBAIPLDZZY5U1HOgrIb88jqMlhOF7xOjw4CWOUYKw9sMDHUWympM1Jis1Jps1JisWO0yEo75faUCIsNCiAkPIVYdSlxEKDHh3sHVRUVF3HPPPcyZM+eU1gIEgnMdoTuCQCDLMkcNRvaWVJNXVofJajv5QaeJMLw+qDZaOVBWw+4CA3X76qi3nPoPEKUKIU2jJi02nPMTIpg2bRrr16/n1Vdf5Z133jkDUgsEwcHesjoK6g0AKCQIVSoID1WgDlUSEaYkNjwUpaL5IwahO4IzSbXRyo7CKvYWV1NtsgbkmsLwOjFabOwrqWFfaQ3FlSZkZGprjURGKk/rfDUmK/tKavjT5RdgtZjd22fNmsWsWbNQq9XuChqeiN69oL2zq6Qec6n/zqokSWjCQ4mPDCMpOozOceEkRasaxUFqNBqMRqP7/yfTHYGguciyTIG+ntXbDjBt8kPc/cx/iIlvXpSBvqyEO/7+CPPmzTvtb3Rw+Wi3AVX1Fn7bX85nGwr4bX85RZVGZFovzu2Zz1aQOfAmQlVqAEJVaq4bdhtrddt8tvfs3QsEZwtV+jLenzyOKn0ZsixTUWfmYFkNG/L0zM85xkdrD7NsRwn7S2uw2uwA7N69m9GjRxMeHg5AeHg4Y8aMYc+ePW15K4J2jN0us6e4mrmbjrJkWxGfvDud/J05rJz7gbuN57vqi3kfv8OmTZta9I0+Zw1vZb2FFbtK+d/GI2w7WsnxspImH7YnJ/thPPfHxHdAFRGJ1WwiJEyF1WzCrFDxyxELi3ILOXS8DnD07tVqNbNmzcJut7t79hqNpjVvWyBoE1bO/YBDTXzgTFYbB8tq+HFnCZ+sL+C3feUoo+KJiYnBZDKhVqsxmUzExMSImSDBKWOzy2w7Vsn/Nh7h592lPJDVg8nDerBh6TxkWWbD0nlMHtaDZ27N9PmuAjxzayaTh/XgxwVzkGW5Rd/oc87wmq121h08zpyNR9hbUu0OJPf3sH3RsG1DQ+zav+yTN3l/8jgqSou4YvhoHnlrLlcMH011xXEAjhnq+X5bEd/kFvLrxi2idy8463B9rE7lA2ey2th2rJJ5uqNs2V/A2L/8ldWrVzN+/HhKSkoAx5JMVlYWxcXFbXFbgnaCLDtGuF9sPMJv+8qpNloA3zORkqTAajE3elefvqUv708ex8S3viJz4E2o1I5jWvKNbtEar1arHQVcqtPpnmrJeQLFrqIq1h/Uu52lqvRlvPznQV5ZXDYsnceGpfMICQ3j2a/Weh3/zK2ZXuu1rraSpABkXr5nMLLd7t6fk70EAElS8May7QCMnDi1kVyFhnoKDaC3hGAUvXtBG3Em9PmZz1bww0f/ZseGbCwmI6EqNVaL2f2BA2+dm7Yk1+v4O596E4CD9nCefvF1UmMdHz3hcCU4GYeP17HuoJ7jtaZG+3zNRCakdaH82GEUCiV2u41QlZqeVw5BoVCyedX3bFw+n7DwSEzGRMLC1JhMeaf9jW7RiFen0y0ENC05RyCoNlr5bmsR2XvKvDyUV879AFmWSUzr4tXzyRx0E898tsLrHFX6MlIzLiYuuWOj88uyHVmWvYxuw/2uXn5TFBWXMGD4ncyc+z1/u+9+d+9eIAgEZ0qfD2zd6PWB6zfo5kajDV8658nRinoWbj5GdEysWJIRNElVvYWl24tZsq3IbXR9LQ/WGPRcMXw0st3x/S4/dhgAu91hIywmI7m/rCQnuwZZ/hcblo5l4/K3gSPccEO21wzMqXLWezXvLKxi7YHjmG12qvRlzJn2BEf2bvMauZYXFrj/bTEZObB1U6PzeI9mlYAasAI2QsJCsJodvfmYhCSOFxYgSQpk2dHe1XO6+f7JTco6buoMAAyAduwTXHtBwunfuEAQBKyc+wHVFeUkn9eVsU++we/LFlBdUU6UJt7LGKsjoprlVfr0pz+xbPZ/2LE+G5OxnvDwcEaMGMFrr70WgLsRBDNWm53cI5XoDhuwNhgEeS5rjJz4HHDie5s19kGvWRmki+nQaSJ1VZdRW9kTCPc6V0SUlcTEOGbMmHHasgbE8JaWlp60TWuHB5isdtbmV3G44kQ4wvLPZ3JoZw69rh2G3WZjz6ZfsZpNSJKCuJSO3DT+KX78dDplR/JZ/vlMrv7TP3lx9JPYbZcB3wAdgU5AMp6TBVYzQBEWUx76osMkny+jCt9Kwe73kaR6rGYTylA1SlUEtbW1zZK/thbmHzeQkRDOgPOiUYUogjKEQsjUPIJRptPlZPpcX1fHs3dd7tW5LTl8kLce/hPKkFBS0i/iwNaN9L5uOANuvov1P8xly+ofufK2cUTHNZ2LT6mKgJAwzCYjoWFhGI1GQkJCUCgUTcoVjM9fyNQ8miNTaY2ZtflVGOq943Bfvutqn8uDnkuJSlUEFmsXLKbHgTtB7kvZEc+z5CIp1iLbf6fnVXFMfulBMhNlSks1p31PATG8zU0m3VqJsEurTfyys4RKs5LIyMhGa7Pbflvm/rcyNAybxUxFyTG+eHEicCHwPLoVt6Nb0Qv4wccVbEAtkiIU2S4BoUAqkIosX0XJIYC7QZpGcpdKkLLZse5tho2ra3asmIsSI/ycb2LQRYloNK33jFoTIVPzCEaZToem7mPhQgXrVym4dOh+9ur2UlFqQrYrkBRWYhM1KJUWju1fC+zHZj6f1C79iYz8HmNNFeu/+5yRE59zz0z92U9spammiiuGj+aK4Xew6ceFHCvTN+vZBuPzFzI1D38yWW12Nh6qIPeIGVmhIjLSuzSfLx8D1+yjKjSK7euj+X1ZHHnb3/c4ygB8DywHsoFSHJOXEnnbY5Gt96LRJLddkQSnM0aGVqvtp9PpNrfkXK3FjsIqVu8vx2Y/4TDl6+FHRMfSre8VWExGtq1Zjyr8cSTpHuprurmPkyQLKedXUF74FRbTGuAIUECvq3sRFathr24NFaWFyChBTgXSgQygJ3ANyP0oPhwP3AHcwX8eLOC6kWr6DzGg6eDdM2vqY1NrtvLD9mLSwu3cEp9IWMg554wuCACtoc+vvKJk926XH0QX93bZDgb3gHQIADvWO/7AU8ClbFi6ig1LL0OS6gHZa1rQE9cUIcBtDz0LwNLtxWRd3AFV6OklvBG0P0qqTPy8u5SKOnOjfZ7f04ZOVApFZ9Z+14NNP2morXKYQElRR48B1VyaZSIl/RA/fjavkbF2OVl9/dEMBk9vWZ6FFhlepzPGwhZJ0ErY7TK/7S9nR2FVo32+PNiq9GXkrNwIPArMw1gb62xdiaRYjGyfS+bgOFRqmaL8eV5rVL8vm+9ev3VgxWGUjwCr3T/UtjUbsFkvBW4HRlFf3YUfP4cfP48jc1At147U06mbYyrc1xqEC9dLNPLRF6mX1AztnkSH6OYXXRYImkNr6POYMXZ+2VKOFFrDusVv0/miVHpedR3b166mMC8PY50EchegG3AB0B3IdP6ZBNQgywuBz9iwdD4bls5DGRLKa99vafK6eeW1zMsx0z9R5tEH7xOZ385iZFlm85FKfs/TY29QV9j1rYxPTnOHdO7etJp+Q26l11V/54ePJHKyrwMc38/UDCPqiDnk75hEtGYY3a94Dkj0shcOJ6ul7mv8uGAOXRbMaVEGtbNi6GS02FiyraiR0fX0ZHN5sD3y1lz6Z91PbOJ8JKkAeBqIRZLWEBF9H5mD7uexmVFIip/ZnP21O+TBtUb1x4pFPDtnVSOvzIS0LkiS5OUs8q/PF5E5MIRQ1WNAKsqQ29F0WIOkCCH3Fw0zHslg8vAyJg97yGecowuXUf5twWwq6sws2HyMrUcrA/NwBYJT4KmnbNw+vhBj7b+wmN8lqfNWrr5Vw9/fuJW+11WCPAt4FhgD9Cepcz9gCJLiNWAjEAXcC/wKHAAeo/fVI5p17cp6C/94+nnWrROZ385WakxWFm8tYv3B442MLjicYPN36MjJXoIsy+RkL6GuOo6clXfw+Us3UF54AxBKfMpGlMprKMoLJ3/HeKDK69vrZS+GjCA2Mdn9vVep1dx2220tyrHQ7r2aK+rM/LCtGEO9pdE+z1HkuKkzsFlh/Q/x7NwwgvoaJWAHvkEZMgO7bS19rr2TkRNfAODRdxYw61/3YzbWYTWbvNYGfI2g7Tabe93J5bnp3S4Mm2UxXXvbKTnyBp26vc+WX9Mw1l0H/Iok/YQsP02oao/7Og3XpnUrvkG34ht3vOPRinoxvSYIKmJjYzGZTsRNejqzXHzptcQld6TzhT0BOLpvB/W15Qy4qQ9XDL+Q35e9zc7fj1B1fCgwDseyzVvk/lpC7q/TUYZ8zGvfr/V53Ya6IvI6n30c0dfx065Sn0VrGv7+Djrg6OQ9CIQh283AF8Ab6Iv3AbgHSw3Xfz2X+8Y88SrfzHyBjcsXEBKmwmwyERUV1aIZlXZteAsNRn7YXtyofJPvRBdFwCzA8bAuyKxBtj9Oh041XDH8MX5f1pHqinL3MRuXz6euqgLAZ8iDq0fkaWhdyTE8k2Q0bOdaF+58wWSm/O//WLsknpVzo7BZhwJDsZi+RKFYS0x8Bya+9RUfTRmPqb7WbeR7XZXlDktyTa8N65Espp4FQcEdd9zBnDlzGiUhaPgx88XIiVOprngUVXgO1YaV7N8chyw/A1wOvEGY+hV+mZ/Prk33cM+/3vA6X59rbyQne4nXdQcMGson7711Zm9YcEbwLBYjSRJ/HKpgY36F3zz6z3y2gm/fe4mdG1YhyyHAY8BUIBqwIym+JDruHeqqt3oNpFzrticLa/P8ju/59TvKy0+eWrgp2q3hPVhWy0+7SrycqFx4O1MpUCj/jd32EACJHU3ccn8Jl1xegyQ97j7GZSx99ZxcIUeuVI/g7eDhKxtVw3b+sl6FhIZxQebNVJb/jdIjN2Kz3k3ur7fRsWs1pUcWUlupBxzG32YxN3oxKustLNxcyHUXJtA9NebkD04gOAM0rCTkmYSguTG6cEJfvpn5ArAAZeiP2CwDiY6bQXXFRSz79ELgPb597zvGTR3bSK88r2tRqllfJHNTYmBKvQlaD1dmshdfepmr7noEvaVpUxUT34Gyo4eQ5UHAu8Alzj1LgSkgbyM8MoNqvdnLyFZXHG80gPKF5/d+4ICXyUxsfllLX7RLw7vtaCWr9x/32/txTfFaTFrgM+y2rkiSlayxFQweXUZIqM/DgKbdz081FKip84aEqVCFRzLh1Y9Iy7gYgIqSPBZ/mMzODTEsmRWJw9szD9jo0/i7sNrtZO8po7TazDXdEk6p1qlA0Brs3r2bp59+miVLllBfX4+kUHJh/yuJio33+c6ejIYzRRuWdsfhpPgWoGXHei2Th81CoYwkc+D1br1qeN2SaiMLNhcyIC2E4AuSETSkYQfuk9kf88nsj0GSmDrnF5/fYEfnSwPMwOE7ALAXeIQBNyVwxfDn+H3ZAnZsyHa/U6sX/Y+ta35i0gffuc/Z1ADKE4XU8u9ru3Ou2pCn57f95W6jW6Uv453H7mLm43e504HZ7ZC3fTjwG9CViOgCMnpN5IY/N210oYEHdGjYKWXVafZ5nT2u2ko9vy+b724Tl2zh3ueOMnrSdsLUxUAfYB0K5X/oc+1t/HPWD149L5fz2LGDe3h/8jjWbT/At1sKqQ1QMWdB8ygqKuKOO+44qxP6p6amuisJhYSFgWwnPimNMZNe9Xpnm8LTGXLc1BmMnDiVtIyLGTlxKlO/XEXmQCMhYZnANMAMTCAm/hhmYz+3Xvm6brXRwrLdekqqjE1dXhAEuEpBqtWObFGhKjXxqZ2RwGcBG1mGWybkogzZh8Po1qFQTqXPdU8x9cspXu/Qc1/+6v5/mFpNfXVlk0Vx1KFKzk+I5Ir0eIb3TGG0thP3X3U+11/S8i5cuxnxyrLMr/tOhAt5uo0f2euobbvskzcpO1ZHSOg8So+MQ5JkBt1ZzvV31xAS+o9mX8vV2+4z6Ga2/vID1RXlKCSJKFUIkSolUaoQQpUKR35mHCWnas02quot1JpsfkfiNQY9SBJWs2/nE1eCeG2Wkn05j5D7qxZ4ArttEvs230mvawtYNGMcIx58hsX/neZ2mf/qjScpPZLHyrkfEDPxOebnHGNYj2RSnAnlBW3LtGnT3PU7z+aE/qWlpYwfP57YXoNYs/w7v9N2/mgqpM7VcbVZ9ISEvYDV/CURMd9jKEvHUPY8Kef34s7HQ/hjxXyf1zVZ7Xy7pYjhPZMJNVW61w9FyFFwkZKSglFSYTI5OkkWkxF9kSONlOtbqQwJpcvFfRjx4Dv8/GUPdm5wLbGtQBk6Ebv1ABFRdzYxOva95Pfa91voqFGTnhDJeQnhxEWEnbH7bBeG12aXWbG7lAOlNe5trtzJ+TtOtMvJLgUWAPEolBWMeeIImQNP/RbHTZ1BlCqEKMnM326/nrRYNQmRYUjNmGKw2WUMdRaOVdZTaDBSaDBSa7a6z1ulL/M7le3Jod0bgIV07VPEsQNPUl/TmTmvpIF8I29PvBOwue+9pOAg4P0SvfHDFgZf1IGLU6JP+f4FrUPDabOz3dN23jxH6N2Hq3Z7TdtJSCgUoJQkJElqpjNk42pFDaefK48/RGXZAxw7eCvFh+7guw/q+MuUHsT6WdO12Ox8v62YjXP+LSobBSFWm52Ve8o4WFDIFcNH0/OqLL5990WOFx9Fttu9HKJysmHm4z2wWWIIU1tJ7PgOXS76gwE3vdbkWq2vpcRLrxvK/738Mpddch7qAEWIBL3hNVvtLNtZws4Dh5kz7QkK9mzFZm0cOgT3AB/hCIxeh902mq/eKCJz4PZmXys8VEm3pCguTIoiNVZFWVkZSUmxJz/QA6VCIiEqjISoMHp3dBxbVm1if2kN+0procGUc8OiDA0/Qge3zsDhjf0KyI8DLwGDgT8DhV7X9vR6ttllft5dSlmNmasy4lGIdd+A41r3/HbxYsxG4zmT0L9rgprkxAQSIsOIjwxt9DGz2uxUm6zUGK1U1Fv48Pu1fPDGi2xdt7LJzqjnlPUfKxY59eRH4CrgSwr2nMfL98Ty4Ov1dO1d10guEXIUvBw8fJSRo8dy5+Q3vH7nbn2v4LgzjMeRyCIb+A/wMDYLwArMxvE8/u5y4Eag6bVazyW/MJUKi8VMz/OTubb3BWfy9hoR1Gu89WYbi7cWcURf556GiktOA0ChcCmzBLwK/A+H0Z0JDASONbscX1K0imE9kvnblecx8MJE0jTqZo1um0uHaBVXdk1g3BWd+VNmR6ivZMBNjuDs5C7dqNaXudcafBVo7nnllfS8cikK5U1AMTAI2Arc7L6GpFD49HrecsTA99uLMfmIfROcWVJTUzErVFhMDiU/V2os906NpFfHGNI0ap8jiBClgriIMDrHR9C7Yyx3D+xNZtdUrBYzYSr/YR2ea8DeerKOkLCriI7LBZKZ9cx5/PZNPK78ClX6Mj597gF3IfMTiRB8FzIvKioiKyvrrF6TDybKa0z8ffJU9m39o9Gaa41BT/8bRvLIW3PpfvkjKJRbgIdxrPH/E4exLWjWdx4gWhWCylLDffePZ+2aNUxoQWm/lhC0I94aoyNDyQNZPbxL+HnVTIzAERA9EkfaxonAh+62rp7zdX/6G+9PHtcoD3Ln+Aj6d9HQOc677NOZQpIk0jRqfl2+mNjYWNb/8LV7n+f0mvb627ycsKLjEp31fheAlAnypzheuO8JU7+LJuldouM0JHVK9znFUqCvY37OMW7unXJG1y0E3hQajOw/7Jg2++eEe1i5bIn4mPuhtLSUCePHc9999zF79mwOHTlGj7QY9hRXu0MGG64Be84c2SyF9BjwBuFRM/llfiI/fJxC0SE1o/5RyMq5H1Cwewsbl8/3OsZsMmKUVI06Qq5QFjEVfeZpKuHKtCW5jJs6g9raWvb8nsb+LbdhtymA3Thm/BzpxJtTdjVGHcpl58dxYXIU936/yL29JaX9WkJQGt6KOjOLtxRRbbK65+S3r1/p5ZSElIQ64leMtZegUFbT5aKXSE2vpKbyBravXYEkKdw9543L53spbHxkGNd2S6BzfESb3eOePXu8wi88X55v33vZZ2xZ/xtGUlF8lOqKf2CzTOR48cOYjROJivkrY588Rky8f29mQ72FBTmF3NgjiS5teN/nCrUmKz/uLHFPm12UEc4tQ4e0sVTBi2t9GLw/hpedH0dahwTMPhwSJUlqpCd/eqSUzhfW89W/08hZqSFnZS6wApDd6V8lScEjb83l92ULOHD4GBvy9AzIiD/n1uTbmh2FVTz1yQq+/+gNn2GWABazxPcfnkfOzw5PYk2HFVzYby71tYlsX4vXd96XM1V4qJLLzo+jR1pMUIVZBp3hLaky8v22Yq+0YAe2bsRmMSMpFI5i9FIGyMsx1l5IfIqZ+14qIanTvQB8/tKjDLhpDFcMv4O3H7nDrWxwQmGDQZE8wy/UajVms5nUxHhiE5L8Jueora0lMjLS/f+8HUf4clpH8nZE8vbEdP78zDEyejVe23JhstqY88t2lrz9FN8t+Pqsn/JsK2x2meU7S9xOdYLTJ0oVwt69e5j85FN8//0STEbfsfWeetLrqmpyf32W7WsfAQaBtAHk4YSqCr2Ocx2jO1yBXZYbxSKfK2vygUaWZTbk6ckpMBAdn9gozNJqNvH7svlcesNj/PepGMzGJEJC7dz292Iuu7ETkvSk13feFZObNfbv7vdBqZDo10VDv86aoKzmFlQS5ZfXsii3iHqLzb2es+yTN6muKCepSwYZvS4lPvUWlMqNwIWoIw/w8PR8kjqdmIr2jP979gvvYgbh4Y41ndWrVwfFGo4r/GL16tWMHz+eUHMVY7QdSYlpXhhQRs86Hp2ZR9fetVRXhPLhM+exdnEcsuy9HubJz3PfZ6tuIw9Nnuoz65eg5fyyt4yiShEz2lqkpqYSp4nFYjajUqmbHOE8c2smk4f1YPva13CkmtwG8kXA71hMPf0et7nAwP6aUK/O8LmyJh9IrDY7P+0qJafA4N7mO8yykHcezcBs7E5YeAkPv3mIy4cZcLneeH7nG8bknp8Qyd2XdeaK9PigNLoQRCPevWV1bC2rQXZ6RDQMFyo5fJCSw2nAHCCGbn1qGTfVijrS7vecMfEdiIyKxmo2eynSxx9/HBRrOP6m10b1S2N3cTXrDuoxnsQpKibexvhXD7P8syR+W5jI4v+mcvRAOCEh07ym1xt6dC5b8AWRC75ApVJTWWlo9Xs7V9EdNrC7uLqtxTjrcHVSXWvA+w8dJTxU2Shhvne4yFGQriU0bDEW03UolKspyn/e7zW2Hq1kZ94R7r//fu6//35mz57d5p3zswmjxcbSHSUUGuq9tnuGWW5fn43VfC/wDhAK/IS5/i5mPFLRKLzMXxhaMMxonoygMLxHK+pZl19FZGSknyoTAEORFN8h29X0vqaSu54oJCSs6RHb+QmRaKQ6JkxwKOyAAQOYNWuWe3+wruFIkkT31BjSEyJZfaCcfSU1TbZXKuHm+0rp1NXIl6/Hk7NSA/wD+NUr6Dxz4E2N4ofHPPwMZdUmUWShFThYVsvvefq2FuOsxFcntc5sY9XeMvLLa937GmeIq0QVcRdd++jYsymNI/teY+uaY/S5xnfn6PYnptM9NYZeFyW2mePN2UhVvYXvtxWj91G0Hhy/W6g6Bqt5BjDBuXU68BQhYSH0uuqmRs5TDWNy29PSQFCMw03WE6PWhuE0AEijgCXIdjWX3VjB3U8da9LoSpLElV0TuLlXMgsXzGfGjBn07t2bAwcOMHr0aMLDHV7MrqnnltRVPJOEhykZ2j2ZYT2SCW9GYHffgVWMf3UXCuUR4FLgD5ShA8kcdBP/+vznRmsp6ogoFJFxLNxceFLjLmia0moTP+8uRUZuNM1fWR7Cnt3Cm7y1iQhTcnOvFAZd1IFQ5YlPmWct1eQu3aipKCIy9lGuHnEcm1XBl9M6sWGZxu95dxVVsXJPGXaxFNMqlFabWLC50K/RBairVrBzw2RgAspQGzEJ/wSeICQsxGeYJDiMdXRMTKMZzfawNBAUhtcTzx6rJCmAe0D+GggjIXUho/5RhKIJGxQZFsLtfVPp30XTKBa3oUNTe/mhuiVFMfayzpwX1/Ta7zO3ZvLRvy7FbusLZAMp2Cw/kftLCjHxHbw+SFcMH+1OXm+12/lpVwnrDx53T/ULmo+hztGbt9gcHUjPsBezUeLTFzozbmwqq1cHj1fl2UTPtBhGazuSEOmYtRk3dQZ/rFjEWw//iZKCAwDk/LyQtYsTUSieQ5YlFs1M47dF8X7Puae4mp/3lArj20IOHa9jUW4hdX4cDav0Zcz4x1Rm/KMLtZWZRMdZeOjfh+ly0XYG3DTGUYj+hpFUVxxv1KHtkRZDLLVMmHDCT6YtYnJPh6CYavakSl/GtjUr6DfkVuKSnmPl3L6AguvvLuP6u7vTVF6LlBg1w3smE6nyf1sN14rayxpORJiSIRdoKLGoWXvgOFZ747Vt19TLltXLke1DkRRvI9snAh+z+L/H+fO/ZqB0dlp8ZXfJKTBQVmPmxu5JqAKUOq29U2tyxJvXma0+1pzms2HpncAldO5ioXt38RE/U8RFhHFH/zRW7S1jX0mNz2pgjqxuQ9ixvohv30vlh49SMNcryBpb7vO7sq+kBpsdhnZPCqpQlPbCjsIqft1X7rMz78q1Hxp2A0f3fwREkJph5G//V4Cmg9UrsuPm8U8RGRnJondf5NDOHFZ9/V9mvf8eXTtEMnj+iSIz7WlpIOhGvCvnfkB9TSX64ltYObcfoGD430q44c9lTRrdi5Kjub1vapNGFxxrRa6p5xkzZnitHbUHenV09O7jIxtPXbpmC5BlQsJCQP4HGb0+RBkis3ZxArOndqGu2vGT+/N6LtDXMS/nGOU1pkbnF3hjttpZsq2YKqMjhWnDZRKF8nVgJKoIK+/+t4TExDYU9hwgVKlgaPdkrr0gEU1CUoMEGyemK6+8uYLR/zyGpJBZMSeJb9+L5r0nGusCwMGyGpbtKMFq8+/EKfDGFS70y94yvzNoL98zmPwdF7Nv8/NAPLCYorxEXr+vV+O2d13N5GE92LB0HrIss+77r+nROQGNRnNG7+NMctqGV6vVarRa7SitVpul1WqzWiqI98OdSP6O+xwCKiYz6A7/9TwlJAZkxHND9yRClEHXjzgjxEeGMbp/R3qkeRe+95wtcE0nh6q+JKnz/UTGmNmfG8XMx9MpPRLmNR3qOtZliCudyTb2inVfv1hsdn7YXuzVQfFcJlEo78dumwxYGfXoTtK7+sovHjy0tj63JX06xTKiTyp1lSeWVlzTlS6011dy91PHUChlNiztzKGd4/j5y//6PN+h47X8sL0Ys1UY35PhChfSHa7wud8V7iXbp+DIOhiGo8bySDIHDeSZz1Y0OubR978jc+BNqNTtwzenObTEUk0ANut0upXAAy0VxPVwFconcbiSQ6cLPmDKFzf7PUapkBjaIwnteXEtvXy7I0SpYPBFHbj+kiRCFI6f0TVbEKZSu2tQxienUXzoUzJ6/Z0w9T7Kj6n494RENiw1OHqmS+cxeVgPXv7zYC9DbLXbWbGrhN/2l4t43wZYbHaW7ijhWIOwCHA49lxy+f9xInXpRPK2TQukeKdLq+pzW9MpLpyfvl/EXye/SFrGxdw8/iluf/hZr1mer//TFbvtFsAETOT3ZQOZPKynz5y/RyrqWbKtSBjfJqg32/huaxH7S/132J/8+Gfik38GXgRswEQkxRNIkuw3xjqxQxIXdU7CYm5fvjlN0RLDeyngip3QtFSQ6LhEjhePwW573bnlQTpf+JvfAvRhSgW39E7lgqSoll66XXNxSjRP39rXayrGZUw9t21f9wlmYyawEMfPtRR43H0eWbZ7Hzu8p2MEfbSSRbmFVNUH94gtULiM7hG9d4Yw14zB1SNeZdfvjzhzyr4FfMiGpfPofVFGsE+Ntao+BwNxEWHc0a8jaRrHSKnhLI9jaQCUoaOAeuAB4lN+4qlPGo+6AIoqjU5HIVFwpCEVdWYWbD7WZOKY+loF89/qh74kC6jFkWP/PXpdleXl7OlJcrSaET0SqK/SeyUbai9OVP4IiHNVaWlpk/v1FUZ+XZhIwZ5LATvX3bGJ2qpyKsqOU1tb26h9eKiSay7UoLJUU1p65pIVBFNsL/iXZ+3atbz40kv89NMKzE5Hkm6ZA0CGPZt+9WhZB9wJPO/88ybQG2XIP5AUZvd6WExCEvrioyz/fCY3T3iKg7W1zC6r4Jr0GLo08KwOtmcEZ04mq01m5QEDhZWN17+Xfz6T/B3H+OT584E4ouPXUVc9BZvFUa5x0OAhTHv5xZPqQnvgZPcQbO/E3wZe4jcRf59BN2O3fotCeRt22yL0xdez5MPj3PZIntsR0ZPa2lo+W13F0IviiFa1zAEx2J4TnJ5MhZUmfjlY6RUW6qK6opyFb03hhnHTWfze5ZQWRKAM1XNh/9e4btSd6H4OpabiOCP/+k8Ar+99t8RwruwUQk1VJTNnznRvnzJlCnDy9/BM0tLfriWG9w8cq+IG5x+/JCUlNXmiN9+38evXEUiSzKhHi7lsaAyOqYjGxKhDGdEnFU1E6OnIfMqcTPZA40uepKQkUpKT3WXVLGYzcYnJbscGSVIgyw6lCFWpiIj+gIRUNQV7J2E134vNehFwO1CC1WxCX3QEAN2Kb9Ct+MadMWZDoYU6RQRXdU3w8vIMtmcErS+T2Wpn6Y5ijpRWMGfaE+5KVyc8maOA3zAbOwAbqNbfgCQZ3Y49cXFx9OzZs1VlamVaTZ+b2yZQ7Nmzh8cee4wff/qpUa5nz4IkK+b8m90bn2T72gQkKYSxTx5D6eMLaQN+KzAzok8qCVEti88Opufk4lRk2naskmXb8vli2qRG1d+q9GXMemoc1fo0Pn++N2ZjBAlptYRH3cuoRyYRE9+Brj0aT+tLksTVXRPo29lRz7xUIbX759SQlkw1zwJcjhgfnqyxPz78UMHbr0UAdm4Zv5fLhhr8ttWEh/KnzLSAGd32hCtMau2aNYy+56/UGo5TY9Az4KYx9LzK4SvjquTR/fKBFOx5DqtZCxQAAwAdoHW0c64Zh6rUZA66ycvhYevRSubnHENf6z8Y/myj1mRlUW4hRyvq3dOVyz55k/cnj2PiW1/R99pbkaSFQD+QDtDzqplcpO3vFTN9vLxxucYgo1X0ORhJTU0lOjoai7lxvV/PnL8jJw4i+fyJqMKtbFsTyxevdsJq9h1KUWu28k1uIYWGcy8nt6te8bfrd/HbvnJ+nvu+1xQ+OJyoXrp7INX6PsBqzMYE4DeOF3bh2P4fGtXddaEKUTKid4rb6J6tnPaIV6fTGXAoa4sYNszOlGdLqKl+irJjRuA5n+3iI8K4rRnhQucqnmFRn816H32tmaXbizHUW7wqebjKp50ot3gNVvMc4BokxTrOu/h9Du/+p1d2q4br7HkFR7n6oTG8O+tTenY8uxWkos7Mkq3FPHxDT68Y3ZzsJQC8PfEOQsLmIctDgTKQhxGt6c/IiSc8ZEdOnMrQjMDUfD5dWkufg5Xy8nJ3/P7rM95nb/5Rdyypa6S2cu4HFOfPp+eVKezPfYGdG2KY/byFv71QSqiPTHkmq43vthaSdXESFyafO74mL738CuvWrWNdltZRLc6Jawr/BOOAj3DkXP7a+X8zsty47i44BlbnSs3wNrViJ+pfqgATG5bS6McASIhUcVvfVCLCRFKH5hIfGcad/TuyfFep3zKDqohIbJYjKEOHY7P8G9n+IId2PUbK+Rdx5z+NrFv8OVvX/MTlw+5k8X+neX2g8nfk8NLLL/PXx/7F7THxRKvPvg5RUaWRH7YXY7TY6HPtjW5j681MrOZRQDV3PbmDQzu17vrJguBh1qxZ7qnBL2f/lz3F1fz94Uc4tDPHXZDFxfZ104BlwEoObEnk0/9Tce9zRwhTNza+NrvMil2lVBktZ310RcN6xTSI0ZUUSmS7jV5X3UDZsb9QfOhu557XgH8REhoKkqND37B4fae4cIb1SEZ9jiTuadOvpav+5eLFSzAaafRjgMPo3t43lXBhdE8ZVaiSW3ulsC5Pz5Yjhkb7XSkkXSPhw7vrKSl4jOJDw/jho1rikn+kvrqSr954ktIjeY0+UK5e60NhKtbvLqBXWkyjNJ3tle3Hqli9v5ynbunrp2gHwMvAQ4ARuIWv3vgNZUgoXS7uQ5W+zD1TEKMOJSLs3Igxbw+czIA42AoMBLLZn5vMlNvLeWlhB5/V0GQcCSMq660MvDDxrM1y9d2vf/CvZ55m+/oThVZiEpI4XlgAgGy3AaFsXzcGuBuwgfQPkN8nPCqG+poqgEazaT3SYhh4QSKKs/S5+aJNvwaeuZNDQsMaTW26ppeF0T19FAqJa7olkHVx47R33utbU3n8vaE88NohoIi87ZHkrHwOWb6CkoKDyLLsZXTBoUARMXE88tZcfttXzrdbiqhoIhF6e8Bqs7Nydym/7ivDLsuNslFJCiWxHVKAJ4ApgBW4g1DVRjIH3UTf64Z5rXdFhIUwok8qYedIcpf2wO7du72KpYSq1CSkdUGSJJShjmlOhUIJ7ARpEHAUuIZZU85zZ37zxa6iKhblFlJj8p2XuL1ittpZsauU7RUKwsK9C63YbTb6DxnBRdprkKR4YDnwV6COrn1eQ6n8CMBtdAFn3V2J6orjXNk1gcEXdTinjC4EQcrI0tJSRt9zL/dP+8Qrlis2PJQRYnq51bgkNZrb+6YREdb0JEd6j3oee/cQkbE7gI7Ar8DDgMPQuj5QLsWrq6rg92WOfKnHDPXM3XSU9QePt8tEAxV1ZhZuLvSqp9uwzByyHXP9Q8C/nS3uRZKWYTEZyf1lKTnZS7ziof85vLdwBgwyGhZLsZrNyHYbVwwfzT/e/orkLl2x252xuvJu4BogjyN7w3n+ThM1Bv/fpOIqI/N0x84ap6vSahPzdEfZW+LQiYaFVjp2vYQxT7xKeFRvZHkNMAQootdVL/Pga7fzr89/9uq4uhw2n5/7C199/TX9u2ja7N7akjZfmJs3bx4Hy2pZsPEgXZ3rj9HqUG7vm0aUcKRqVVJj1dzZvyNLtxdT1kQu5o5dY+kx4GU2/XgZ8BjwLnAdVvN47DYbSJKz1+qgoaNEToGBvSU1XNU1oV04nciyzJajlfyeV+Gz+ITnlPw3M40U7BkLQFrXNznvYhtXDF/A6kX/48DW36mrrnRPww2/6RZmvPnvRucTtD0Ni6UcKjjG8Eefx2yz06FTOp0u6ElNpZ79m9djtx9yjHzln4A+vPeEifGvFBCf7DupTHFxEVmP3s3MWZ8ypN8FAb2v1kKWZXKPVPLTH3v5n0eokC9/kcO7w9m+7iUgjoTUSs67+FXMpj3An33URzYRGRXNuMF9SI1tutra2UzQWbaIsBBu65N6VjrrBAPR6hD+1C+NlXvKONBEare6qjIG3LSRvO3PUnrkaWT5DlQR1xCf/AIT3xzjVfnF19p8jcnKT7tK2HKkkgEZcXSOjwjE7Z0yhjoL2XvLKPSR/tGF62OT/XUiBXuSnPHmRVw29EbgRgDGPPEq38x8gY3LFzg/MGaSE+PadVq7sxnPKABXVZuSKiOLtxa7f+9vZr7gMfItAK4FfqL8WCbT7o1g0gdGUs53dEA9PaRXzv2A/J05vPzKy9Q+/xpDLurQrqIxKurMZO8po6jSyApnqNCyT95EX1LYKFZ386oYFrydhs2i4ILMGu6ZUkR45N+9zufZcd284hsibNXntNGFIDO8jjSQKWJq7gwTqlRwY/ck/ogMY1N+BTKNnUs8e7Zlx4qZ82pHCvNSOLTrPbatKXGv9ShDQrGYjCgUSp/pPUuqjXy3tYhOceFckR4fNApnttrJKTCw5Uilz1GuJ7IMK77owMqvOiBJMnc8Xsil11c2aldj0DNg+Ggem/gg2d991W5KTgocJMeoua1vKou3FmG02Kgx6Ok/ZITHyLeMkLAbUUespMbQi/cn2xg9aRu/fTOB+OQ08nfoeOnuge7zuWeCwlRsPVhI1w6RbXdzzUCWZXYU17K3opbJN/fxGT738p8H88ay7djt8NP/klg1z1Fya8DNekY8UIwyhEZhWq5vSUqMmufuuVH47BBEhlepkLipVwpJ0aq2FuWcQJIkLjs/joTIMH7eXYrFZm+kMC46dDQz8a1DfD8rmQ1L41n831SiNJPoNziO+roD7NqQTf6OnCavd7SinoUVx0iNVZPZWUNGYkSbeEDb7DI7C6vYdKiCesvJc+7abLBoZiqbfopDUsjc+Xgh2qzGRhccnZXrLkikd6dYbh10RWuLLggASdGOKIrvthR5jXxl2e7MQlZGjyveoKbyPXZuiOGzFy4COpK/Y3Gjc3nOBC3bUcxFydFc1TU+KEe/hQYjq/eXc6ikmsjISHecf+6vS73aybKdycMuR5LmIMvdUShkbn2wmKtuOVGNyDMn9siJjrwMXTtEcv0lSYQKJ0MgCJyrABQSXJcRi9JoICsrS4wUAkjXDpGM6pdGjDq0URJ5T0LDZEZOLOaeKUcIj7JRY7iMnOwX2bXB4RmqLznG5GE9ePqWvj7r/LooqjSybEcxczYeIfeIgdoAeYAaLTZ0hw18sfEI32/czfTH/uxXRhdmo8TnL3Zm009xhKrs3Dv1iF+jC3B5ejy9O53dCUXOBRKjVIzMPOGI2NChqLaqmD1/pAAfAxHAIhy+EN5YTEa2/Lrc3Yn9Y3cel109kJWb92Fv44pfruxTBw8f46ddJXyTe8zL78Nd2xtHxjsXIaG9UUXsQpZHoI6wkNRlIr2u2gecKPnXsFjLlBH9GNYjWRhdD4LiSZwXH8H58WqmTZvG+vXrefXVV9tapHOKbp2SeWDQhY0Uxld5tN5XVzPpg4Ok96gAkoDvgU8ICUvxGU7jD0O9hbUHjvPphgIWby1ib0kNpmaMQE8Fu13mmKGe7D1lfLq+gA15x6k2WprsYLio0it574mO7N4UTXiUhQnTDtP9Cv9r4n07a7js/LM7gcK5RHxkGCMzU4lShXiF3WWNfZDaKgOPvD2Hvtd9h0L5PI7P6FvATMAxjRoSGkZcckfUkVHuDt7KuR+QtyOHV155ha90RymoMPotFH+medGZfWrCpCns81N3u2HKWbgdq2UNprrOJJ9n5OLLJlNy+IMG1Z68PZiH3jqSffv2njXx/a1FUMx5xMfHeQW0z5o1i1mzZqFWq4OygsfZhiuRyZIlS6ivr/fpLOVJbKKVB98o4v1J33N4z53AX7Gah5L7y98Bx1qQr5RwvpBlmQJ9HQX6OiRJIjlaRae4cDrHhZMYFXZKmWxkWabKaCX/uJEt5aUc1tdh9DDmJwoa0KSMh3aF88UrnajShwKHuKj/DM6/ZLzf6/ZIi+GabgnNllPQPoiLCGNkZhrfbimi2ujwYHZ12jYun486MhLZ/hJwAPgEmAhkAPdgs1YQplJjKC30m3gmJDSMj1bt4tLzNHTtEBkQ4xQbG+tVqWndD1+z7oevfeqpa6r9sxf+SVrX7yg8OMK5Zz4lh/9GyeHaRvejvf42Dw9mM+elJAoHQx8ExYh39+7djBgxwh3QHh4ezpgxY9izZ08bS3Zu4CuuMSIy2m8tZACFAqLjv6Dn1U+SmnEcSAMWIynmA6k+Cyw0xFXD1jUikGWZ4iojusMVfLulkI/WHmL2usN8t6WQVXvLWPL7Li6/eiCrcvez5Uglmw5VsObAcVbuKWVBzjFmrTnE/34v4JeDBvaWVHsZXfDdI/eUUZZh3fdxvDepo9Po/gZczpbf3vY7A3BRcjSDLkw8lcctaEfEOguz/GtEv0bTqL8vm48MPP7en+gxYBqSdBwYDmxElnv6TTzjeu8eff87ymtMLN9ZwpebjpJTYDgjyTdsdpn9pTUs2VbEU580rQPVFeVeOqkvCaWmchGFB0cgSTLxKTN4dGYImQMH+jxHjUHPVbeMYfGPK5kwof3XzT1TBIXhdVUPcX34TSYTMTExoqcUQFxxjatXr2bChPFEyzXEqL29yxsaynFTZzDqsXt47J0SRjxYjEJpRLbfAezFYnqYMLVjvdPXmm+Vvoy3H7nD55Sv53XqzFaOVNSzs7CKt/79OttyNvLSyy+z5kA5G/MdqTB3F1VTXGXE7HQQ+/S5B3yu3/qKKVRHOOKMZ/7zYf73SgLfvZ8KhAFvIklDgVK/nYiuHaLIuriDmEY7y4lWh7Bl+04uz7rFK4MZQP/Bt5CWcTH3PncnT3+qp2PXeqAbCsUmlCF/AbwzY3m+d9Fxjg5blb6MVx6+ix//2MNnGwr4bmsRu4qqMNT5jhNuDvVmG/tKasjeU8Yn6w/z484SDh+vIzo+0acOuDrZvy2Y7dbJbWujeXtiBod3RxCbaKH7FS9QUfI4m36c7/ccj037gIWfzWLIlZcyY8YMr7AtwQmCYqoZvKuHzJ49WzhYBRhfcY0mi40Vu8s4dNwxpeTLWxFAoYSrR+jZvekJyo49SkXJpcB/2PTTHo4XfU3+jhPHVOnLePnPg7zWthpO+XqW3tOXFHJk77ZmTRG7ZCzYvaWRjC4a5qeurihnwdsbKNj9BZAM1AL3AfPcKXwtJmOjKk0ZiZEM7Z50zqW6O1fpel4nundJYqPJsSQmO+N7c7KXkJO9xP0uPjz9EN+8m0rOSg3YP0dSXInF9E/sNluj985FQ706oq/jiL4OgChVCKmxahKjwohShbj/SBLYZRlZBqtdpqregsH5p6zazPEas88wQfDWgdWL/sfWNT/xx4pvPXRMw4alQ9mwtDMAkvQ9leX3UlmuB3BXIJIkBY+8Ndd9P2macG7qee4UOmgJ0ple3H/66aflF1544aTtSktLg67YcbDJ1BbyyLJMrEaD2dQ401VIaBjPfrWWyEjv+MQ9ukhmTzUBFzq3/Ao8A/zu9zqSJKEMCfVbkCAyNh5TfW2jyiYuY9hw/dZTRn9rzE/fMhCb9VXgfueWDcC9gMNLU1IoubD/lUTFxmOqq3WveXVLimLoJadmdJvz2z3//PO89tprQW3Jm6PPwaY30DoyjR49msioaHbmHWH7xrXY7Taf76Isw1sPr6T40EPIcgjqyEI6dn2FB19/2Ot8p/POtjaL3n2R35fNp9/gW7DbbGxdG4rd+i7QEUlh5oa7C7h06DGWfuw7YY7rni9OiWbwRR3OSIGI9vo+NaXPQTHVLAheJEli75493DpylN91IU+euTWT2VPPB3rhKCSgx1HlZQOwGLjM53X6Db7FvQbri9pKvd/pMWi8fhsSpiIyNp5H3v6q0RS51QJrF8cTqjqCw+gagcnA1biMLjhGNfFJaYyZ9Krb6F6YfOpGV3B2MG/ePD6Z/THa7he443p9vYuSBP98P4tHZx4m+Twjxto08ne+y8qvErF5LOE++v53Ta63nkkahv7kZG8k99c7sVu/xZGjfR19rnmKrLFmYhN8L9HExHdAQmJARgLXX9K4CIvAP8LwCk5KamoqyQlx2CxmQv18bFycMIAKYDqQjqN8Xg1wK7ARWAPchuv1Sz6vK6b6Or+xg554VjbxxGv91lnpqrZSz+/L5run8n7+8r9sXRPNfx7oxuL/pmCsVeFwoOqHJL0JOJxgJIWSiy69hv5ZI7yuc3FKNNdfLIzuuU55eRnjx4/nrS+WeBV2aUjHriYefSefy288it0m8dP/knjr4QwObD2RPvXA1o1eBk2hUDJn2hMcO7inyXj4luLS05CwaBwd5H3AX5AUJq770zauGP4eVssOd/uGsczVFccJVSoY1jMZ7XmaMyLj2UzQrPEKghuX89Xf/vY3Xn/nA/bnH/XZztMASpICWa6i9zVrCFUdZef6nhjr7sYxsryaMHURSZ3XER69hHFTnwZOxA7WVOrZvta75x+qUnNR/6uo0pdx+8PPAt7p6WoMekcBhwbrwaABHuP3ZQ/z+7LOzj27gadwxCGfKMkqSQqQ7cQnpXmtEffrouHKjHjhSCXw8ofo2asXusMVftuGhskolE8A5agjvqCkIJUPnz6fPtdWIvM21RXlJJ/XlbFPvsHvyxawV7eGitJCdw1sf74K/rLM+cJX2yhNB2qrBmI1v4lrSSguKYc/Tw2lS7cQ4FmvczQsjhCtDuWmnsl0EJkGTwtheAXNwvNjM3f2fymrNrFiVym1tbWN2vpyYBoz6R98XvcoEdE6IqIfZ8PSVEz1qRzdPwoYxbuT6uh9dRVD7/kvyeeZ+N/Lj3oZYElSYDWbKDt6yOuD5OmYMm7qDKr0ZXz7wTT2bPoDq/kaHCPrOwDHSDpMXczVI4rQFz3Pzo0/YzE5DHpEdCzd+l7Btbf/xcv5RULimgsS6CMyUgl8MCAjntjwEH7ZW469gb9MwzVcY106MAmYwtbVscA7QG9KDr/OWw//yevYkoKDgH9HQn+Ojp64DG58cprbWfF4cTH9Bn3K+h8yKD70BABxydUkpL5JYd7bhKnnAuc1ec+d48IZ2j1Z5FxuAcLwCk6LDtEqRms7snyzmcM1kpcHpa/SYQ23D7u3nAPb6snJjmXHuhgO74rg8C7HFFyUxsSF/RfQoaORovy36D+kK7m/fozdZmn0QXIgsWHpJjYsfRqFoi9Rca9iNXcHThRkkKRsZPkd+g2OYNi9U/lmZrjXFF/3ywe6P2AumZUKicx4O5Pu/RNz5swR4W0Cn3RPjSFaFcLynaWYrB4JW5z5jk84JUn0vHIrmYN+Y9HMKAxl1+BIujGBuOSfSUhdxKFdX3qV3AwJU9Hrqiyu+9PfeH/yOAr2bMVmPRFm1JSHvytxR/4OgARyss8DXufQzosBiE2wMPCOcq4YXsGSD3dzcGslvy2YzejHX/R7r/27aLgiPV4st7SQFhlerVY7CrhUp9M91UryCNoRIUoFl58XTf+wGFbtLaOizrdHsi8USrgws5YLM2spH7udT57/EZt1CPriS6gxpLI5W8XmbIDXObQLQlUzUSorsFhKQTaApAQ5FogG4gFH8hW7HaqOA9hJTT9O6ZGZ2Kxf0evq84mMiXePZH2Nyj2JUoUwrEcyr0590p3G9J133mn5QwtihD6fPp3jIxjVL40fthdTWe8wjP7ixi+59DwuvvQFfl/2EJJiKrJ9FBUlN1FRchPwEEhzQJ6LpCjHZjGjjohi4/L5HNqZQ0JaF8qPHUahUDbyqnZxYqQdCQwG/oLDvyLM2SIfmEZN5VcsnW1m8X9P6K1uxTfoVnzTyJCHKRUMubgD3ZKCv752e6BFhlen0y3UarXXt5YwgvZJmkbNwI4KRo65j9v++bo7MUBzeWP8QGd2H1eO7h7AdUiKPpx/yZ8pPqyivkYJJDr/QOMQxSJiE8vofU1HEjtVsuSDaynKP+beu22Nw1s5JNTx8fE3KgfoqAnnrwO783fTuZXGVOhzy4iPDGN0/45ese++OnhV+jK2rVlBn4FXMfiOEH6Z/wS7N2Viqr8Z6A9yf2A6oWEHMBtXsGHpb0AOEEr5scMA7jrBrhjzyNgOHN5TzYK355CQ+hElBecDAwBXEhwbsAz4lJCwH+l11XXcfP8PAF6jctcI29OQJ0erGdojidhwUa61tRBTzYJW4Y03XmfH5k30+fl/3DhhSpOF5V34i2OEnWQOOp+b77+TmPhDyDKY6hR8Me1F9m3eAXIMjg9JFVDp/FOLxRzLwFHfoVRF0POKLxtM8zWdf9pFZmeHE9WAPd75q8PDwxkxYgSvvfbaqT8cwTmDKlTJzb2S2XSogj8OGXx28Ba9+yL1NZWEhqmZ+dhdHjqgAm7GMUK9EbPxQhyOTxM9rlAMHAUsqMJjsZit6FbGs+mnjtisEnC5R1sbjtj579B0yMZQpkOSFNgssldUgueo3DXCdoUKZXaJ5Yr0eBEq1Mqc1PA6p5+80Ol0C0/lIqWlpSdtE4yjiGCTKdjkAejWrRtm8wnj+cWns/ni09mEhal4acF6d83b6opyFr41hVGPv+IeET/6/nf89PkM9mz61b2uJUkKZGSUoWqUqggv562xzzxOdUW51zEhYSouuXwQCoWCbauXs/zzmWT95VHU6ggUzuk9V3iRr3O6iAxTcnV6DB1jbJSXl6FUKgkNDcVoNKJSqTAajYSEhKBQKJr1PjckWH67QOhzsNyrJ4GWKT0SQjqGsiavEqPVEab28l1Xe3U0dSu+AZzJY5zvaEjYD1xyeS2DxqioKj+PFf/bQuHBeOBiHPnQU5x/wOTs29rdy8pHgf3ADiAb+JVumT2IjI3j4NZDaIf+Ce31t6P7+VsqykrcemAoL6X/DSPRXn87G5cvpKKsBCxGrkqPoWO0jePlZyakqbmcje/TSQ3vqSqlL5qbdSTYspNA8MkUbPKsW7eO//znPz5HhvGJSeQeMbDlSCU/fTqdgt1bWP/d524npsjISKJiYrFZzM7QIzu9rs5yr8U2zIjV8BjXutn2NT+693uuUV186bWNpvl8nfOS1Giu7ZZIWIh37HB1dTUTJkzwSmPakucfDL9doPQ5GO61IYGWKSkJLjk/jVV7y8gvr23kbOWa1lUolGxe9b17xBkVE8v7j93qYzZIicPodsIRA291/qkHDjv/diApFCDLdEjr3MjruUNqJ+ZMewKbyRE7f98L77n3pZx/IZddkMrVXRMa6UNbcra9T63hXJWh1Wr76XS6zS05l6B9kpyc7FXZqGGBi6GZGV4lHxt6YfpaA2u45tqQhrlmD2z9nbrqSveU8sWXDeS2B5/2im/0dc74iDCu7pbAeQkRjfaB7/zVZzNCn1ufiDAlN/dKYWdhFWsOKHxO61ZXHG+kA409or2XShrui0lIQl90BJAcHdirsrycCT3xF4oUow7lmk5x9L2gA0VFRdxzzz3Cm/8M0WLnKqDFPWhB+8aVXMNXgYuGtX5V6nB6XpnF8PsmAU07OfnD85gxT7zKNzNfYOPyBe4RsCo8ssmkAhFhIVyeHkf3lGgRFuGB0OczR4+0GLrER/BVfaXbyK5ZMpfqiuN+daCpKkIN9/kqwtBQn/zWow5T8ev2Q/TtFIv+uMNQT5s27Zzx5m8LhHOVoMU0NTJsWOvXbDbRNyOFmy+/hG1HKzHUn37pMxcNR80VZb4rW0WrQujZMZbeHWOCahpNcG4QrQ7ht+WLKdDXsXr/cWKTO/tc+nDRVMhbw1mf3X/8RtbYB4mJ7+C3A+trFH3l4Bt5f8Z0unaJA+CCCy7A5FEQ5Vzw5m8LhOEVnHF8jYj7dIqlT6dYiiuN7C2pYV9pTaPC9c2l4YihofNUp7hweneMRW2uZNy4UWL6TNCmdImP4K5Lw1m3y0ZetYJqo+/OZ1OzQZ77wtRq6qsrm8xiBd5xxaFhKqwWMxd1TqJrl47uNmvXrvXrsyFoPYThFZxxmhoRp8SqSYlVc3W3BIqrjBTo6zmir6e02uS3nujJCFFInJ8QyfkJEaQnRBCldrzm//jHM2L6TBAUKBUSFydFcHX3DuwpqSbnsOGUZ3/8Th37KSsYqlSgMFYx7m/38dADE3zWPT+Zz4agdRCGVxAUKBUSHTXhdNSEMyADzFY7x2vNHK81o681U2W0YrTYMFntWKx2kCTClArCQiRUIUriIkJJiAwjPjIMW52BtJRk97k1Go2Xg5eYPhMECwqFRPfUGC5JieaYwciekmoOltZittlPeuzJHLDAkWs8TaPmkpRounaI5MFl37n3+XMYbMpnQ9A6CMMrCErCQhSkxqpJjVWfvHEDSo3eDlMNHbzE9Jkg2JAkiU5x4XSKC+e6C+wU6Os5aqin0GDkeI3Z5+yPv5SUaampdIoLp7PzfFGqU/vMn2ve/G2BMLyCs56GDl5i+kwQzIQqFXTtEEnXDg7HK5PFhr7OQrXRSo3J8UeWQZJgsamaW8f8hbv+fC8/LJiD4XgZf72y6epCgrZHGF7BOYGYPhO0V1ShSlJjlaT6qEy5cum37n/fNmRAAKUStARheAXnBGL6TCAQBAsimFEgEAgEggAiDK9AIBAIBAFEGF6BQCAQCAKIMLwCgUAgEAQQYXgFAoFAIAggwvAKBAKBQBBAhOEVCAQCgSCACMMrEAgEAkEAEYZXIBAIBIIAIgyvQCAQCAQBRBhegUAgEAgCiDC8AoFAIBAEkNMukqDVajVABqAF9DqdbmFrCSUQCAKL0GeBIHC0ZMR7J2DQ6XSzgNdbSR6BQNA2CH0WCALEaY94nQrq6ilvbi2BBAJB4BH6LBAEjpMaXq1WO6rhtgbTUM8A45s6R2lp6UkFMRgMJ20TaIJNpmCTB4RMzSVYZAqEPgfLvXoiZGoeQqbm0VKZTmp4m1rrcSrxNCAe8CtJUlJSs4RpbrtAEmwyBZs8IGRqLsEgU6D0ORjutSFCpuYhZGoeLZGpJc5VWTh6x6Odm+44bSkEAkGbIvRZIAgcLVnjXQn0b0VZBAJBGyH0WSAIHCKOVyAQCASCACIMr0AgEAgEAUQYXoFAIBAIAogwvAKBQCAQBBBheAUCgUAgCCDC8AoEAoFAEECE4RUIBAKBIIAIwysQCAQCQQARhlcgEAgEggAiDK9AIBAIBAFEGF6BQCAQCAKIMLwCgUAgEAQQYXgFAoFAIAggwvAKBAKBQBBAhOEVCAQCgSCACMMrEAgEAkEAEYZXIBAIBIIAIgyvQCAQCAQBRBhegUAgEAgCiDC8AoFAIBAEEGF4BQKBQCAIICEtOVir1WYBGsCg0+lWtopEAoGgTRD6LBAEhtMe8Wq12n6AXqfTLQQeaD2RBAJBoBH6LBAEjtMe8ep0us1arVaj1WqfBKY11fb5558/3csIBIIAIPRZIAgckizLTTbQarWjGm5z9opd+zOAp3Q6neglCwRBjtBngaDtOanh9YezZzxLp9MZtFpthU6ni2td0QQCQaAQ+iwQBI6WOFetBDK0Wq0WeKqV5BEIBG2D0GeBIECc9ohXIBAIBALBqSPieAUCgUAgCCDC8AoEAoFAEEBalEDjdNBqtRogCzAAuAL1/W0PApkyAC0nYhzbVCaP/aOAPJ1OtzkYZNJqtRMAHaDV6XSzgkSmLEAPZATyt3NeexRwqU6ne8pjm19Z2yvBps9Cl1suk9DlRnK1ui63xYh3ArDZKegDzdjeljLdiSOLzyzg9SCRyfWjXwrEB4NMLqVwfjgCbUyakgmnTIF+Tq4QHU2DzW35jp8pgk2fhS63QCahy405E7rcFob3Uhw9F/C+GX/bA4HPa+t0ulk6nS7PqRwB6402JZMTLfBHQKVx4E+mO3B4xGYB/YJBJqdCvK7VahcQ+A+IP9ryHT9TBJs+C11uHkKXW0aL3m+xxts8ngHGt7UQ4Oj9BekU5UqnXM+0tSDgToE4DcdHLdAjHEHwInT55AhdPsO0heH9gxPTBYZmbA8Efq/tnN+fRuCnOPzJpHf2Rq8n8D1SfzLl0AZTQE78yZSl0+kW6nS6N4C8gEvlm7Z8x88UwabPQpdbJpPQ5ebRove7LQzvLCDL+cJ9CG6FaLS9rWVy/v8Z4CMC39PyKZNznUPnbJMQJDLNAvo5tzeZ5zdQMgErnb9hP+DnAMvkkiHDef1geMfPFMGmz0KXWyaT0OUGnAldFgk0BAKBQCAIIGKNVyAQCASCACIMr0AgEAgEAUQYXoFAIBAIAogwvAKBQCAQBBBheAUCgUAgCCDC8AoEAoFAEECE4RUIBAKBIID8PynPwRNvTeuVAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 576x216 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Set into eval mode\n",
    "model.eval()\n",
    "likelihood.eval()\n",
    "\n",
    "# Initialize plots\n",
    "f, (y1_ax, y2_ax) = plt.subplots(1, 2, figsize=(8, 3))\n",
    "\n",
    "# Make predictions\n",
    "with torch.no_grad(), gpytorch.settings.fast_pred_var():\n",
    "    test_x = torch.linspace(0, 1, 51)\n",
    "    predictions = likelihood(model(test_x))\n",
    "    mean = predictions.mean\n",
    "    lower, upper = predictions.confidence_region()\n",
    "    \n",
    "# This contains predictions for both tasks, flattened out\n",
    "# The first half of the predictions is for the first task\n",
    "# The second half is for the second task\n",
    "\n",
    "# Plot training data as black stars\n",
    "y1_ax.plot(train_x.detach().numpy(), train_y[:, 0].detach().numpy(), 'k*')\n",
    "# Predictive mean as blue line\n",
    "y1_ax.plot(test_x.numpy(), mean[:, 0].numpy(), 'b')\n",
    "# Shade in confidence \n",
    "y1_ax.fill_between(test_x.numpy(), lower[:, 0].numpy(), upper[:, 0].numpy(), alpha=0.5)\n",
    "y1_ax.set_ylim([-3, 3])\n",
    "y1_ax.legend(['Observed Data', 'Mean', 'Confidence'])\n",
    "y1_ax.set_title('Observed Values (Likelihood)')\n",
    "\n",
    "# Plot training data as black stars\n",
    "y2_ax.plot(train_x.detach().numpy(), train_y[:, 1].detach().numpy(), 'k*')\n",
    "# Predictive mean as blue line\n",
    "y2_ax.plot(test_x.numpy(), mean[:, 1].numpy(), 'b')\n",
    "# Shade in confidence \n",
    "y2_ax.fill_between(test_x.numpy(), lower[:, 1].numpy(), upper[:, 1].numpy(), alpha=0.5)\n",
    "y2_ax.set_ylim([-3, 3])\n",
    "y2_ax.legend(['Observed Data', 'Mean', 'Confidence'])\n",
    "y2_ax.set_title('Observed Values (Likelihood)')\n",
    "\n",
    "None"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}