tf_prob_pca.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "Probabilistic_PCA.ipynb",
      "version": "0.3.2",
      "provenance": [],
      "collapsed_sections": [],
      "toc_visible": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    }
  },
  "cells": [
    {
      "metadata": {
        "id": "P20Z1RhWhPzX",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "##### Copyright 2018 The TensorFlow Authors.\n",
        "\n",
        "Licensed under the Apache License, Version 2.0 (the \"License\");"
      ]
    },
    {
      "metadata": {
        "id": "qS8MroChhSxR",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "#@title Licensed under the Apache License, Version 2.0 (the \"License\"); { display-mode: \"form\" }\n",
        "# you may not use this file except in compliance with the License.\n",
        "# You may obtain a copy of the License at\n",
        "#\n",
        "# https://www.apache.org/licenses/LICENSE-2.0\n",
        "#\n",
        "# Unless required by applicable law or agreed to in writing, software\n",
        "# distributed under the License is distributed on an \"AS IS\" BASIS,\n",
        "# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.\n",
        "# See the License for the specific language governing permissions and\n",
        "# limitations under the License."
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "metadata": {
        "id": "0ufLLTrrPIhi",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "# Probabilistic PCA\n",
        "\n",
        "<table class=\"tfo-notebook-buttons\" align=\"left\">\n",
        "  <td>\n",
        "    <a target=\"_blank\" href=\"https://colab.research.google.com/github/tensorflow/probability/blob/master/tensorflow_probability/examples/jupyter_notebooks/Probabilistic_PCA.ipynb\"><img src=\"https://www.tensorflow.org/images/colab_logo_32px.png\" />Run in Google Colab</a>\n",
        "  </td>\n",
        "  <td>\n",
        "    <a target=\"_blank\" href=\"https://github.com/tensorflow/probability/blob/master/tensorflow_probability/examples/jupyter_notebooks/Probabilistic_PCA.ipynb\"><img src=\"https://www.tensorflow.org/images/GitHub-Mark-32px.png\" />View source on GitHub</a>\n",
        "  </td>\n",
        "</table>\n"
      ]
    },
    {
      "metadata": {
        "id": "Py-bl6M32ZXY",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "Probabilistic principal components analysis (PCA) is a\n",
        "dimensionality reduction technique that\n",
        "analyzes data via a lower dimensional latent space\n",
        "([Tipping and Bishop 1999](#1)). It is often\n",
        "used when there are missing values in the data or for multidimensional\n",
        "scaling."
      ]
    },
    {
      "metadata": {
        "id": "WHNNrlLNPbpA",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "## Imports"
      ]
    },
    {
      "metadata": {
        "id": "a0_ac6wS2gti",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "!pip install -q tensorflow-probability"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "metadata": {
        "id": "mbM5dCFdUior",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "import matplotlib.pyplot as plt\n",
        "import numpy as np\n",
        "import seaborn as sns\n",
        "\n",
        "import tensorflow as tf\n",
        "import tensorflow_probability as tfp\n",
        "from tensorflow_probability import edward2 as ed\n",
        "import warnings\n",
        "\n",
        "plt.style.use(\"ggplot\")\n",
        "warnings.filterwarnings('ignore')"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "metadata": {
        "id": "QjIXsU55eq-d",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "## The Model\n",
        "\n",
        "Consider a data set $\\mathbf{X} = \\{\\mathbf{x}_n\\}$ of $N$ data\n",
        "points, where each data point is $D$-dimensional, $\\mathbf{x}_n \\in\n",
        "\\mathbb{R}^D$. We aim to represent each $\\mathbf{x}_n$ under a latent\n",
        "variable $\\mathbf{z}_n \\in \\mathbb{R}^K$ with lower dimension, $K <\n",
        "D$. The set of principal axes $\\mathbf{W}$ relates the latent variables to\n",
        "the data.\n",
        "\n",
        "Specifically, we assume that each latent variable is normally distributed,\n",
        "\n",
        "\\begin{equation*}\n",
        "\\mathbf{z}_n \\sim N(\\mathbf{0}, \\mathbf{I}).\n",
        "\\end{equation*}\n",
        "\n",
        "The corresponding data point is generated via a projection,\n",
        "\n",
        "\\begin{equation*}\n",
        "\\mathbf{x}_n \\mid \\mathbf{z}_n\n",
        "\\sim N(\\mathbf{W}\\mathbf{z}_n, \\sigma^2\\mathbf{I}),\n",
        "\\end{equation*}\n",
        "\n",
        "where the matrix $\\mathbf{W}\\in\\mathbb{R}^{D\\times K}$ are known as\n",
        "the principal axes. In probabilistic PCA, we are typically interested in\n",
        "estimating the principal axes $\\mathbf{W}$ and the noise term\n",
        "$\\sigma^2$.\n",
        "\n",
        "Probabilistic PCA generalizes classical PCA. Marginalizing out the the\n",
        "latent variable, the distribution of each data point is\n",
        "\n",
        "\\begin{equation*}\n",
        "\\mathbf{x}_n \\sim N(\\mathbf{0}, \\mathbf{W}\\mathbf{W}^\\top + \\sigma^2\\mathbf{I}).\n",
        "\\end{equation*}\n",
        "\n",
        "Classical PCA is the specific case of probabilistic PCA when the\n",
        "covariance of the noise becomes infinitesimally small, $\\sigma^2 \\to 0$.\n",
        "\n",
        "We set up our model below. In our analysis, we assume $\\sigma$ is known, and\n",
        "instead of point estimating $\\mathbf{W}$ as a model parameter, we\n",
        "place a prior over it in order to infer a distribution over principal\n",
        "axes."
      ]
    },
    {
      "metadata": {
        "id": "Gwc6maZifKnb",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "def probabilistic_pca(data_dim, latent_dim, num_datapoints, stddv_datapoints): # (unmodeled) data\n",
        "  w = ed.Normal(loc=tf.zeros([data_dim, latent_dim]),\n",
        "                scale=2.0 * tf.ones([data_dim, latent_dim]),\n",
        "                name=\"w\")  # parameter\n",
        "  z = ed.Normal(loc=tf.zeros([latent_dim, num_datapoints]),\n",
        "                scale=tf.ones([latent_dim, num_datapoints]), \n",
        "                name=\"z\")  # parameter\n",
        "  x = ed.Normal(loc=tf.matmul(w, z),\n",
        "                scale=stddv_datapoints * tf.ones([data_dim, num_datapoints]),\n",
        "                name=\"x\")  # (modeled) data\n",
        "  return x, (w, z)\n",
        "\n",
        "log_joint = ed.make_log_joint_fn(probabilistic_pca)"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "metadata": {
        "id": "x2zF-wrTVSK4",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "## The Data"
      ]
    },
    {
      "metadata": {
        "id": "n4HYhXlTVZ1_",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "We can use the Edward2 model to generate data."
      ]
    },
    {
      "metadata": {
        "id": "b23iIkX8VVyn",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 68
        },
        "outputId": "86fc2c3d-61f8-4625-e272-5d3c9b836aeb"
      },
      "cell_type": "code",
      "source": [
        "num_datapoints = 5000\n",
        "data_dim = 2\n",
        "latent_dim = 1\n",
        "stddv_datapoints = 0.5\n",
        "\n",
        "model = probabilistic_pca(data_dim=data_dim,\n",
        "                          latent_dim=latent_dim,\n",
        "                          num_datapoints=num_datapoints,\n",
        "                          stddv_datapoints=stddv_datapoints)\n",
        "\n",
        "with tf.Session() as sess:\n",
        "  x_train, (actual_w, actual_z) = sess.run(model)\n",
        "\n",
        "print(\"Principal axes:\")\n",
        "print(actual_w)"
      ],
      "execution_count": 5,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Principal axes:\n",
            "[[-0.10648464]\n",
            " [ 1.3358579 ]]\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "metadata": {
        "id": "O2ZdIFz7VuSO",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "We visualize the dataset."
      ]
    },
    {
      "metadata": {
        "id": "ubJJvk0KVyVW",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 364
        },
        "outputId": "5f6effe3-f36d-474d-eb08-49f061c8b012"
      },
      "cell_type": "code",
      "source": [
        "plt.scatter(x_train[0, :], x_train[1, :], color='blue', alpha=0.1)\n",
        "plt.axis([-20, 20, -20, 20])\n",
        "plt.title(\"Data set\")\n",
        "plt.show()"
      ],
      "execution_count": 6,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeoAAAFbCAYAAAAa+83qAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3XmQnMV9xvGn33dmdnZXe4x2V8fo\nsBAywpjVbRtxOogyKlRRuRwrVXEAqxBOLBxIEROMyiJghEI5HIEAwTYOpyBlUSFAcKDKIGycklW2\nkWSEDg6vZB3LHtr7nON9O3+87EpCEkK7s7vv7H4/VSqk2Zl3uhnEM91v96+NtdYKAACEkjPSDQAA\nACdHUAMAEGIENQAAIUZQAwAQYgQ1AAAhRlADABBiBDUAACFGUAMhc+mll+rzn/+8qqurVV1drQsv\nvFB/93d/py1btpzWddrb2/XMM88MUStP7fnnn1d9ff2IvT8wWhDUQAj9wz/8g3bs2KHt27fr2Wef\n1bx58/S3f/u3euKJJz71NX7zm9/oP//zP4eukZ/A8zzdddddamhoGJH3B0YTghoIMdd1NX36dF17\n7bVav3697rnnHh04cECSlE6ntW7dOl188cWaP3++li9frl/96leSpBdeeEE33nijPvjgA1VXV2vr\n1q2y1urBBx/UpZdeqvnz5+vyyy/Xf//3f5/0vZubm3XDDTfoS1/6kubPn6+/+Iu/OGZUv23bNl15\n5ZVatGiRvvSlL+nWW29Vd3e3JGnu3Llqb2/XX/3VX+m2224bwn9DwBhgAYTKn/3Zn9mf/vSnJ/zZ\nRRdd1P+zRx55xC5ZssQ2NDTYbDZr/+M//sPOmzfPtre3W2ut/bd/+ze7bNmy/te+9NJLduHChbam\npsb6vm//93//15599tm2pqbmhO+1du1ae80119iOjg6byWTs008/bS+88EKbyWRsfX29XbBggX36\n6adtOp22hw4dsl//+tftHXfcYa219sCBA/ass86yb7/9di7/1QBjEiNqII/MnDmzf0R97bXX6oUX\nXlBVVZVc19WyZcvU3d2tP/7xjyd87RVXXKE33nhDZ5xxhowxWrp0qVzX1a5du074/Pb2dkWjUcXj\ncUUiEV155ZV68803FYlE9POf/1wTJ07UlVdeqWg0qmQyqdWrV+v5558fsr4DY1VkpBsA4NPzPE+u\n60qSWlpa9M///M/asmWLOjo6ZIyRJKVSqRO+NpVK6V/+5V/0xhtvqK2tTZKUyWRO+vy/+Zu/0XXX\nXaeLL75Y559/vr785S9r6dKlikQi2rt3r/bu3avq6urj2tfc3Jyr7gIQQQ3kjXQ6rXfffVdLly6V\nFCw4y2Qy2rhxo6ZOnaqmpiZdcMEFJ339HXfcoW3btumxxx7TrFmz5DiO5s6de9Lnf/7zn9drr72m\nzZs361e/+pXuvPNOPfPMM9qwYYPi8bgWLlyoDRs2nPC1ffeqAQweU99Annj66aeVSqV0+eWXS5K2\nb9+uFStWaNq0aTLGaOfOnZ/4+u3bt2vZsmU666yz5DiOPvjgA/X29p70+e3t7ZKkSy65RP/0T/+k\n5557Tlu3btWePXv0mc98Ru+//74ymUz/8zs6OvpH6gByh6AGQq6lpUWPPfaYHnjgAd16662qrKyU\nJE2dOlXbt29XJpPR9u3b9V//9V9yHKd/73JBQYGamprU3Nys3t5eTZ06VTt27FAqldL777+v+++/\nXxUVFSfd6/yXf/mXuv/++9Xd3S3f9/WHP/xBsVhMyWRSf/7nfy7f93Xvvfeqs7NTzc3Nuvnmm7Vm\nzRpJUjwelyTt27dPnZ2dw/BvCRi9jLXWjnQjABxx6aWXqr6+Xo4TfI8uKCjQvHnzdO211+q8887r\nf95vfvMb3XbbbWpoaNDcuXN111136YEHHtArr7yie++9V7NmzdKqVat0+PBh3XfffZo6daq+973v\nad++ffrsZz+rO++8U6+++qoee+wx3XTTTbr66quPace7776rO++8s3+kfsYZZ+iGG27QJZdcIkl6\n6623dPfdd2vXrl0aN26cLrroIq1Zs0bl5eWSpL//+7/X66+/rssuu0z333//cPyrA0YlghoAgBBj\n6hsAgBAb1KrvDRs2aPfu3fJ9X1/96ld15pln6qGHHpLv+yovL9f111+vaDSaq7YCADDmDHjq+513\n3tH//M//aM2aNero6NDNN9+s6upqzZ8/X4sXL9azzz6ryspKfeUrX8l1mwEAGDMGPPV9zjnn6MYb\nb5QkFRcXK5VKaefOnVq0aJEkadGiRXr77bdz00oAAMaoAQe14zj9WzA2bdqk+fPnK5VK9U91l5aW\nqrW1NTetBABgjBr0YrLf/e532rRpk1atWpWL9gAAgKMMajHZ9u3b9fzzz+v73/++ioqKFI/HlU6n\nFYvF1NzcrEQiccpr1NbWDqYJoZZMJulfHqN/+Ws0902if/kumUye1vMHPKLu7u7Whg0bdMstt2jc\nuHGSpOrq6v7zards2aJ58+YN9PIAAECDGFFv3rxZHR0d+td//df+x77zne/oRz/6kV577TVVVlb2\nVzACAAADM+Cgvuyyy3TZZZcd9/itt946qAYBAIAjqEwGAECIEdQAAIQYQQ0AQIgR1AAAhBhBDQBA\niBHUAACEGEENAECIEdQAAIQYQQ0AQIgR1AAAhBhBDQBAiBHUAACEGEENAECIEdQAAIQYQQ0AQIgR\n1AAAhBhBDQBAiBHUAACEGEENAECIEdQAAIQYQQ0AQIgR1AAAhBhBDQBAiBHUAACEGEENAECIEdQA\nAIQYQQ0AQIgR1AAAhBhBDQBAiEUG8+L9+/fr7rvv1rJly7R06VI9/PDDqqmpUUlJiSRp+fLlWrBg\nQU4aCgDAWDTgoO7t7dXjjz+uc88995jHv/GNb2jhwoWDbhgAABjE1Hc0GtWaNWuUSCRy2R4AAHAU\nY621g7nAxo0bVVpa2j/13draqmw2q7KyMl1zzTUqLS3NVVsBABhzBnWP+uMuvvhilZSUaMaMGXrh\nhRf03HPPadWqVZ/4mtra2lw2IVSSyST9y2P0L3+N5r5J9C/fJZPJ03p+Tld9V1dXa8aMGZKkRYsW\naf/+/bm8PAAAY05Og/qee+5RfX29JGnnzp2aNm1aLi8PAMCYM+Cp75qaGj311FNqbGyU67rasmWL\nli5dqvvvv1+xWEzxeFzXXXddLtsKAMCYM+Cgnjlzpm6//fbjHj/vvPMG0x4AAHAUKpMBABBiBDUA\nACFGUAMAEGIENQAAIUZQAwAQYgQ1AAAhRlADABBiBDUAACFGUAMAEGIENQAAIUZQAwAQYgQ1AAAh\nRlADABBiBDUAACFGUAMAEGIENQAAIUZQAwAQYgQ1AAAhRlADABBiBDUAACFGUAMAEGIENQAAIUZQ\nAwAQYgQ1AAAhRlADABBiBDUAACFGUAMAEGIENQAAIRYZzIv379+vu+++W8uWLdPSpUt1+PBhPfTQ\nQ/J9X+Xl5br++usVjUZz1VYAAMacAY+oe3t79fjjj+vcc8/tf2zjxo26/PLLdccdd2jSpEl64403\nctJIAADGqgEHdTQa1Zo1a5RIJPof27lzpxYtWiRJWrRokd5+++3BtxAAgDFswFPfruvKdd1jHkul\nUv1T3aWlpWptbR1c6wAAGOMGdY86F5LJ5Eg3YUjRv/xG//LXaO6bRP/GkpwGdTweVzqdViwWU3Nz\n8zHT4idTW1ubyyaESjKZpH95jP7lr9HcN4n+5bvT/RKS0+1Z1dXV2rJliyRpy5YtmjdvXi4vDwDA\nmDPgEXVNTY2eeuopNTY2ynVdbdmyRTfccIMefvhhvfbaa6qsrNQll1ySy7YCADDmDDioZ86cqdtv\nv/24x2+99dbBtAcAAByFymQAAIQYQQ0AQIgR1AAAhBhBDQBAiBHUAACEGEENAECIEdQAAIQYQQ0A\nQIgR1AAAhBhBDQBAiBHUAACEGEENAECIEdQAAIQYQQ0AQIgR1AAAhBhBDQBAiBHUAACEGEENAECI\nEdQAAIQYQQ0AQIgR1MAY5XlSOh38E0B4RUa6AQCGl7VSY6NROi3V1TkyxqigwFdVlZUxI906AB/H\niBoYYxobjbJZo1hMisWkaNQqmzVqbCSlgTAiqIExxPOkVMqR87G/+Y4TPM40OBA+BDUwhnieZK09\n4c+stQQ1EEIENTCGuK5kTnIj2hgj1x3mBgE4JYIaGENcVyoo8OX7xz7u+8HjBDUQPgQ1MMZUVVlF\nIlbpdLA9K5MxikSsqqpOPCUOYGTldHvWzp07dd9992natGmSpOnTp+uaa67J5VsAGCRjpAkTrCZO\nlKRgFM1IGgivnO+jPuecc/Td734315cFkGOuG2zPAhBuTH0DABBiOR9RHzx4UD/84Q/V2dmpFStW\naM6cObl+C2DM8rzgF9PVwNhh7Mk2VQ5Ac3Oz9uzZo8WLF6u+vl4/+MEP9OCDDyoSoVIpMBjWSg0N\nUk9P8HtjpMJCacIEnbTs5+mEet/Csr5qZQDCI6cJOn78eJ1//vmSpEmTJqm8vFzNzc2aMGHCSV9T\nW1ubyyaESjKZpH95LEz9a2gIyn4eXVHM96UPP7SaMOHY79p9tbxTKUfW2pPW8k4mkzp4sFbvvuuo\no8OR7wcVykpKfM2e7cva/B29h+mzGwr0L78lk8nTen5Og/rXv/61WlpatHz5crW2tqqtrU3jx4/P\n5VsAY05f2c9o9NhAPlL205PrHhlBNzUZWWuOen5fLW8dF+rvvuuot9dRLKb+oO7pcbRli9G0afrE\noAcwPHIa1IsWLdIDDzyg3//+98pms7r22muZ9gYG6VRlP7PZIJxTKUeplFV9vaOSEntMsFordXY6\nKinx+kfInZ1Sfb0j3zfyPEeeZ2WtUW+vVUGBr6lT/Y+mwU8e9ACGXk5TtLCwULfccksuLwmMeUfK\nfh4JSc8LRsC+b9TcbOT7Ru3tUleXUXOzo7Y2q5YWq5kzfTU1GWUy0uHDjj780KigwKq5Obju9u0R\nFRQYGWNVWSlZG1wvEnE0e3a6/371x0fvAIYPw10g5PrKfmazRsZIzc1GqVQQ0q5rZYyrSMQqlTJq\nbTVqbw9Ox9q71+q99xxNnSp1dEjZrKNYzNM777hqb3eUSEh/+pNRRYVkTDCirqqychyjdNqqtdWo\nuPjIl4O+QzsIamB4EdRAHqiqsmpslGprXWUykusaxeO+SkutamuNamuDKexUylFdndTebuS6wWKy\nnh4r15V6e40+/DCiQ4eCEXhnp9Tc7KqnRxo3zqq52ZG1nrq6pOJiKZWSenulaPTIqJ6QBoYfQQ3k\nAWOkioogdB3HynHUv4Csrc1Re7tRV5dRe7tRZ6fR4cOOjLFyXaOODqmiwpcx0ocfGtXUuPI8RxMm\nSF1djpqafJWVBc8vKwu2Z2WzRnv3unIcX7GYVTQqTZjAoR3ASCCogTzheZIxQUgfffpVW5v03nuO\nOjtdtbZKjhOsIEunjbJZX4WFRu+956ilxVVPj6O2tr4pdCkSMbLWUVGRFI9LM2YEU92trUZtbY4k\nX5Ija/3jGwRgWBDUQJ5wnL6tV8GUtrVGLS1Wvb1GqZSRtVaO48iYoDiK77vKZIzSaUctLVJPj5Hv\nS9lscL3GRsl1HUUiVgUFVomEr+Zmo8JCq0mTrLJZXxUVvoqKgtF7Os1iMmAkENRAnghWY0ueZ9XW\nZtTT4+j9940OHAiqoEQiVtlssO+5q8vI2mBPdCYjpVJGPT3BFLrjBKNzKbgHHYkYtbZaFRQY7dlj\nVFrqqLzc1/Tpjjo6jMaNC0bZLCYDRgZBDeSBvqInlZVWH3wQrPpuaQmCNxq1/YHc2mrU0OCouzsY\nWff0HLmXfXRRlL5rSsEIu6XFkbVSeXlwP9v3g/CWjJqapMpKy2IyYIQQ1EAe6Ct64vtScbFRUZFV\nJhPcP66rc3T4sFVnZ0QtLY56eoIgD0bfQRBnMp98fWuDqe26OquJE4PFYz09Rp4n9fY66u72VFzM\nYjJgJHDMJZAH+rZH+X5fWU8pFjPq7Q1GvI7jqrXVUXd3UNyk75CNdPrUIS0Fz+npkRobHf3xj0Hx\nlIYGozffdLVjR1AopbMz+GfujvEB8GkwogbyQF/Rk1TKyBgjx7HyvGBvdTbrqKkp2JaVTgf7n7PZ\nI4vGPq1USpKMGhpiSqU8FRV5ikYdNTZaZTKuenqsysul1lars8/2qfsNDBOCGsgTfUVPjAmqlDU3\nSy0twT/b2pz+UXE6PfD3SKWkw4eltjZXpaWOPC+rZFI6eNBXIuH1lxRtaAimyAEMPaa+gTxhTHAo\nRnW1p6oqT9msVUODo6YmR6mU+WjB2eDfp7c3KDka1A131dZm1NTkaM8eR93dwdR7V5fTvxgNwNBi\nRA3kmb4D6TwvWKldVGTU0hIEbC51dQX7rGtrHaXTnjo6rKZNczR+vC9jfLZqAcOEETWQZzxPOnTI\nVU+PUTQa7KvOZgc35X0imUwwsm5rk9rbHTU2RrR5s1FHR1CalJAGhgdBDeQRa6X9+4127zbats3R\nO+9E9Kc/BaVDh0I6Hdz37u42ymalzs6IamsddXZaghoYJkx9A3nCWmnPHkc1NRHt2uWqoyPSf296\nKPX2Sq2tUjweLPNOpYL710x9A8ODETWQJ+rqjDo6jLJZq/p6R42Njrq6jLq7h/Z9nY/+L2GMVTwe\n3BtPpZycT7UDODGCGsgDnhfU7bbW6MABV729jrq6nJys8j6VTCb41dkZnKLV2iplMr7SabHyGxgG\nBDWQB/qOuPS84EzqTMaqu/v0i5oMVHe31Nvra9y4YFR/4ICr5mZHtbWuGhqoVgYMJYIayAOuG5xy\nFY2q/1jLTGb4RrSOI2WzRgcPBtXPEgn70UlcVtmsUWMjZcqAocJiMiAP9JUQdRxHhYVWkYgd1nvE\nvh+s/D582FF5uad9+xxFo1JVla94XCopMZxVDQwRRtRAnkgkrNrbJd83chyjceOG9/2Dk7isiosd\nSY46O4Nf2WxwMAj3q4GhQVADeeLDD42iUaNZs6wmTszK2uFNRs8LaoAH510HZ18HC9ykTIapb2Co\nENRAyFkr1dcbHTzoat8+R7t2OUqnHcViwx+O7e1GBw4EB4EEtcCDBW2RCKvJgKHCPWog5Bobg3On\nu7qMYjGjnh4jz3NG5H6wtcFJXVu3RnXokK9p06TJkz3NnEnxE2CoMKIGQiw4EcuR40ieZ1RaamWt\n1eHDVum0MyJnQhtjJBnFYsEBIYcPO+rupqQoMFQIaiDEPC84VlKSXNeqpcWopcWop8cdtj3UH5dO\nB6Pq3l6jdNrR4cNSNDr0pUyBsYqgBkLMdSVrgxD0faP6eqPu7mBE6zhHynsOp0hEisd9FRb6Gj/e\nV2Gho7o6Vn0DQ4WgBkLKWqmpKRhB19c72rs3uD89YYKV6/ojVmvb9yUpaEtxsRSJ+JJY9Q0MlZwv\nJnviiSf0/vvvyxijlStXatasWbl+C2BMaGw0ymaNJk70VVsbVASrrzdqbHRlbbDaeqRGsa4rGSN1\ndlpNmyYVF7PqGxgqOR1R79q1S3V1dVq/fr2+/e1v6/HHH8/l5YEx4+hFZMboo3KdRhMnWlVW+iou\n9kdkIZkUHHvZ2RnUGrdW+sxnPFVUsOobGCo5DeodO3boC1/4giRp6tSp6urqUvdQn8EHjEJHLyLz\nvOD0qmjUKhaTGhqkTMZVPK4RCetsVioslCoqPI0bZxWPW0WjviRxohYwBHI69d3a2qqZM2f2/7m0\ntFStra0qKio66WuSyWQumxA69C+/jVT/+sIuFgvCr7dXmjlT2rs3COfSUqmrS2pvD0J84ApP+xWO\nE9yn9v1iffih1NgojRsn1ddL48cfafeECSPzRaIP/23mt9Hev9MxpAVP7Kc4+662tnYomzCikskk\n/ctjI92/1tbgHrW1UktLUOSkt9dVPO5+dDiHkeMMZr65UFLPab/KmODIy9ZWK2s9/f73WTU2Wp1x\nhq+WFquKCivflz780GrChJG5dz3Sn91Qo3/57XS/hOR06juRSKi1tbX/zy0tLUokErl8C2DMqKoK\nTsnyfdO/TSuRsJo5M7hP3dNjVFAw/O0KpuKDfdTFxZIxjnp7HR086CiVCn7uOME9dqbBgcHLaVDP\nnTtXW7ZskSTV1NQokUiosPD0p9YABCPXCROskklP55zjac6cjKqqPBUV+ers1Ijto5aC6XbflzIZ\nK9cN7p13dEiplPlo+1Ywo0ZQA4OX06nv2bNna+bMmVq7dq2MMVq1alUuLw+MSa4b/JoyxcpxPJWU\nWHV2SgcOuP2hONyCoLbq7nbkukFYS448z+v/8mCMYSU4kAM5v0f913/917m+JICPVFVZHTpkNGOG\nr/fey6qhwVV7+/C3w/elbNYok5GMCe5JW+uruDj4UuH7UkGBT1ADOUBlMiCP+L5UXu4rErGaNs0q\nm7UjMv0djUrl5VaFhVIqJfX2WhUXS+PHW2UyRpGIVVUVRVCAXOCYSyCP9C0qk4wqKoI9zI5jhn0K\nPBKRPC9Y7FZQYPW5z2VUXe3L2iNT9QByg6AGQs7zgl99AVhQEByE0dlpZO3wh7QUjKgjEekzn/E1\nebKvz33OVzQ6/O0AxgKCGggpa4N636mUI2utjDEqKPA1eXJwJnVTU7BYq6REamsbvnYFp2dZjR9v\nVVholU6zaAwYStyjBkKq71COvtKh0ahVNmt0+LBRMhlszSostMNe/cv3pVTKqrPTl+8bFRX5bMMC\nhhBBDYTQ0YdyHM1xpJ4eR5GIr1jMauJET9LwLtoqLpaqqqSCAlc9PUHBEwBDh79hQAgdfSjH8axa\nW52P6m1bFRQM75DacaSyMl8lJb6yWSPPYxsWMJQIaiCEgvOeTxzA7e2OSkut4nEpnY4oEhneVdZF\nRVJbm6OWFqNUSvJ9R42NQU1yALlHUAMh1Le6++MrujOZYKQ9ZYrVmWdmFY16ikaPf95QtiuVknp6\njIyxmjrVatIkT75v1Ng4gkdlAaMYQQ2EVN+hHJmMUTotZTJGjhOstnZd6ZxzfJ1xhqe5cz1Nnjw8\nq7n6zsb2fSPPk6JRX0VFwXYtDuEAhgbbs4CQ6juUw/O8/n3UklRb60gKKpNNm2Z16JCv2bODs6k7\nO4d2Djwe72uHL9c1qqsLTsxqajIaN84e004AucGIGgg515VisWMLnvh+8PsvfzmjsjJfxcVWkyYF\nzxsqjhPcny4r8zVpkqfPftbvH/Fns0YtLQ4hDQwBRtRAnqmqsmpsDKaaS0qsFizIqqTEyveDEe3e\nvUOTln3HWhYVWVVUBCN+ax1ls8HPT75KHcBgENRAnjl6SrynJ3hs1ixf6bTU3W3U3u6qqSn375tI\nBIduzJrlq6zMKpORsllfkYiR6/oqL2fqGxgKTH0Decp1pcJCKRIxKiyUzj7bavLk4JCMoqLcv5/v\nS55n1N0dTIOXlFideaavadM8VVYGh4MQ0kDuMaIG8ljfPets1qi42Fci4auy0pe1rjKZYIV2rgRT\n28FKdGOsiout5swJ9oVx/jQwdAhqIM9VVVk1NAQrsquqrMrLPbW1GZWVOWpry01YR6NSImF11lkZ\nnXeer3HjrBzHVyZj5PvBYSGcPw0MDYIayHPGBCG6f7/VnDmeslkpErH605+iikQctbdL3d0Dv37f\nau+CAikeD/Z0V1RYdXVJEyZkNXWqZSQNDCGCGhglrDUqLpZmzw6mwnt7s2pri8jzgu1T6fTArltY\nGCwkKywMjtNsbXW0d69UXOzL9x01NQWj6eE+xQsYK1hMBowS0aiV7wcrs6uqfJ15pq8pU7IqLw/2\nWEdO82u56wYj6XHjpFjMqrTU9p9FXVoqTZoU7N/OZikfCgwlRtTAKOC6UkWF1N5udfiwUTYbBKy1\nvurqgspmxcXBPetPEon0nXN95LzrsjLbf1JXPB7c8zbGqrLyyJR3UD7UYwocGAIENTAKuK4Uj/tq\na3NUXi55nlVnp6+uLkelpVJvb1CkJBIJSo0evcDMmL7TuoLRcmGhUSKRVVGRVF5uFY0GwTxjhi/J\nqKgoCOijF49Zyx5qYKgQ1MAoMX681aFDwT1ka6WKCiPft4rFPFkrdXYaGROU+cxmjyww6ytPWlAQ\nHLJRXGw1e7anWEzq7Q1WdUciwc9931c8HgT50fekjWEPNTBUCGpglPB9fVQdzMpxrBIJX47jqrvb\nqLJSchxf1lqlUq6iUam4ODiyUuo7bMPTuHG+Kir6TsYKKp1Ze2QqPBIJ7lcbY/rrjbOHGhhaBDUw\nSrhusPLbmGC62nGCAzTmzg22UpWVSd3djt57z1N3d+SjKW1fiYRUVOSpuNgqmZTS6aC6WVWV5Hm+\niouD7VgVFVbRqK8JE6y6u/tCnj3UwFAjqIFRwnWDLVPNzRFJtn86O52W5s711N7uqK0tq3jc0f79\nGXmeVFjoKZGwisV8JRLBkZVFRUfCOZPpOzHLatIkX9ZKpaW2/8+FhdyXBoYaQQ2MIhMmWDU3e+rq\nMjLGqKjIV0eHo1jMKp02stZqxgyr887Lqr5e+tOfPHmekeMYFRQEW7BcNwj3YOW4lbXB1LYUjLBj\nMavychHSwDAhqIFRxBjp7LN91dUZ9fQ4ammRJk6UXDe4Z93W5sh1fdXVOTrzzODYytpayXEcSb4y\nmaDKWVtbsDfa2mCRWmmpkeQrmQy2ZUUiVCMDhkvOgvqXv/ylfvazn2nixImSpDlz5uhrX/tari4P\n4FMyRpo82Sqd9iS5Kijw5PtSba2V6wYVxNrbpSlTghDv6ZEyGV+lpVJXl1Vhoa+iIqOGBqOqKk9V\nVUGgS8Eis0jEck8aGEY5HVEvXrxYV199dS4vCWAQjh75um4wHZ7NWvm+079qu6QkqDRWUGDV1WU0\naVKw/7qiwtNZZ/lyHMna4AuA5xlVVFAuFBhOTH0Do1RQxMSob2FZPO4rEpF838hxggCWgsVihYVW\n48dbFRVZTZ7sq60tuD9dUHDsNSlsAgy/nAb17t27tX79enmep6uuukpnnHFGLi8P4DQcfVa14wSr\nuCVftbWOCgocdXcHI+Xp0z05TlAApaAgqEzW1RWMtD8eyhQ2AYafscFp8Kfl9ddf16ZNm4557IIL\nLtCkSZO0YMECvffee/rxj38pjXnGAAAM/ElEQVSse++9N2cNBXD6rJUaGoL70H3T17GYVFoq7dlz\npDqZ4wS1wWfPllpapEOHgjKjfVu8qqqC10ejwX1tAMNnQEH9aXzrW9/Sj3/8449Wk55cbW3tULx9\nKCSTSfqXx0ZT/zxP/aPjvhFxMpnUvn21SqeD8I7FpIaG4EhMY6SmJqN02pHnBQVUkkkvb46zHE2f\n3YnQv/yWTCZP6/k5m/p+8cUXVVFRoQsvvFD79+9XaWnpKUMawPA4OqCP1hfQUhDkqZSjaDT47l5Z\naeV53kcnZ7GIDBgpOQvqCy+8UA899JB+8YtfyPd9rV69OleXBjAMPC9YLHa0voBPp1lEBoyUnAV1\nRUWFbrvttlxdDsAwO3qV+MexiAwYOcxNA5B0ZJW47x/7OKdjASOLoAbQr6oqKA+ayRil01ImY6hE\nBowwCp4A6GdMcLCH53nHrRIHMDIIagDHIaCB8GDqGwCAECOoAQAIMYIaAIAQI6gBAAgxghoAgBAj\nqAEACDGCGgCAECOoAQAIMYIaAIAQI6gBAAgxghoAgBAjqAEACDGCGgCAECOoAQAIMYIaAIAQI6gB\nAAgxghoAgBAjqAEACDGCGgCAECOoAQAIMYIaAIAQI6gBAAgxghoAgBAjqAEACLEBB/WuXbt07bXX\n6q233up/bN++fVq7dq1uvfVWPfroozlpIAAAY9mAgrqurk4vv/yyZs+efczjTz75pFauXKl169ap\nu7tb27Zty0kjAQAYqwYU1IlEQjfddJOKior6H8tms2poaNCsWbMkSQsXLtSOHTty00oAAMaoyEBe\nVFBQcNxj7e3tKi4u7v9zWVmZWlpaBt4yAABw6qB+/fXXtWnTpmMeW7FihebNm/eJr7PWfqoGJJPJ\nT/W8fEX/8hv9y1+juW8S/RtLThnUS5Ys0ZIlS055odLSUnV0dPT/ubm5WYlE4pSvq62tPeVz8lUy\nmaR/eYz+5a/R3DeJ/uW70/0SkrPtWZFIRFOmTNGePXskSb/97W9POeoGAACfbED3qLdu3aqXXnpJ\nhw4dUk1NjV555RWtXbtWK1eu1E9+8hNZazVr1izNmTMn1+0FAGBMGVBQL1iwQAsWLDju8alTp+qO\nO+4YdKMAAECAymQAAIQYQQ0AQIgR1AAAhBhBDQBAiBHUAACEGEENAECIEdQAAIQYQQ0AQIgR1AAA\nhBhBDQBAiBHUAACEGEENAECIEdQAAIQYQQ0AQIgR1AAAhBhBDQBAiBHUAACEGEENAECIEdQAAIQY\nQQ0AQIgR1AAAhBhBDQBAiBHUAACEGEENAECIEdQAAIQYQQ0AQIgR1AAAhBhBDQBAiBHUAACEWGSg\nL9y1a5fuu+8+rV69WgsXLpQk3X777UqlUiooKJAkXX311Zo5c2ZuWgoAwBg0oKCuq6vTyy+/rNmz\nZx/3s9WrV2v69OmDbhgAABjg1HcikdBNN92koqKiXLcHAAAcZUAj6r6p7RPZuHGjOjo6NGXKFK1c\nuVKxWOwTr5VMJgfShLxB//Ib/ctfo7lvEv0bS04Z1K+//ro2bdp0zGMrVqzQvHnzjnvuFVdcoenT\np2vSpEl69NFH9eqrr2r58uWfeP3a2trTbHL+SCaT9C+P0b/8NZr7JtG/fHe6X0JOGdRLlizRkiVL\nPtXFvvjFL/b/fuHChdq8efNpNQYAABwrZ9uzrLVat26durq6JAWrwllUBgDA4AzoHvXWrVv10ksv\n6dChQ6qpqdErr7yitWvXasmSJbrjjjsUj8eVSCS0YsWKXLcXAIAxZUBBvWDBAi1YsOC4x88//3yd\nf/75g24UAAAIUJkMAIAQI6gBAAgxghoAgBAjqAEACDGCGgCAECOoAQAIMYIaAIAQI6gBAAgxghoA\ngBAjqAEACDGCGgCAECOoAQAIMYIaAIAQI6gBAAgxghoAgBAjqAEACDGCGgCAECOoAQAIMYIaAIAQ\nI6gBAAgxghoAgBAjqAEACDGCGgCAECOoAQAIMYIaAIAQI6gBAAgxghoAgBAjqAEACLHIQF7keZ4e\neeQR1dfXy/d9XXXVVTr77LO1b98+/fSnP5UxRtOnT9e3vvWtXLcXAIAxZUAj6jfffFPxeFzr1q3T\nt7/9bT355JOSpCeffFIrV67UunXr1N3drW3btuW0sQAAjDUDCuqLLrpIV199tSSptLRUnZ2dymaz\namho0KxZsyRJCxcu1I4dO3LXUgAAxqABTX1HIkde9vOf/1wXXHCB2tvbVVxc3P94WVmZWlpaBt9C\nAADGsFMG9euvv65NmzYd89iKFSs0b948vfrqq9q7d6++973vqb29/ZjnWGs/VQOSyeRpNDf/0L/8\nRv/y12jum0T/xpJTBvWSJUu0ZMmS4x7ftGmT3nrrLf3jP/6jIpGISktL1dHR0f/z5uZmJRKJ3LYW\nAIAxZkD3qOvr6/WLX/xCN910k2KxmKRgOnzKlCnas2ePJOm3v/2t5s2bl7uWAgAwBhn7aeeoj/Ls\ns89q8+bNqqys7H9s7dq1qqur009+8hNZazVr1ix985vfzGljAQAYawYU1AAAYHhQmQwAgBAjqAEA\nCLEB7aMejLFQfnTXrl267777tHr1ai1cuFCSdPvttyuVSqmgoECSdPXVV2vmzJkj2cwBO1H/RtPn\nJ0m//OUv9bOf/UwTJ06UJM2ZM0df+9rXRrhVufHEE0/o/ffflzFGK1eu7C9SNBrs3LlT9913n6ZN\nmyZJmj59uq655poRbtXg7d+/X3fffbeWLVumpUuX6vDhw3rooYfk+77Ky8t1/fXXKxqNjnQzB+zj\n/Xv44YdVU1OjkpISSdLy5cu1YMGCEW7lwG3YsEG7d++W7/v66le/qjPPPPO0Pr9hD+qjy48eOHBA\n//7v/6677rqrv/zorFmz9MADD2jbtm2aP3/+cDdv0Orq6vTyyy9r9uzZx/1s9erVmj59+gi0KndO\n1r/R8vkdbfHixf0V+EaLXbt2qa6uTuvXr9fBgwf1yCOPaP369SPdrJw655xz9N3vfnekm5Ezvb29\nevzxx3Xuuef2P7Zx40ZdfvnlWrx4sZ599lm98cYb+spXvjKCrRy4E/VPkr7xjW/0DwTy2TvvvKMD\nBw5o/fr16ujo0M0336zq6urT+vyGfep7tJcfTSQSuummm1RUVDTSTRkSJ+rfaPr8RrsdO3boC1/4\ngiRp6tSp6urqUnd39wi3Cp8kGo1qzZo1x9Sl2LlzpxYtWiRJWrRokd5+++2Rat6gnah/o8k555yj\nG2+8UZJUXFysVCp12p/fsI+oR3v50b6p7RPZuHGjOjo6NGXKFK1cubJ/D3o+OVH/RtPnd7Tdu3dr\n/fr18jxPV111lc4444yRbtKgtba2HnPLpbS0VK2traPqi+XBgwf1wx/+UJ2dnVqxYoXmzJkz0k0a\nFNd15bruMY+lUqn+qdK+zzBfnah/kvTqq6/q5ZdfVllZma655hqVlpaOQOsGz3EcxeNxSUGhsPnz\n5+sPf/jDaX1+QxrUQ11+dKR9Uv8+7oorrtD06dM1adIkPfroo3r11Ve1fPny4WrqgJxO/46WL59f\nnxP184ILLtCKFSu0YMECvffee3rooYd07733jlALh06+fVanMnnyZK1YsUKLFy9WfX29fvCDH+jB\nBx88ZoCA8Lv44otVUlKiGTNm6IUXXtBzzz2nVatWjXSzBuV3v/udNm3apLVr1+qGG244rdcO6X+9\no7386Mn6dyJf/OIX+3+/cOFCbd68eaialTOftn/5+vn1OVU/zzrrLLW3t8v3fTlOfm+USCQSx3x7\nb2lpyavP6lTGjx+v888/X5I0adIklZeXq7m5WRMmTBjhluVWPB5XOp1WLBbLu79vn0Z1dXX/7xct\nWqRHH310BFszeNu3b9fzzz+v73//+yoqKjrtz2/Y/68zFsuPWmu1bt06dXV1SQoW9OT7orKjjcbP\n78UXX9T//d//SQpWpJaWluZ9SEvS3LlztWXLFklSTU2NEomECgsLR7hVufPrX/9aL730kqRgmr+t\nrU3jx48f4VblXnV1df/nuGXLlrz/+/Zx99xzj+rr6yUF9+P7VvHno+7ubm3YsEG33HKLxo0bJ+n0\nP79hr0w22suPbt26VS+99JIOHTqk0tJSJRIJrV27Vps3b9aLL76oeDyuRCKh1atXf+L97LA6Wf8O\nHjw4Kj6/Pk1NTf3bJ3zf1ze/+c1Rs43pmWee0e7du2WM0apVqzRjxoyRblLO9PT06IEHHlB3d7ey\n2ay+/vWv5/W2Hin4QvXUU0+psbFRrutq/PjxuuGGG/Twww8rk8mosrJS1113Xd5O75+of0uXLtWL\nL76oWCymeDyu6667TmVlZSPd1AF57bXX9Nxzz2ny5Mn9j33nO9/Rj370o0/9+VFCFACAEMv/uTwA\nAEYxghoAgBAjqAEACDGCGgCAECOoAQAIMYIaAIAQI6gBAAgxghoAgBD7f1kLGUIGmWmMAAAAAElF\nTkSuQmCC\n",
            "text/plain": [
              "<matplotlib.figure.Figure at 0x7fc4297919e8>"
            ]
          },
          "metadata": {
            "tags": []
          }
        }
      ]
    },
    {
      "metadata": {
        "id": "qIZepLMgItfq",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "## Maximum a Posteriori Inference\n",
        "\n",
        "We first search for the point estimate of latent variables that maximizes the posterior probability density. This is known as maximum a posteriori (MAP) inference, and is done by calculating the values of $\\mathbf{W}$ and $\\mathbf{Z}$ that maximise the posterior density $p(\\mathbf{W}, \\mathbf{Z} \\mid \\mathbf{X}) \\propto p(\\mathbf{W}, \\mathbf{Z}, \\mathbf{X})$.\n"
      ]
    },
    {
      "metadata": {
        "id": "s2AvQAYqIh6K",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "tf.reset_default_graph()\n",
        "\n",
        "w = tf.Variable(np.ones([data_dim, latent_dim]), dtype=tf.float32)\n",
        "z = tf.Variable(np.ones([latent_dim, num_datapoints]), dtype=tf.float32)\n",
        "\n",
        "\n",
        "def target(w, z):\n",
        "  \"\"\"Unnormalized target density as a function of the parameters.\"\"\"\n",
        "  return log_joint(data_dim=data_dim,\n",
        "                   latent_dim=latent_dim,\n",
        "                   num_datapoints=num_datapoints,\n",
        "                   stddv_datapoints=stddv_datapoints,\n",
        "                   w=w, z=z, x=x_train)\n",
        "\n",
        "\n",
        "energy = -target(w, z)\n",
        "\n",
        "optimizer = tf.train.AdamOptimizer(learning_rate=0.05)\n",
        "train = optimizer.minimize(energy)\n",
        "\n",
        "init = tf.global_variables_initializer()\n",
        "\n",
        "t = []\n",
        "\n",
        "num_epochs = 200\n",
        "\n",
        "with tf.Session() as sess:\n",
        "  sess.run(init)\n",
        "\n",
        "  for i in range(num_epochs):\n",
        "    sess.run(train)\n",
        "    if i % 5 == 0:\n",
        "      cE, cw, cz = sess.run([energy, w, z])\n",
        "      t.append(cE)\n",
        "\n",
        "  w_inferred_map = sess.run(w)\n",
        "  z_inferred_map = sess.run(z)"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "metadata": {
        "id": "UL1PmOM34brs",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 364
        },
        "outputId": "3929bdbb-0666-444d-a9c5-414bbc2acefc"
      },
      "cell_type": "code",
      "source": [
        "x = range(1, num_epochs, 5)\n",
        "plt.plot(x, t)"
      ],
      "execution_count": 8,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "[<matplotlib.lines.Line2D at 0x7fc422fcc898>]"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 8
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAe8AAAFKCAYAAAA0WNeQAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3X98leV9//HXdeecEE7Ij4NIMEDE\nJjEiBIJJU2QrWsaErR1aRbvKOlOVtU2H/VK7VpxYNVLaaa1MQIdsLZYyC9RRq5R2DNvG2gxBBBJN\nVqAUgQqRcCAkJOHkXN8/TnJMIJAEEu5zOO/n45FHzrnPde7z+XDIeZ/7t7HWWkRERCRmOG4XICIi\nIr2j8BYREYkxCm8REZEYo/AWERGJMQpvERGRGKPwFhERiTEKbxERkRjjcbuA3jh48GCfzCczM7PP\n5hUN1E90Uz/RTf1Et3jvJzMzs8vp3YZ3VVUVTz31FCNHjgQgKyuLGTNmsHjxYkKhEOnp6cyZMwev\n10t5eTnr16/HGMPUqVOZMmUKwWCQpUuXUltbi+M4lJaWkpGRwd69e1m+fDnGGLKyspg9e3aPmxER\nEYlnPVryvvbaa7n//vsj95cuXcq0adO4/vrrWbVqFa+99hqTJ09m7dq1LFy4EI/Hw7x58yguLmbL\nli34fD7KysrYvn07q1atYu7cuaxYsYKSkhJycnJYtGgR27ZtY8KECf3WqIiIyKXivLZ5V1VVUVRU\nBEBRURE7duxg165dZGdn4/P5SExMJC8vj+rqaiorKykuLgYgPz+fmpoagsEghw8fJicnB4DCwkJ2\n7tzZRy2JiIhc2nq05L1//36+853vcOLECW6//Xaam5vxer0ApKamEggECAQCpKamRp7T1XTHcTDG\nEAgESE5OjoxNS0vj6NGjfdmXiIjIJavb8L7iiiu4/fbbuf766zl06BCPPvoora2t5/2CXV0HpafX\nRjnbhvvz0ZfzigbqJ7qpn+imfqKb+jlTt+E9ePBgJk2aBMCwYcNIT09n9+7dtLS0kJiYSF1dHX6/\nH7/fTyAQiDyvrq6O3NzcTtODwSDWWtLT06mvr+801u/3d1us9jbvmvqJbuonuqmf6Bbv/Zwt6Lvd\n5l1eXs7LL78MQCAQ4NixY9x4441UVFQAUFFRQUFBAbm5uezevZuGhgaampqoqalh9OjRjB8/PjJ2\n69atjBkzBo/Hw/Dhw6murgZg8+bNFBQU9LgZERGReNbtkndRURGLFi1iy5YtBINB7r33Xq666ioW\nL17Mxo0bGTJkCDfccAMej4dZs2axYMECjDHMnDkTn8/HpEmT2LFjB/Pnz8fr9VJaWgpASUkJy5Yt\nw1pLTk4O48aN6/dmRURELgXG9nSDcxTQavOuqZ/opn6im/qJbvHez3mvNhcREZHoEpfhbVuDhCp+\nhW2o736wiIhIlInL8G7ZXYP996ewr613uxQREZFei8vwdpJTwjdq33e3EBERkfMQl+GdMGQoAPbo\nBy5XIiIi0ntxGd7OgCRISYM6hbeIiMSeuAxvAPxD4Ghtj0/NKiIiEi3iN7wHD4GWFjihPc5FRCS2\nxG14m8GXh28crXW3EBERkV6K2/Bm8JDwb233FhGRGBO/4e0Ph7et05K3iIjElrgNb6MlbxERiVFx\nG95EtnkrvEVEJLbEb3inDQbjYLXkLSIiMSZuw9skJED6YNA2bxERiTFxG95AeI/zwBFsqNXtSkRE\nRHosrsPbDL4cQiE4FnC7FBERkR6L6/BuP1xMq85FRCSWxHd4D24/1ls7rYmISOyI6/DWKVJFRCQW\nxXV46xSpIiISi+I7vHWKVBERiUHxHd4paeDxaslbRERiSlyHtzEmvOpcp0gVEZEYEtfhDYRXnR8P\nYE+dcrsSERGRHon78I5cXUxL3yIiEiM8PRnU0tLC/fffz2233UZVVRV79uwhJSUFgBkzZnDddddR\nXl7O+vXrMcYwdepUpkyZQjAYZOnSpdTW1uI4DqWlpWRkZLB3716WL1+OMYasrCxmz57dr02eU8er\niw29wr06REREeqhH4f2Tn/yEQYMGRe7feeedFBYWRu43NTWxdu1aFi5ciMfjYd68eRQXF7NlyxZ8\nPh9lZWVs376dVatWMXfuXFasWEFJSQk5OTksWrSIbdu2MWHChL7vric6nKjFuFOBiIhIr3S72vzA\ngQPs37//nOG6a9cusrOz8fl8JCYmkpeXR3V1NZWVlRQXFwOQn59PTU0NwWCQw4cPk5OTA0BhYSE7\nd+7so3Z6z/jblrx1uJiIiMSIbsP7hRde4K677uo0bcOGDTz66KM8/fTTHD9+nEAgQGpqauTx1NRU\nAoFAp+mO42CMIRAIkJycHBmblpbG0aNH+6qf3tOJWkREJMacc7X5r3/9a66++mqGDh0amTZ58mRS\nUlIYNWoU69atY82aNeTl5fXoxay1PZp2NpmZmT0e29N5hdJTOQAMaKzn8j6c/8XWl/820UD9RDf1\nE93UT3Tri37OGd5vvfUWhw8f5q233uLIkSN4vV5mz57NqFGjACgqKuL5559n4sSJBAIfXlazrq6O\n3Nxc/H5/ZHowGMRaS3p6OvX19Z3G+v3+HhV78ODB3vbXpczMzM7zGuij6U/7+2z+F9sZ/cQ49RPd\n1E90Uz/Rrbf9nC3ozxnec+fOjdxevXo1Q4cO5Ze//CUZGRlkZGRQVVXFyJEjyc3N5bnnnqOhoYGE\nhARqamooKSnh5MmTVFRUUFBQwNatWxkzZgwej4fhw4dTXV3NNddcw+bNm5k+fXqPG+kX/iFabS4i\nIjGjR3ubdzR9+nSefvppEhMTSUpKorS0lMTERGbNmsWCBQswxjBz5kx8Ph+TJk1ix44dzJ8/H6/X\nS2lpKQAlJSUsW7YMay05OTmMGzeuzxvrlcGXw8F92KZGTJLP3VpERES60ePwvuOOOyK3Fy5ceMbj\nEydOZOLEiZ2mtR/bfboRI0bw2GOP9abOfmUGD8FCeOk7M8vtckRERM4p7s+wBkSuLqbDxUREJBYo\nvKHTiVpERESincIbMB1PkSoiIhLlFN6gE7WIiEhMUXhDZJu31TZvERGJAQpvwHgTISVNS94iIhIT\nFN7tBl8ORz/o1elaRURE3KDwbucfAqda4MRxtysRERE5J4V3G6Od1kREJEYovNtFDhfTTmsiIhLd\nFN7tdKIWERGJEQrvNkanSBURkRih8G7XvtpcS94iIhLlFN7t0v3gOFidIlVERKKcwruNcRIgfbBW\nm4uISNRTeHfkHwKBOmyo1e1KREREzkrh3YEZfDmEQhA46nYpIiIiZ6Xw7qj9RC3a7i0iIlFM4d2R\nP7zHua4uJiIi0Uzh3YFOkSoiIrFA4d1R5BSpCm8REYleCu+O2k+RekSrzUVEJHopvDsalAreRC15\ni4hIVFN4d2CMCR/rrR3WREQkiim8Tzd4CNQfw55qcbsSERGRLim8TxO5uphWnYuISJTy9GRQS0sL\n999/P7fddhtjx45l8eLFhEIh0tPTmTNnDl6vl/LyctavX48xhqlTpzJlyhSCwSBLly6ltrYWx3Eo\nLS0lIyODvXv3snz5cowxZGVlMXv27P7us+c6Hi42NNPdWkRERLrQoyXvn/zkJwwaNAiA1atXM23a\nNB577DGGDRvGa6+9RlNTE2vXrmX+/Pk88sgjvPrqq5w4cYLXX38dn89HWVkZt956K6tWrQJgxYoV\nlJSUUFZWRmNjI9u2beu/DntrcPuJWrTkLSIi0anb8D5w4AD79+9nwoQJAFRVVVFUVARAUVERO3bs\nYNeuXWRnZ+Pz+UhMTCQvL4/q6moqKyspLi4GID8/n5qaGoLBIIcPHyYnJweAwsJCdu7c2V/99dqH\nJ2rRTmsiIhKdug3vF154gbvuuityv7m5Ga/XC0BqaiqBQIBAIEBqampkTFfTHcfBGEMgECA5OTky\nNi0tjaNHo+hCIH6dqEVERKLbObd5//rXv+bqq69m6NChffJi1toeTTubzMy+2wZ9tnmF0lM5AAxo\nrOfyPny9/taX/zbRQP1EN/UT3dRPdOuLfs4Z3m+99RaHDx/mrbfe4siRI3i9XpKSkmhpaSExMZG6\nujr8fj9+v59AIBB5Xl1dHbm5uZ2mB4NBrLWkp6dTX1/faazf7+9RsQcPHjyfHs+QmZl57nkNTKbp\nTwf67PX6W7f9xBj1E93UT3RTP9Gtt/2cLejPGd5z586N3F69ejVDhw6lpqaGiooKJk+eTEVFBQUF\nBeTm5vLcc8/R0NBAQkICNTU1lJSUcPLkyciYrVu3MmbMGDweD8OHD6e6upprrrmGzZs3M3369B43\nclEM1olaREQkevXoULGO7rjjDhYvXszGjRsZMmQIN9xwAx6Ph1mzZrFgwQKMMcycOROfz8ekSZPY\nsWMH8+fPx+v1UlpaCkBJSQnLli3DWktOTg7jxo3r88YuiH8IHPgj9mQjZqDP7WpEREQ66XF433HH\nHZHb8+fPP+PxiRMnMnHixE7T2o/tPt2IESN47LHHelPnRWUGX46F8LHew7PcLkdERKQTnWGtKzpc\nTEREopjCuyttp0i1RxXeIiISfRTeXTAdT5EqIiISZRTeXWk7RarCW0REopHCuyvtq821zVtERKKQ\nwrsLxuuFlDSdIlVERKKSwvtsBl8OR4/06vStIiIiF4PC+2wGD4FTLXDiuNuViIiIdKLwPgsT2WlN\n271FRCS6KLzPxq/DxUREJDopvM+mbcnbKrxFRCTKKLzPwugUqSIiEqUU3mfTvtpch4uJiEiUUXif\nTbofHEcnahERkaij8D4L4yRA+mVa8hYRkaij8D6XwUPgaB22tdXtSkRERCIU3udg/EPAhuBYndul\niIiIRCi8z0VXFxMRkSik8D6Xwbq6mIiIRB+F9zlEjvXWTmsiIhJFFN7notXmIiIShRTe5+LXKVJF\nRCT6KLzPZVAKeBN1ilQREYkqCu9zMMbAkAw4dBAbDLpdjoiICKDw7pa5Jh+aT8LuardLERERARTe\n3TJjCwGwlVtdrkRERCTM092A5uZmlixZwrFjxzh16hS33XYbFRUV7Nmzh5SUFABmzJjBddddR3l5\nOevXr8cYw9SpU5kyZQrBYJClS5dSW1uL4ziUlpaSkZHB3r17Wb58OcYYsrKymD17dr83e17yxoHH\nGw7v2+5yuxoREZHuw3vr1q1kZ2dz8803U1tby+OPP87VV1/NnXfeSWFhYWRcU1MTa9euZeHChXg8\nHubNm0dxcTFbtmzB5/NRVlbG9u3bWbVqFXPnzmXFihWUlJSQk5PDokWL2LZtGxMmTOjXZs+HGTAA\n8sZC1Tbs0SMY/2VulyQiInGu29XmkyZN4uabbwbgyJEjDB48uMtxu3btIjs7G5/PR2JiInl5eVRX\nV1NZWUlxcTEA+fn51NTUEAwGOXz4MDk5OQAUFhayc+fOvuqpz2nVuYiIRJMeb/N+6KGHWLRoESUl\nJQBs2LCBRx99lKeffprjx48TCARITU2NjE9NTSUQCHSa7jgOxhgCgQDJycmRsWlpaRw9erSPWup7\nH4b3Wy5XIiIi0oPV5u0ef/xx9u7dyzPPPMNdd91FSkoKo0aNYt26daxZs4a8vLwezcda26NpXcnM\nzOxpuX06L3vFFfxp2HBC1du5YuhQjKfH/2wXTV/+20QD9RPd1E90Uz/RrS/66TaF9uzZQ2pqKkOG\nDGHUqFG0traSlZVFWloaAEVFRTz//PNMnDiRQCAQeV5dXR25ubn4/f7I9GAwiLWW9PR06uvrO431\n+/3dFnvw4MFeN9iVzMzMXs8rNHo89rX1HCzfhMkb2yd19JXz6SeaqZ/opn6im/qJbr3t52xB3+1q\n83feeYdXXnkFgEAgQFNTE8uWLePQoUMAVFVVMXLkSHJzc9m9ezcNDQ00NTVRU1PD6NGjGT9+PBUV\nFUB457cxY8bg8XgYPnw41dXhY6c3b95MQUFBj5txg7Z7i4hItOh2yfumm27i2Wef5eGHH6alpYV7\n7rmHpKQknn76aRITE0lKSqK0tJTExERmzZrFggULMMYwc+ZMfD4fkyZNYseOHcyfPx+v10tpaSkA\nJSUlLFu2DGstOTk5jBs3rt+bvSA6ZExERKJEt+GdmJjIV77ylTOmL1y48IxpEydOZOLEiZ2mtR/b\nfboRI0bw2GOP9aZWV5kBA+DqsfCODhkTERF36QxrvWDyrwO06lxERNyl8O4FHTImIiLRQOHdGxnD\nw1cZe/dtXWVMRERco/DuBWNMeOn7ZCPs0VXGRETEHQrvXtIhYyIi4jaFd29dkw8eD3antnuLiIg7\nFN69ZAYkhQ8Z2/8HbOCI2+WIiEgcUnifB+11LiIiblJ4nwdt9xYRETcpvM/HsOFw2VB4Zzu2tdXt\nakREJM4ovM+DMQaTXwgnG2C3DhkTEZGLS+F9nrTqXERE3KLwPl/XjAsfMqbwFhGRi0zhfZ4ih4y9\n9wdsoM7tckREJI4ovC9AZNV5lQ4ZExGRi0fhfQHaw5udWnUuIiIXj8L7QrQdMmbffVuHjImIyEWj\n8L4AkUPGGhtgT43b5YiISJxQeF8gHTImIiIXm8L7QumQMRERucgU3hfIDEiC3DGwbw/22FG3yxER\nkTig8O4DusqYiIhcTArvPmDy2w4Z06pzERG5CBTefWHYiPAhY+9s0yFjIiLS7xTefcAYgxl7XfiQ\nMS19i4hIP1N49xEzeRokeAiteAYbOOJ2OSIicgnzdDegubmZJUuWcOzYMU6dOsVtt93GlVdeyeLF\niwmFQqSnpzNnzhy8Xi/l5eWsX78eYwxTp05lypQpBINBli5dSm1tLY7jUFpaSkZGBnv37mX58uUY\nY8jKymL27NkXo99+Y7KyMbd/Hvvi84SefxLnq49jEhLcLktERC5B3S55b926lezsbB599FHmzp3L\nCy+8wOrVq5k2bRqPPfYYw4YN47XXXqOpqYm1a9cyf/58HnnkEV599VVOnDjB66+/js/no6ysjFtv\nvZVVq1YBsGLFCkpKSigrK6OxsZFt27b1e7P9zUz5FKbwz+D/qrDrVrpdjoiIXKK6De9JkyZx8803\nA3DkyBEGDx5MVVUVRUVFABQVFbFjxw527dpFdnY2Pp+PxMRE8vLyqK6uprKykuLiYgDy8/Opqakh\nGAxy+PBhcnJyACgsLGTnzp391eNFY4zB3DUHhl6B3fAT7PY33S5JREQuQd2uNm/30EMPceTIER54\n4AHKysrwer0ApKamEggECAQCpKamRsZ3Nd1xHIwxBAIBkpOTI2PT0tI4erT7E5xkZmb2uLGLOa/T\ntcz/Lofv/zz8YBFD/3Ulnoz+e612/dmPG9RPdFM/0U39RLe+6KfH4f3444+zd+9ennnmGay15/2C\nXT23p/M7ePDgeb9uR5mZmX02ry4lDYK/nU3ohcX86bH7cb7+bUzbl53+0O/9XGTqJ7qpn+imfqJb\nb/s5W9B3u9p8z549fPDBBwCMGjWK1tZWBg4cSEtLCwB1dXX4/X78fj+BQCDyvK6mB4NBrLWkp6dT\nX19/xthLifnzv8RcPwX2/h675j/cLkdERC4h3Yb3O++8wyuvvAJAIBCgqamJ/Px8KioqAKioqKCg\noIDc3Fx2795NQ0MDTU1N1NTUMHr0aMaPHx8Zu3XrVsaMGYPH42H48OFUV1cDsHnzZgoKCvqrR1cY\nYzCzvgiZWdjXXiX0ZrnbJYmIyCWi29XmN910E88++ywPP/wwLS0t3HPPPWRnZ7N48WI2btzIkCFD\nuOGGG/B4PMyaNYsFCxZgjGHmzJn4fD4mTZrEjh07mD9/Pl6vl9LSUgBKSkpYtmwZ1lpycnIYN25c\nvzd7sZkBSThffIDQgvuxKxZjR16FGTbC7bJERCTGGXshG7AvspjZ5n2a0P/+Grv8uzD8Spx5T2IG\nDOjT+cf7NqFop36im/qJbvHez3lv85YL53zsBsyNfwUH/oj9z+fcLkdERGKcwvsiMXfcC1fmYH/7\nP4R+u9HtckREJIYpvC8S4/XifOHrMDAZ+6PnsPv/4HZJIiISoxTeF5G5fBjO3V+BUy2EnnoYu7va\n7ZJERCQGKbwvMlMwETPrS3CintB3H8Jued3tkkREJMYovF3g3PhXOHPmg5NA6N/+hdD6NRd01joR\nEYkvCm+XmPxCnAe+DYOHYP/rh9gVz2CDp9wuS0REYoDC20VmxFU4855s2wt9I6FFj2IbTrhdloiI\nRDmFt8tM+mCcf/oWFEyE6h2Evv11bO37bpclIiJRTOEdBcyAJJwvfQNz0y3w/n5CC/9Je6KLiMhZ\nKbyjhHEScG6/O7wnekM9oSf/WRczERGRLim8o0x4T/SHwePBLnuC0Po1bpckIiJRRuEdhczY63C+\n8Z3Inuih19a7XZKIiEQRhXeUMiNG4Xz925CShn1xGfbd7W6XJCIiUULhHcXMZUNxSueBcQg99x3s\noUvnsngiInL+FN5RzuRci/ncl6HxBKHFZdhGHQcuIhLvFN4xwPmzv2g7jOwAoWVPYFtb3S5JRERc\npPCOEea2uyC/CKq2Ydf+wO1yRETERQrvGGGcBJzZX4MrRmI3/pRQ+S/dLklERFyi8I4hZqAP5x8f\nguQU7I+ew/5fpdsliYiICxTeMcYMvQLnSw8AltCz39Z50EVE4pDCOwaZvHzMnV+AE8cJLVlAqLHB\n7ZJEROQiUnjHKGfydMyUT8GBP3LkiYewIe2BLiISLxTeMczccQ+MHk/T5nLsf610uxwREblIFN4x\nzCQk4HzhG3iGZ2E3/ITQm6+7XZKIiFwECu8YZ5IHMeThp8CbiP3xcmxTo9sliYhIP/P0ZNDKlSt5\n9913CYVC3HLLLWzZsoU9e/aQkpICwIwZM7juuusoLy9n/fr1GGOYOnUqU6ZMIRgMsnTpUmpra3Ec\nh9LSUjIyMti7dy/Lly/HGENWVhazZ8/u10YvZd4RozDTb8X+7EXs+jWYW+9yuyQREelH3YZ3ZWUl\n7733HgsWLKC+vp6vf/3rjB07ljvvvJPCwsLIuKamJtauXcvChQvxeDzMmzeP4uJitmzZgs/no6ys\njO3bt7Nq1Srmzp3LihUrKCkpIScnh0WLFrFt2zYmTJjQr81eysy027C/3Yj9759i//wvMUMz3S5J\nRET6Sberza+99lrmzp0LQHJyMs3NzYRCoTPG7dq1i+zsbHw+H4mJieTl5VFdXU1lZSXFxcUA5Ofn\nU1NTQzAY5PDhw+Tk5ABQWFjIzp07+7KvuGMGDMDMvBuCQUKr/8PtckREpB91u+TtOA5JSUkAbNq0\niQkTJuA4Dhs2bOCVV14hLS2Nu+++m0AgQGpqauR5qampBAKBTtMdx8EYQyAQIDk5OTI2LS2No0eP\ndltsZmbfLU325byiQWZmJnbG7dT+biPN2zcz+P19JF030e2yztul+P5cStRPdFM/0a0v+unRNm+A\nN998k02bNvHQQw+xe/duUlJSGDVqFOvWrWPNmjXk5eX1aD7W2h5N68rBg31zPevMzMw+m1c06NiP\n/fRdUPk2tUu/g/PwIoynx29x1LiU359LgfqJbuonuvW2n7MFfY/2Nn/77bd56aWXePDBB/H5fOTn\n5zNq1CgAioqK2LdvH36/n0AgEHlOXV0dfr+/0/RgMIi1lvT0dOrr688YKxfOjLwKM/km+NN72F+9\n6nY5IiLSD7oN78bGRlauXMkDDzzAoEGDAHjyySc5dOgQAFVVVYwcOZLc3Fx2795NQ0MDTU1N1NTU\nMHr0aMaPH09FRQUAW7duZcyYMXg8HoYPH051dTUAmzdvpqCgoL96jDvm5r8DXzL25Rex9cfcLkdE\nRPpYt+tU33jjDerr6/ne974XmXbjjTfy9NNPk5iYSFJSEqWlpSQmJjJr1iwWLFiAMYaZM2fi8/mY\nNGkSO3bsYP78+Xi9XkpLSwEoKSlh2bJlWGvJyclh3Lhx/ddlnDEpqZgZs7AvLsP+1w8xf/+Pbpck\nIiJ9yNiebnCOAtrm3bWu+rGtrYQe+wr86T2cf34Kc2W2S9X1Xjy8P7FM/UQ39RPdLuo2b4k9JiEB\n529ng7WEXlzW450CRUQk+im8L2Fm9HiYMBF2vYvd/Bu3yxERkT6i8L7EObffDR4vdu0PsM1Nbpcj\nIiJ9QOF9iTOXD8Pc9GkIHMH+fK3b5YiISB9QeMcB89czIf0y7C/+C1v7vtvliIjIBVJ4xwEzIAkz\nswSCpwit/b7b5YiIyAVSeMcJUzwZckbDW7/Dvrvd7XJEROQCKLzjhDEG52//AYwh9OPl2C6uDCci\nIrFB4R1HzJXZmImfgAN/xG79rdvliIjIeVJ4xxnzqc+A42B/9iI21Op2OSIich4U3nHGDL0Cc/0n\nwlcde/N1t8sREZHzoPCOQ+aTn4GEBOwrP9bSt4hIDFJ4xyFz+TDM9VPg/f1a+hYRiUEK7zhl/vr2\ntqVvbfsWEYk1Cu84ZS4fhpn0F/D+AV20REQkxii849iHS9+rsa1a+hYRiRUK7zhmhmRg/mwqHNLS\nt4hILFF4x7nw0rcnvO1bS98iIjFB4R3nzGVDw0vfh/+E/d9fuV2OiIj0gMJbOix9/1hL3yIiMUDh\nLZjLLsd8/C+h9n1sxa/cLkdERLqh8BYAzF/NBI8H++qPscGg2+WIiMg5KLwFADP4csyf39S29P2a\n2+WIiMg5KLwl4sOl79Va+hYRiWIKb4kwg4dgPj4NPjiE/d0mt8sREZGzUHhLJ+Glb2/b0vcpt8sR\nEZEueHoyaOXKlbz77ruEQiFuueUWsrOzWbx4MaFQiPT0dObMmYPX66W8vJz169djjGHq1KlMmTKF\nYDDI0qVLqa2txXEcSktLycjIYO/evSxfvhxjDFlZWcyePbu/e5UeMP7LMJOnYTe9gn1jE2byNLdL\nEhGR03S75F1ZWcl7773HggULePDBB/nBD37A6tWrmTZtGo899hjDhg3jtddeo6mpibVr1zJ//nwe\neeQRXn31VU6cOMHrr7+Oz+ejrKyMW2+9lVWrVgGwYsUKSkpKKCsro7GxkW3btvV7s9Iz5q9uA28i\ndv0aLX2LiEShbsP72muvZe7cuQAkJyfT3NxMVVUVRUVFABQVFbFjxw527dpFdnY2Pp+PxMRE8vLy\nqK6uprKykuLiYgDy8/Opqamlx2v5AAAZQ0lEQVQhGAxy+PBhcnJyACgsLGTnzp391aP0kkkPL31z\n5DD29Y1ulyMiIqfpdrW54zgkJSUBsGnTJiZMmMD27dvxer0ApKamEggECAQCpKamRp7X1XTHcTDG\nEAgESE5OjoxNS0vj6NGj3RabmZnZu+4u0ryiQV/301ryZf70+i9xNqxl2G2zcAYk9en8u6P3J7qp\nn+imfqJbX/TTo23eAG+++SabNm3ioYce4r777jvvF7TW9mhaVw4ePHjer9tRZmZmn80rGvRbP1M+\nRevPf8LB//x3nJs+3ffzPwu9P9FN/UQ39RPdetvP2YK+R3ubv/3227z00ks8+OCD+Hw+kpKSaGlp\nAaCurg6/34/f7ycQCESe09X0YDCItZb09HTq6+vPGCvRxUy7FQYmY3++Fnuy0e1yRESkTbfh3djY\nyMqVK3nggQcYNGgQEN52XVFRAUBFRQUFBQXk5uaye/duGhoaaGpqoqamhtGjRzN+/PjI2K1btzJm\nzBg8Hg/Dhw+nuroagM2bN1NQUNBfPcp5MskpmGmfhhP12P/+qdvliIhIm25Xm7/xxhvU19fzve99\nLzLty1/+Ms899xwbN25kyJAh3HDDDXg8HmbNmsWCBQswxjBz5kx8Ph+TJk1ix44dzJ8/H6/XS2lp\nKQAlJSUsW7YMay05OTmMGzeu/7qU82b+4m+w//Mz7H+vw075JGZQavdPEhGRfmVsTzc4RwFt8+5a\nf/cT2vhT7I//HTPt0zgzP99vr9NO7090Uz/RTf1Et4u6zVvim7nhr8A/BLvpVWzgiNvliIjEPYW3\ndMt4EzF/87dwqgX76mq3yxERiXsKb+kRc/0UGHoFtvyX2Nr33S5HRCSuKbylR4zHg5lxJ7S2Yn/2\nn26XIyIS1xTe0mPmox+H4VdiK36FPbDP7XJEROKWwlt6zDgOzqc/B9YSevlHbpcjIhK3FN7SO+M+\nClddDW/9Drv3925XIyISlxTe0ivGmPDSNxBat9LlakRE4pPCW3rNjB4Po8dD1TZsTaXb5YiIxB2F\nt5wX55a/AyC07oc9viqciIj0DYW3nBfzkTwYXwy73oXKt9wuR0Qkrii85bw5t8wCY8JL36GQ2+WI\niMQNhbecNzPiKsxHJ8O+PfDWG26XIyISNxTeckHMzZ8FxyH00gvYlma3yxERiQsKb7kgZmgmZuoM\nqH0fu36N2+WIiMQFhbdcMPM3nw1fMnTDS9g/7Xe7HBGRS57CWy6YSRqI89l/gNYgoR89q0PHRET6\nmcJb+kbBx8KHjtXsxFb8yu1qREQuaQpv6RPGmPDSd+IA7Jr/wDaccLskEZFLlsJb+oy5bCjmb/4W\n6o9hX3rB7XJERC5ZCm/pU2bqzZCZhf3NBuzuarfLERG5JCm8pU8Zjwdn1pcACK1cim1tdbkiEZFL\nj8Jb+py5egzmz6bC/r3Y//mZ2+WIiFxyFN7SL8xtJTAoBfvyKmxdrdvliIhcUhTe0i9MSipm5ueh\nuYnQj5e7XY6IyCVF4S39xlw/BXKvhbd+h93+ptvliIhcMjw9GbRv3z6eeOIJPvnJTzJ9+nSWLFnC\nnj17SElJAWDGjBlcd911lJeXs379eowxTJ06lSlTphAMBlm6dCm1tbU4jkNpaSkZGRns3buX5cuX\nY4whKyuL2bNn92ujcvEZx8GZVUqo7CuE/vPfcK4ZhxkwwO2yRERiXrfh3dTUxPe//33Gjh3bafqd\nd95JYWFhp3Fr165l4cKFeDwe5s2bR3FxMVu2bMHn81FWVsb27dtZtWoVc+fOZcWKFZSUlJCTk8Oi\nRYvYtm0bEyZM6PsOxVVmeBbmL2/BbvgJ9tUXMbfe5XZJIiIxr9vV5l6vl3nz5uH3+885bteuXWRn\nZ+Pz+UhMTCQvL4/q6moqKyspLi4GID8/n5qaGoLBIIcPHyYnJweAwsJCdu7c2QftSDQyn/oMXDYU\n+8t12AP73C5HRCTmdbvknZCQQEJCwhnTN2zYwCuvvEJaWhp33303gUCA1NTUyOOpqakEAoFO0x3H\nwRhDIBAgOTk5MjYtLY2jR492W2xmZmaPmuqJvpxXNIj2fk7+4zw+eHQu3jX/zuXf/jeMc+7vjdHe\nT2+pn+imfqKb+jlTj7Z5n27y5MmkpKQwatQo1q1bx5o1a8jLy+vRc7u64lRPr0J18ODBXtV5NpmZ\nmX02r2gQE/2MyIYJE2neVsGBZ5/E+fTfnXVoTPTTC+onuqmf6Bbv/Zwt6M9rb/P8/HxGjRoFQFFR\nEfv27cPv9xMIBCJj6urq8Pv9naYHg0GstaSnp1NfX3/GWLm0OX//j3D5MOz61YT+99dulyMiErPO\nK7yffPJJDh06BEBVVRUjR44kNzeX3bt309DQQFNTEzU1NYwePZrx48dTUVEBwNatWxkzZgwej4fh\nw4dTXR0+9/XmzZspKCjoo5YkWplBqThzHoaBydgf/Ct2T43bJYmIxKRuV5vv2bOHF154gdraWhIS\nEqioqGD69Ok8/fTTJCYmkpSURGlpKYmJicyaNYsFCxZgjGHmzJn4fD4mTZrEjh07mD9/Pl6vl9LS\nUgBKSkpYtmwZ1lpycnIYN25cvzcr7jNXjMD5wtcJLXqU0JIFOP/8Xczgy90uS0Qkphjb0w3OUUDb\nvLsWi/2E/udn2Befh5FX4XzjO5gBSZHHYrGfc1E/0U39RLd476dPt3mLXCgz5VOYydPgvT8Q+o/v\nYUMht0sSEYkZCm9xhTEG89kvQF5++PSpL69yuyQRkZih8BbXGI8H54vfCO+B/qr2QBcR6SmFt7gq\nvAf6fBjo0x7oIiI9pPAW15krRuL8wz9Bayuhpd8iWPu+2yWJiEQ1hbdEBTO2EHPH3XDsKB+U3Y9t\nbnK7JBGRqKXwlqhh/uJvMB+/iVO7awj9x9PaA11E5CwU3hI1jDGYO7/AgPzr4K03sN9/Gtt00u2y\nRESijsJboorxeLnswX+BK3OwFb8i9PhXsft2u12WiEhUUXhL1ElITcd54DuYmz4Nhw4QWvhPhDb+\ntMdXnxMRudQpvCUqGY8X5/bP43zlm+ELmfz43wk9U4atP+Z2aSIirlN4S1QzYwtxvvmvcG0B7NxC\n6NH7sO9ud7ssERFXKbwl6pk0P85XHsHMLIETxwl972FCL63ABoNulyYi4gqFt8QE4zg4027F+ca/\nwJAM7M9/QuhfHsDqhC4iEocU3hJTzFW5OPOfxnzsBvjD/xEq+3+EfvMLHVImInFF4S0xxwz0Ye75\nKubz/w9CIewPlxC6/3OE/u1fsG/9Dnuqxe0SRUT6lcftAkTOhzEGM2kK9uox2N/+D3bzb7BbXsdu\neR0G+jAFEzHFH4drxmM8+m8uIpcWfapJTDNDMjA334md8VnYtwf75m+wb5Zjf7cJ+7tNMCgVUzgJ\n89HJkHstxtHKJhGJfQpvuSQYY+DKbMyV2dhb74I9NeEQ3/I69tcbsL/eAGl+GHkVJmM4ZAzHZGRC\nxnDwX6ZQF5GYovCWS45xHMgZjckZjb3jHvi/ynCQ73gTKt/CVr4FQOR8bYmJMDQTMjLbgj0TMyQD\n0gZD+mDMgCTXehER6YrCWy5pJiEBRo/HjB4PgG1sgEMHsYcOwKGDcOgA9tDB8O39eyOB3ulErAN9\nHwZ5+uAzb6elQ6ofkzTwYrcnInFK4S1xxfiS4apczFW5naZba+FY3YfBXvcBBOqwx+og0Pbz/v5O\noX7GmdYHJEFqevgnzY9J9UNqOieyRmGtgZQ0SE6BlFQYmKxV9SJy3hTeIrRtM0+/DNIvw+TldznG\nnjoFx49GwtwG6uDYETgewB4LhB87HoA//F/4ELa25x3t8gUdGJQSDvNBqeEd61JSw/eTB8HA5PBe\n822/ifz2wYAkBb9InFN4i/SQ8XrhsqHhH8CcZZwNhaChHo6FwzzdgaN//AOcOA4njmNP1Educ+J4\neJW9DZ25JE8XS/cAxkBSW5AnDQwv8ScNhAEDMUlJ4fsDOkxPSoLEpPC2+wEDILF9TOfbxknoo38p\nEelvCm+RPmYcJ7yKPCUNgOTMTI4dPHjW8TYUgpMNUN8W5g0nsCcb4GRjeHqH37bTtMbwF4Tmk9B2\nnvdzXTS12wuqerzhIE8c0PaT2OH2AEzbtKP+wYROBcHbNsab2PZ7AKb9vrftuZHHOv54IcETXtsh\nIuelR+G9b98+nnjiCT75yU8yffp0PvjgAxYvXkwoFCI9PZ05c+bg9XopLy9n/fr1GGOYOnUqU6ZM\nIRgMsnTpUmpra3Ech9LSUjIyMti7dy/Lly/HGENWVhazZ8/u715FopJxnLbV5SnA8PC0Xs7DBk9B\ncxM0NYXDvLkJmk5C88nwqWObm8PTWpqhpanTfdvc1OGxtp+mxvBmgOYmCIXCr9H2WifOVUfPmw6H\neHuYtwe7x/vh/bbbxtNhnKfDYx5P2xcBL3g9bdPax3sj9yNjPW1jEjyRLxB4veEvTyIxptvwbmpq\n4vvf/z5jx46NTFu9ejXTpk3j+uuvZ9WqVbz22mtMnjyZtWvXsnDhQjweD/PmzaO4uJgtW7bg8/ko\nKytj+/btrFq1irlz57JixQpKSkrIyclh0aJFbNu2jQkTJvRrsyKXKtMeUskpZz52gfO2wWCnYL88\nPY3a/e/BqRZoaYFTzdiWlk73idxvDv8+dQpOtYS/ZETGtUD7/VOnwmsUTgUh2AKtrR++fm/r7eX4\n/QAJCeEw93jafns73O7wuz38PR5ISAj/u5/xXE/btA6PdRrz4TTTcZpz2tiEhA4/HaY7CZ0fcxK0\nD0Qc6ja8vV4v8+bNY926dZFpVVVVkSXloqIiXn75ZTIzM8nOzsbn8wGQl5dHdXU1lZWVTJ48GYD8\n/HyeffZZgsEghw8fJicnB4DCwkJ27typ8BaJQqZ9qdWXDEBiZiYm0dd5TB+/pg21fhjkbcFPMNgW\n9G23g6ciP/bUqU73w89tvx8847dt/fD+AMehubEhfL812Pl308m2223jT1tK7+0XhTP6vMDnRxgH\nEhxI8LDf48Eap1O4dwr8c912nPCXgYQPb384zuk0pvNznM73T5tfZIyTgElwTpvexbiEhHBPjsOp\nYDP2gw86zePD207X041zyW+W6Ta8ExISSEjovCNLc3MzXq8XgNTUVAKBAIFAgNTU1MiYrqY7Tvgf\nNBAIkJycHBmblpbG0aNd7pMrInHIOAkwICG8U11Pxl/Aaw3NzOTgOfZJ6MiGWsNrBSIh39r2+1SH\n2x2+ALQGPxzf2hpei3HatMjYjvNube3wc9o8Qh0eC505zpPgcKqp6cPHQ63hLzynjw+1nvFlJNLn\nBfx7dvtv2Mvx533RX9NFwJ9rWnfjI+NM5y8MHe6b5BTMzM+HD0ntZxd9hzVrz3zruprWlczMzD6r\noy/nFQ3UT3RTP9HtUuunp2wodMYXgvbbkd/t09u/RITapoU6PvbhFwPbGl5D8eG0UHhekeeFoDX4\n4XNDocgXich82+rq+nYIGwpCayj8nPYvIe01h06bHmrFto2NzKP99UK2w+1TXT+3wxi6ySrr8ZDx\n2XvwdvP/qS/+v51XeCclJdHS0kJiYiJ1dXX4/X78fj+BQCAypq6ujtzc3E7Tg8Eg1lrS09Opr6/v\nNNbv93f7uj39dtydzF58044F6ie6qZ/opn66075qPbEP59lzF+v9MWe53ZG1FmwoEvJY2+F2CDxe\nar0D4Rz19rafswX9ee3lkJ+fT0VFBQAVFRUUFBSQm5vL7t27aWhooKmpiZqaGkaPHs348eMjY7du\n3cqYMWPweDwMHz6c6upqADZv3kxBQcH5lCIiInJRGGMwTnhHRZM4ADMgCTPQh0kehBmUelFPkdzt\nkveePXt44YUXqK2tJSEhgYqKCu677z6WLFnCxo0bGTJkCDfccAMej4dZs2axYMECjDHMnDkTn8/H\npEmT2LFjB/Pnz8fr9VJaWgpASUkJy5Ytw1pLTk4O48aN6/dmRURELgXG9nSDcxTQavOuqZ/opn6i\nm/qJbvHeT5+uNhcRERH3KLxFRERijMJbREQkxii8RUREYozCW0REJMYovEVERGKMwltERCTGKLxF\nRERijMJbREQkxsTUGdZERERES94iIiIxR+EtIiISYxTeIiIiMUbhLSIiEmMU3iIiIjFG4S0iIhJj\nPG4XcDH94Ac/4Pe//z3GGEpKSsjJyXG7pPOycuVK3n33XUKhELfccgtbtmxhz549pKSkADBjxgyu\nu+46l6vsXlVVFU899RQjR44EICsrixkzZrB48WJCoRDp6enMmTMHr9frcqU9t2nTJn7zm99E7u/e\nvZvs7Gyam5sZMGAAAH//93/PRz7yEbdK7JF9+/bxxBNP8MlPfpLp06fzwQcfdPm+lJeXs379eowx\nTJ06lSlTprhdepe66ufZZ58lGAzi8XiYM2cO6enpfPaznyUvLy/yvIcffhjHib5lnNP7WbJkSZef\nAbH6/jz11FMcP34cgBMnTpCbm8unP/1pvva1r0X+dlJTU/nqV7/qZtlndfpndHZ2dt///dg4UVVV\nZRcuXGittfa9996zDz74oMsVnZ+dO3fab33rW9Zaa48fP26/+MUv2sWLF9stW7a4XFnvVVZW2ief\nfLLTtCVLltg33njDWmvtj370I/uLX/zCjdL6RFVVlX3++eftN7/5TfvHP/7R7XJ67OTJk/aRRx6x\nzz33nP35z39ure36fTl58qS97777bENDg21ubrZf/epXbX19vZuld6mrfp555hn729/+1lpr7c9/\n/nP7wx/+0Fpr7d133+1anT3VVT9dfQbE8vvT0ZIlS+zvf/97e+jQIfuNb3zDhQp7p6vP6P74+4m+\nr5T9ZOfOnXz0ox8FYMSIETQ0NNDY2OhyVb137bXXMnfuXACSk5Npbm4mFAq5XFXfqaqqoqioCICi\noiJ27NjhckXnb+3atcycOdPtMnrN6/Uyb948/H5/ZFpX78uuXbvIzs7G5/ORmJhIXl4e1dXVbpV9\nVl31c++99zJx4kQgvARXX1/vVnm91lU/XYnl96fdwYMHaWxsjKm1pF19RvfH30/crDYPBAKdVlWm\npqYSCATw+XwuVtV7juOQlJQEhFfRTpgwAcdx2LBhA6+88gppaWncfffdpKamulxpz+zfv5/vfOc7\nnDhxgttvv53m5ubIavL29ygW7dq1i8suu4z09HQAVq9eTX19PcOHD6ekpITExESXKzy7hIQEEhIS\nOk3r6n0JBAKd/p9F6/vVVT/tf0OhUIhf/OIXkS9ZLS0tLFq0iA8++ICPfexjfOpTn7ro9Xanq36A\nMz4DYvn9abd+/XqmT58euR8IBPjud7/L0aNHmTZtGh//+McvVpk91tVn9Pbt2/v87yduwvt0NsbP\nCvvmm2+yadMmHnroIXbv3k1KSgqjRo1i3bp1rFmzhnvuucftErt1xRVXcPvtt3P99ddz6NAhHn30\nUVpbW90uq09s2rSJG2+8EYC//uu/Jisri2HDhvH888+zYcMGZsyY4W6BQigU4plnnmHs2LHk5+cD\n8LnPfY7JkycD8M1vfpPRo0eTnZ3tZpk9Mnny5DM+Azpuu49FwWCQ6upq7r33XgBSUlL4zGc+w8c/\n/nEaGxt58MEHGTt2bLdrINzS8TP6vvvu6/P5x81qc7/f3+lbzdGjR6P2Te/O22+/zUsvvcSDDz6I\nz+cjPz+fUaNGAeFVMvv27XO3wB4aPHgwkyZNwhjDsGHDSE9Pp6GhgZaWFgDq6upi9j2qqqqKfHgW\nFxczbNgwAAoLC2Pm/ekoKSnpjPfl9L+pWHu/li5dGvkC2e6mm24iKSmJpKQk8vPzY+a96uozINbf\nn3feeafT6vKBAwfyiU98Ao/HQ2pqKh/5yEc4cOCAixWe3emf0f3x9xM34T1+/HgqKioA2LNnD36/\nn4EDB7pcVe81NjaycuVKHnjgAQYNGgTAk08+yaFDh4BwaLTvvR3tysvLefnll4Hw6rBjx45x4403\nRt6niooKCgoK3CzxvNTV1ZGUlITH48FaS1lZGQ0NDUD4AykrK8vlCnsvPz//jPclNzeX3bt309DQ\nQFNTEzU1NYwePdrlSnumvLwcj8fDHXfcEZl28OBBFi1ahLWW1tZWampqYuZvqavPgFh+fyC86enK\nK6+M3K+srGTFihUANDU18cc//pHMzEy3yjurrj6j++PvJ66uKvajH/2Id999F2MM99xzT+SbaizZ\nuHEja9as4YorrohMu/HGG/nFL35BYmIiSUlJlJaWkpaW5mKVPXPy5EkWLVpEY2MjwWCQmTNnctVV\nV7F48WJOnTrFkCFDKC0txeOJra07e/bs4cUXX+TBBx8E4I033uCnP/0pSUlJ+P1+vvSlL0UOG4tG\ne/bs4YUXXqC2tpaEhAQGDx7Mfffdx5IlS854XyoqKnj55ZcxxjB9+vSo3AbZVT/Hjh0jMTEx8gV+\nxIgR3HvvvaxcuZKqqiqMMRQVFXHrrbe6XP2Zuupn+vTp/PSnPz3jMyBW35+vfe1rrF69mmuuuYZJ\nkyYB0NraynPPPcfBgwcJhULcdNNNfOITn3C5+jN19Rn95S9/meeee65P/37iKrxFREQuBXGz2lxE\nRORSofAWERGJMQpvERGRGKPwFhERiTEKbxERkRij8BYREYkxCm8REZEYo/AWERGJMf8fNbjDOVS5\ndOkAAAAASUVORK5CYII=\n",
            "text/plain": [
              "<matplotlib.figure.Figure at 0x7fc429791198>"
            ]
          },
          "metadata": {
            "tags": []
          }
        }
      ]
    },
    {
      "metadata": {
        "id": "TqTuhkTsP90b",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "We can use the Edward2 model to sample data for the inferred values for $\\mathbf{W}$ and $\\mathbf{Z}$, and compare to the actual dataset we conditioned on. "
      ]
    },
    {
      "metadata": {
        "id": "T3O6PHe3XX8a",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 402
        },
        "outputId": "56daaadd-1752-482a-f108-fdb3e22f6520"
      },
      "cell_type": "code",
      "source": [
        "print(\"MAP-estimated axes:\")\n",
        "print(w_inferred_map)\n",
        "\n",
        "def replace_latents(w=actual_w, z=actual_z):\n",
        "\n",
        "  def interceptor(rv_constructor, *rv_args, **rv_kwargs):\n",
        "    \"\"\"Replaces the priors with actual values to generate samples from.\"\"\"\n",
        "    name = rv_kwargs.pop(\"name\")\n",
        "    if name == \"w\":\n",
        "      rv_kwargs[\"value\"] = w\n",
        "    elif name == \"z\":\n",
        "      rv_kwargs[\"value\"] = z\n",
        "    return rv_constructor(*rv_args, **rv_kwargs)\n",
        "\n",
        "  return interceptor\n",
        "\n",
        "with ed.interception(replace_latents(w_inferred_map, z_inferred_map)):\n",
        "  generate = probabilistic_pca(\n",
        "      data_dim=data_dim, latent_dim=latent_dim,\n",
        "      num_datapoints=num_datapoints, stddv_datapoints=stddv_datapoints)\n",
        "\n",
        "with tf.Session() as sess:\n",
        "  x_generated, _ = sess.run(generate)\n",
        "\n",
        "plt.scatter(x_train[0, :], x_train[1, :], color='blue', alpha=0.1, label='Actual data')\n",
        "plt.scatter(x_generated[0, :], x_generated[1, :], color='red', alpha=0.1, label='Simulated data (MAP)')\n",
        "plt.legend()\n",
        "plt.axis([-20, 20, -20, 20])\n",
        "plt.show()"
      ],
      "execution_count": 9,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "MAP-estimated axes:\n",
            "[[-0.20085734]\n",
            " [ 2.7221382 ]]\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeoAAAFOCAYAAABJ4F7CAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3Xl0HNWB/v1vVXVXd2tpqbXa8m4c\nTDysxiEsmRB+5hASJkw2M5nMCTAQcmLCmkBYYgYIQxZ2xxCSwAyEYUgC2WDIwLzENluMB4YlgG2W\nYMCrpLb23ruW94+2FRvvUksqSc/nHB9L1V1V96osP31v3brX8H3fR0RERALJHOkCiIiIyO4pqEVE\nRAJMQS0iIhJgCmoREZEAU1CLiIgEmIJaREQkwEKD2fn+++9nzZo1eJ7HZz/7WQ444ABuv/12PM+j\ntraW888/n3A4XK6yioiIjDvGQJ+jfv311/mv//ovrrjiCvr6+vj2t7/NIYccwhFHHMExxxzDAw88\nQENDAyeddFK5yywiIjJuDLjre86cOVx88cUAVFZWks/nWbVqFfPmzQNg3rx5vPrqq+UppYiIyDg1\n4KA2TZNoNArAsmXLOOKII8jn8/1d3fF4nO7u7vKUUkREZJwa9GCyF154gWXLlnH22WeXozwiIiKy\nnUENJnvllVf47W9/y3e+8x0qKiqIRqMUCgVs26azs5NEIrHXY2zatGkwRQi0lpYW1W8UU/1Gr7Fc\nN1D9RruWlpb9ev+AW9SZTIb777+fyy+/nKqqKgAOOeQQVq5cCcDKlSs5/PDDB3p4ERERYRAt6hUr\nVtDX18ett97av+0b3/gGP/nJT/jjH/9IQ0MDxx9/fFkKKSIiMl4NOKhPPPFETjzxxJ22X3XVVYMq\nkIiIiPyVZiYTEREJMAW1iIhIgCmoRUREAkxBLSIiEmAKahGRMWTp0qWceOKJ9PT07PF977zzDuvX\nr9/v43/pS18im83u9vU777yTxx9/fLev//nPf6arq2u/zzueKahFREaI60KhUPq7XJYuXUpLSwtP\nPfXUHt/3zDPPsGHDhvKdeB899thjml56Pw1qZjIREdl/vg/JpEE+b+L7PoZhEIl4NDb6GMbAj9vb\n28sbb7zBpZdeyi9/+UtOPfVUAN5++21uu+02DMPg4IMP5qSTTuKRRx6htraW2tparr32Wu655x5i\nsRh33nknM2bM4G//9m/513/9V3K5HLlcjgsuuIAPf/jDuzzvE088wS9+8QsaGxuJRCLMmDGDdDq9\n0/7pdJpnn32W9957j2uvvZannnqKp556Cs/zOProoznjjDMGXvkxTC1qEZFhlkwaOI5BOOxj2xAO\n+ziOQTI5iJQGnnzySY4++miOOuooNm7cSDKZBGDJkiV885vf5Pbbb6ezs5OKigqOOuoozjnnnN2G\nb2dnJ6eccgq33nor55xzDr/4xS92+T7f97n77ru5+eabuf7669m4ceNu9583bx6zZs3isssuo7m5\nGYAf/ehH/PjHP+bxxx8nnU4Pqv5jlVrUIiLDyHUhnzcJh/0dtptmabvruljWwI69dOlSTj/9dCzL\n4uMf/zjLly/ntNNOY/369RxwwAEAXHnllft0rEQiwX/8x3/wq1/9imKx2L9a4gf19vYSi8X613Y4\n+OCD93n/SCTCRRddhGVZ9PT00NfXR2Vl5UCqPqYpqEVEhpHrllqhu+L7Pq7LgII6mUyyZs0afvzj\nH2MYBvl8nqqqKk477TSMvfSnb/+64zgA/OY3v6GhoYErr7ySN998kzvvvHO3ZTbNv3bOep63T/u3\ntrby0EMPcddddxGLxfjnf/7n/a/0OKGubxGRYWRZ7DY4DcMYVGv6s5/9LP/2b//G3XffzX333Udv\nby8bN25k+vTprF69GoAbbriB999/H8MwcLeOYquoqKCjowPXdfvf19PT07/K0zPPPNMf4B8Uj8dJ\np9OkUikcx+H111/f4/6mWeo16Onpoba2llgsxltvvUVbWxvFYnFglR/jFNQiIsPIsiAS8dja8Ozn\neaXtAw3qZcuW8alPfar/e8Mw+OQnP8ny5cs577zzuPPOOzn//POprq5m2rRpHHrooSxZsoQXX3yR\nz33uc1x55ZX8y7/8C9OnTwfgpJNO4qGHHuLSSy/lwx/+MJ2dnTz22GM7ndc0Tc444wwuvPBCrr76\nambMmLHH/Q877DCuvvpqbNsmFotx3nnnsXz5cj7zmc9w2223DazyY5zh764PZpiM9TVHVb/RS/Ub\nvYJet8GO+g56/QZrPNRvf+getYjIMDMMaGrycV23/570QFvSMvYpqEVERogCWvaF7lGLiIgEmIJa\nREQkwBTUIiIiAaagFhERCTAFtYjIGPG73/2Oc889lwsvvJCFCxfy4osvAnD77bezefPmAR/3oosu\n4t13393t63/605/2ebKSd999l4suumiP7/noRz+6x9f3tjLYB33ve99j1apV3HvvvZx++uk7leeE\nE07glVde6d+2atUqTjjhBP7yl7/0b3v88cc57bTTuOiii7jwwgu5+OKLee+990gmk1x66aW7nRCm\nHBTUIiIjpYzrXLa2tvKHP/yBH/3oRyxevJjvfOc73HfffQCcd955TJw4cdDn2J2HHnpoSINqe8Vi\nkYceemif3//cc88RiUT4m7/5GwBc1+Xtt9/uf33ZsmU7Pde8dOlSpkyZwrJly3bYfsIJJ3Dbbbex\nePFizjjjDJYsWUJjYyNHHXUUv/nNbwZRqz3T41kiIsPN9zGSScx8vn/CEy8SwW9sZKDrXKZSKQqF\nAsVikVAoxOTJk1m8eDFAfyvwqaeeoqenh40bN7J582bOOussHnvsMVpbW/nBD35Ae3s7v/vd77j2\n2msB+Pu//3sefvjh/nMkk0m+973vAaU5wS+//HJWrVrF6tWrueyyy7j55pt59NFHWbp0KaZp8rGP\nfYzTTjuNZDLJNddcQzgc7l8cZHuu63LdddeRTCaZPXt2//YXX3yRf//3fycUClFdXc3VV1/NHXfc\nwdq1a7n11lv52te+ttelOH/961+zcOHC/u8/+tGPsnTpUj70oQ8B8Pzzz++wj+u6PP3001x11VX8\n4Ac/4Gtf+9ouf95z5szpXyns7/7u7/jqV7/KP/zDP+z7BdsPalGLiAwzI5nEcBz8cBhsGz8cxnAc\njK3LUg7ErFmzOOigg/jyl7/MD37wA5YvX94/l/f2ent7ueGGGzj++OP5n//5H2644QY+8YlPsGLF\nir2eo6Ojg9NPP51bb72VT33qUzz88MOcdNJJ1NXV8cMf/pAtW7bw9NNPs2TJEhYvXszTTz9NW1sb\nv/nNb/pbow0NDTsd94UXXsB1Xe644w5OPPFEuru7Aejr62PRokUsXryYiooKXnjhBb70pS8xZcoU\nLr744r0uxek4Du++++4OHw6OOuoonnvuOXzf580332Tq1KmEQn9ts7744otMnTqVww47jHg8zqpV\nq3b5s3juuec46KCDAPpXD9uwYcNef4YDoRa1iMhwct1SSzoc3nG7aWLm86VwHeAsKFdeeSXvv/8+\nL7zwAr/85S955JFHuOWWW3Z4z7bWY319ff/iIIlEgt7e3r0ev66ujiVLlnDvvffS19fHgQceuMPr\na9asYcOGDVx88cUAZDIZWltbef/99/nEJz4BwGGHHcb//u//7rDf+++/3981PWfOnP4lMWtra7nx\nxhtxXZfNmzczd+7cHfbb21Kavb29xOPxHRZBiUajzJgxg9dee40//elPHH/88Tz77LP9ry9dupT5\n8+cDMH/+fJYuXdpftuXLl/Pmm2/i+z719fWcd955/fs1NDTQ3t7O5MmT9/pz3F8KahGR4eS6e1zm\ncqDrXPq+T7FYZNq0aUybNo3Pfe5znHHGGbS1te3wPmu7Y2//9a7K9MEW+T333MNHPvIRTj31VJ56\n6imee+65HV4Ph8McffTRfOtb39ph+y9+8Yv+sNzVeXa3VOYNN9zA97//faZNm9bfjb+9fV2K84M+\n8YlP8OSTT/Lyyy9z9tln9wd1oVBgxYoVvPXWW/zud7/DcRxSqVR/IJ9wwgk7dKMPF3V9i4gMJ8va\n4zKXA21N//d//zc33XRTfxCm02l83yeRSOzzMSorK+no6ADgnXfeIZPJ7PD6tqUrfd/fYaT3tqUr\nDzzwQF555RVyuRy+77NkyRLy+TxTpkzhzTffBODll1/e6bzbv/76669TKBT669DU1EQqleLll1+m\nWCzusDzn3pbijMfj9Pb27vTh4Oijj+bZZ59l+vTp2Lbdv33FihUcccQR3HPPPdx9993ce++9TJ06\ndZdl/qAtW7bQ2Ni41/cNhIJaRGQ4WRZeJMKu1rn0IpEBB/XJJ59MIpHg3HPP5Zvf/CaLFi3ivPPO\nIxKJ7PMxDjjggP6lJ5944gkmTJiww+uf+cxn+NGPfsRll13GCSecwKuvvsoLL7zAYYcdxvnnn080\nGuULX/gCF154Ieeeey51dXVEIhG+8IUv8Nhjj3HppZeSSqV2Ou9HP/pR8vk8F154IcuWLaO5uRko\nDWY7//zzuemmm/jSl77EAw88gGEYOI7DNddcs9elOEOhENOnT2ft2rU7nC8ajTJnzhyOP/74HbYv\nXbp0h6VCt/1cly9fvsefWy6Xo6uriylTpuz9hzwAWuZyCI2HpdpUv9FrLNcv8HUb5KjvwNdvkMpZ\nvz/96U88//zz/ffNh8Kvf/1risUi//iP/7hP79/fZS7VohYRGW6Ggd/UhNvSgjdhAm5LC35T04Af\nzZLdO+6448hms6xevXpIjp9MJlm5ciVf/OIXh+T4MMjBZOvWrePGG2/klFNO4eSTT+5/vq26uhqA\nU089dadReiIispXWuRwWV1555ZAdu7GxkZtuumnIjg+DCOpcLsc999zDwQcfvMP2L3/5yxx55JGD\nLpiIiIgMous7HA5zxRVX7NeIQhEREdk/gx5M9uCDDxKPx/u7vru7u3Ech5qaGs466yzi8Xi5yioi\nIjLulHXCk49//ONUV1czffp0fv/73/PQQw9x9tln73EfjVwcvVS/0W0s128s1w1Uv9FuREd9H3LI\nIUyfPh2AefPmsW7dunIeXkREZNwpa1DfdNNN/dPVrVq1asge/hYRERkvBtz1vXbtWu677z6SySSW\nZbFy5UpOPvlkbrvtNmzbJhqNcu6555azrCIiIuPOgIN65syZXHPNNTttP/roowdTHhEREdmOZiYT\nEREJMAW1iIhIgCmoRUREAkxBLSIiEmAKahERkQBTUIuIiASYglpERCTAFNQiIiIBpqAWEREJMAW1\niIhIgCmoRUREAkxBLSIiEmAKahERkQBTUIuIiASYglpERCTAFNQiIiIBpqAWEREJMAW1iIhIgCmo\nRUREAkxBLSIiEmAKahERkQBTUIuIiASYglpERCTAFNQiIiIBpqAWEREJMAW1iIhIgCmoRUREAkxB\nLSIiEmChwey8bt06brzxRk455RROPvlktmzZwu23347nedTW1nL++ecTDofLVVYREZFxZ8At6lwu\nxz333MPBBx/cv+3BBx/kk5/8JN/97neZMGECy5cvL0shRURExqsBB3U4HOaKK64gkUj0b1u1ahXz\n5s0DYN68ebz66quDL6GIiMg4NuCub8uysCxrh235fL6/qzsej9Pd3T240omIiIxzg7pHXQ4tLS0j\nXYQhpfqNbqrf6DWW6waq33hS1qCORqMUCgVs26azs3OHbvHd2bRpUzmLECgtLS2q3yim+o1eY7lu\noPqNdvv7IaSsj2cdcsghrFy5EoCVK1dy+OGHl/PwIiIi486AW9Rr167lvvvuI5lMYlkWK1eu5IIL\nLuCOO+7gj3/8Iw0NDRx//PHlLKuIiMi4M+CgnjlzJtdcc81O26+66qrBlEdERES2o5nJREREAkxB\nLSIiEmAKahERkQBTUIuIiASYglpERCTAFNQiIiIBpqAWEREJMAW1iIhIgCmoRUREAkxBLSIiEmAK\nahERkQBTUIuIiASYglpERCTAFNQiIiIBpqAWEREJMAW1iIhIgCmoRUREAkxBLSIiEmAKahERkQBT\nUIuIiASYglpkPHJd3GyBQtbFdUe6MCKyJ6GRLoCIDCPfh/YknZsKuHGHzi1JvEiM0MR6GpvAMEa6\ngCLyQWpRi4wjRjJJZ5uLa9qEKyOEK2wiRhG/fQvJpFJaJIgU1CLjheviZ/MUnBDm9r/5pkmomCOf\n8dUNLhJACmqR8cJ1cR0f3/d3esnwAddRUIsEkIJaZLywLKyQgbGLG9G+AVghLGv4iyUie6agFhkv\nLAsjFsEOOXjedts9DyccJVJhKKhFAkhBLTKO+I2N1DVbWF6BYjpPMVMg74cxmhpobNy5S1xERl5Z\nH89atWoVt9xyC1OmTAFg6tSpnHXWWeU8hYgMhmFAcxN1DS7NDc24m0JYtqWWtEiAlf056jlz5vCt\nb32r3IcVkXJw3dIfy8KK2dgxa7evK71FgkETnoiMB76P0d6OmU7jGwZGKAS2XZoAxTBKryeTmPk8\nvu9jGAZeJILf2KhZUERGWNmDesOGDfzwhz8klUqxYMECDj300HKfQmTcGlCD1/cx33gDa2tIm4aB\nZ9tQKGB0dOA3NWEkkxiOgx8Ol3YBDMeBZBK/qWnI6iMie2f4u3qocoA6Ozt54403OOaYY2hra+Pa\na69lyZIlhEJquIsMhu9Deztks39tBMdi0LSHaT/7Q33Teqx175Z2APC80k6RCNTUQEsLhfc3UcQm\nHC41tPsVCjBlirrBRUZQWRO0rq6OY489FoAJEyZQW1tLZ2cnTXv4RL5p06ZyFiFQWlpaVL9RLEj1\na283cBxjhxnFPA82b/Zpatrxs7bvQ3JjkWKqiNHZQax9A5F0Jwk7C/j48TiGaVLf2Mim3hRvvuZQ\n2NSNa9mYJlRWehxwgIfvg5crgAtWzGY0CdK1Gwqq3+jW0tKyX+8va1A/88wzdHV1ceqpp9Ld3U1P\nTw91dXXlPIXIuOO6kM+bhMM7BrJplra7rouFi1twcT2Dvpfewe7uprJnCzgOZtHF9316CyY1cR+y\nWfzaWigW2fByJ5kp04jZpeOZJuRyJi+9ZNDSAhRNipZNpMKgsdHX7WqREVDWoJ43bx6LFy/m//7v\n/3Ach69+9avq9hYZJNdll9N+Avieh9faTk97gWLWw1j9OulN3dh1CRrMPnwTzO4u/FyBXPNUojET\ns1CAWIGMH2VLp0VP3CDdFcZzXHzDopB1sKMGfrNPqCoKUZOtt6t3ar2LyNAra4rGYjEuv/zych5S\nZNyzLLZO+/nXkHTdUte31bGFHsPB6OzGePUNYq+8gJk1cUNheqoj1MxsJteege4Okl3V9NXYRPJ9\n9Bh5vIoQGzemcf7yZzKxBNN4H9Nx6XXiFC0TvymGO+kg4AOtd92uFhlWau6KBJxlQSTi4TgGhgGd\nnQb5PHhFj8quHOGeDmrWr8FsbcdL5ykWIti5NnKtId7vNLEn1lO1+S9YqY3EOvN09oZoi0bJt9SQ\n2ryeUB6qihvpTVRSlQjhhivoijXS6bk0d23BrS+NMfF9v3/EuYgMHwW1yCjQWOewpc1j42abomdi\nWQaxUI4JHW+QeXM9zhtvguMR3bSBEBWQyxLDxezuIdfZBNkMOTdPKJWkWAwRCxWp2/g2aXcSG71D\naLYKbHFqaTJ93FwH7pQmCgUfp6MPr6IWK2ZjGJoLXGQkKKhFgmzrRCRWPk+D4+MULfxoBL+hkdhf\n3iLS3UpmSw9WTwq7K4nd04rpWVD0cQ2wMIilkxQra5jc838Y2QwFr5G4t5GKUI5k1MRIpchWNdFT\nWU+sIkSMNLGuVjpf96hqzuN2+ZjVFdTObsRSUosMOwW1SIBtPxGJBxA2CXkOzvp1hNa9h5Er4PRk\nMdvaMTLd5PIWMSdDmBwGLh4GhUIEI5UCitjkaSZLAZuQCxPyb5GPHoBphai1O2hpqCPW18HGYjWp\nvkqsmhRmn0Uo202EVoxwi2YrExlmCmqRoHLd0pSeW2cLM03o6oRQTw/hjlbst9+h06vAcXzCRY+Q\nlwd8wmQwcLAAE5cwKWyigIFFkQrAxcQAXApEcr20eTMJUaDh1RwVfpbIhzyMLWka6qO4B9ZjhQ2c\nnIOXzWNqtjKRYaWgFgkq193hsayuLoNYphPPc0gXQoRT0J6L0NtTzxS7lrrsJmynjyhpHAx8TMIU\nCOERJYVJady41/+3SS0+Ng65QoxYqkCBStxomOI76whPmUjGrSC+6V2caTPB8/EcCOXzuBpVJjJs\ntB61SFBZ1tbHskqPYxWyPolYlkzWpODZ9PgV4DkY0TBgkTMrAChi4WETwsHEwwBK7yh9Mre3fm/i\nEaNAiCJNbKSQzZPxK3CicexCjmymACGLbHsGCg6YBqZtlj48uO6I/EhExiMFtUhQWRZeJAKeh+dR\namF7PrGIT9WECopTDyASDxNyi/jZDOm8RQ4Thwje1pZyhF13m5mUAtvCJUaOWrppZjOhvh6KDuAV\niXa04hULFPM+ha4UVjSEFd764UGtaZFho65vkQDzGxshmcQq5jE8H99x8SyLYk09TqNBX69JZ+Y9\nas0ohldBLT4h0kS23qPeE4NSWBvkcICprCWe78V9N4ZfVU1bbDLppWnClRHcWB3VsUq2rElSN6dJ\nQS0yjBTUIkFmGKWBW66LGfLIRSMU2w0iqQ4yxSKFth5q05vxcy4FwviAhUeIUhDvi23d4iYFbDbh\nUEE21QupFK2RuXQVKym8k6W3yqIiY9Fpmxywh1W7RKS81PUtMhpYFg0tYUKTmglnu/DSWfq6XfK9\nefqKMZy8RzObsLZOM7q/GWoCka1/bDKESRFnC0Z7Kz3FavraciTdBnqjzeR7HNo3e2WuoIjsjlrU\nIqOEYUBTg0vjERVsbqsm3Zoj2+1QkeqjSIhKMlSSIoQ34E/g2zq0K8niUCBRaCeUfIUeZzJvvHog\n0cRUKkI+mT5NJyoyXNSiFhlNXBfL9ME0CaV6aUqvZZKznim8TzXdxMgTGeQpto0OD+NSR5JQppto\n9ybcN98m/ZckPgZYpgZ+iwwTBbXIaGJZeL5B5xsd2MlN1Bfaqc9tppFNRMnvdQDZPp+GUlhX0Uu9\nv5lwNkN0Swc9z75Krr0HK2ypNS0yTBTUIqOIb5hsXtWN88dnKL7+DoVkLwk2UkkWA79sQQ2l/xyi\nFKmlk7i7hXhhC4WswaauSvIb2hXUIsNEQS0ySviez7rH3qT9+U2k1vUR6+2g3k8SY9vc3kNwTsDB\nwsenKb8BPCi4JvnuPG5Bfd8iw0FBLTJKJF9rI5dM4Xo+Xk+GWG8rPgXC5DEpPWZVbial+cItPKJG\nlhqjC9cLU8iZFLIKapHhoFHfIqOAW3ApducwHI/IKy9S27MKO99LjD5snCE7rwlUkKGBNlKFOGEn\nRWcyhTM9TMG1sDXyW2TIKahFRgE3W8D0HCrffoV83xYswwU8LPyt04UOjW3d6TZZmthMr7MJc+1z\npCPT6JsVpbepmUjUp7HR1wQoIkNEXd8iQeb7GO3t2B3tRDo2E29fC7kMRTeERxgfFx8ff+9HGrAQ\nEKFAlF4SXeuJpTuJNlYR2tJOuLsdxzFIJpXSIkNFLWqRADOSSQzHwYjZRAp9mI5D2PIJO91E6CNK\nDhsfD8o64nt72/6TsPHxsjkKPVl4cz3r7WlUOhl8PKprLVzXVTe4yBBQi1okqFwXM58H08To6KCy\npZqsGaEYiuJhbg3n0ixkQ52PPqVBZRWkmJl7AyubIt9bINsHbsGjo0MrX4oMFQW1SFC5bv/az2ah\nwJZsNV6skuoaiIQdQjjYFIe023ub0oeBAlX0Udu7kcSGV8l0F8kWQ3iGRbGorm+RoaKgFgkqq7T2\ns++WWqyd76Z4N3Qgve0OMSdDjALm1vvTQ92YtSitXR2mSKiQor5tLfH3XyeVzOJ4JqHQcHxcEBmf\nFNQiQWVZeJEIHZ0GhZxPIVXEitqkQglS4Tg91JCigiLD84scAmxyhMhQ6ErT85cu3ludobChnfp6\nPaYlMlQU1CIB5tQ1kiOGYRmQy1EVzpCK1dNdrMSkQJQCUfZ/WcuBMAADhxTVbIlMwg65RI0chbWt\nZPocBbXIENGob5EAcz0Dp64RtypO/h2HXK+D15cC16GOLiroG5IZyXaldJ86h4VLRb4TJ1NJMW+Q\nSuaZQhHXtRXWIkNALWqRALMs8H0DNxQh3XgAqbYM0WwP3aEJOFjD0pLexgRi+NTQQYIuqv0UU1lP\nJApt7Rr1LTJUFNQiAeX70NFh0NVl0NZqsO5dl6quddSG00TMPMPT4b0jG6ihjygZsk6ICdl3MSwD\nPzTYVbBFZHfK3vV977338vbbb2MYBmeeeSazZs0q9ylExoVk0sBxDJqbPbasSmJm0rQWJ5D284SN\nVjwMPErPOA93ZJsUiZCDvi6s+hoqqk2Gfuy5yPhU1hb16tWraW1t5frrr+frX/8699xzTzkPLzJu\nuC7k8yamCYbnEirmIGRTmTAJVYbJxWpJU4W/deKT4WRQ+nAQLaRxCVN/UI1GfYsMobIG9WuvvcZH\nPvIRACZPnkw6nSaTyZTzFCLjgutSmuwE8HIFnEwBIxLCChl0ZCvJ+RX0UEtxyOck25m59U8s5BA1\n8kTiUcLh0seFQkH3qkXKraxd393d3cycObP/+3g8Tnd3NxUVFbvdp6WlpZxFCBzVb3Qbqfq5LuD7\nRHracTyPkOlRiGdJx5vptWNMiKbJ9tSQxyY2AnFtAb4Xos2ZRmxDlMoJE2lrg7q60uu2DU1NjOiK\nWvq3ObqN9frtjyF9PGtbi2BPNm3aNJRFGFEtLS2q3yg20vXrW5ukJ+fiGyadxSKelyJrQKEqSps5\nhebOV8kWK6gmPazl8gCTPF4hQ0+6QM//9zQbCv+P6R+y6Oryqa/38TzYvNmnqWlkZiwb6Ws31FS/\n0W1/P4SUtes7kUjQ3d3d/31XVxeJRKKcpxAZH1yXxsocVtjA8wy8RB1eKEJVlUdDs0l1LI9VzGHh\nDHvRPMDAI0sVVXaR6kySynfXsGGDST5f6g0wzdI9dnWDiwxeWYP6sMMOY+XKlQCsXbuWRCJBLBYr\n5ylExodS3zcNDT7NzS4HzvZMJQYBAAAZRklEQVSYcXQCa/pEIuE8tV3rSdCDjYsBwzagbFv7OEU1\nnXYzvudR4aSoS20g3VUgnzfwthbG930FtUgZlLXre/bs2cycOZNFixZhGAZnn312OQ8vMn5sW5Cj\n9CWWBROaXCor3yMbep2U/x71JImS6x/cNRy23XJ2CNHgttKWjlOTyuP2RbDyOVy3EnNrYQzD0Ehw\nkTIo+z3qf/qnfyr3IUXGn60LchiOw7bkM99+m0R6I5mGRkITe3A32VgUhn0gmQlU002uUEeWGMVI\nhEi6h2img8rKSiwLPA8iEU9BLVIGmplMJKD8xkb8UAijWIRUCqOjA9/3CTfVYnkF/GiUAvawTnay\nrevbwyYccoibKTzHJE+EyqhLXY1DsWgQCvk0NmrpS5Fy0KIcIkFlGPhNTbiuCz09mLW1mJkcdPhY\nroMfrSCfq8ShgD1cRdr6d4QsNcVO3PQb9LpVVB3xN8w5vojTVMSKqctbpJzUohYJOBeLgl2FH4lB\n1CbsFbBsi87KyfRRS4HhG0zmAUVCuNgUCUG8hliDzdRpBuFYCDtmKaRFykwtapGA8v3SfN/5vInv\nR4kU6qgyssTrN5HBIOVGgSrSVBIiTXQYyuQCBUIUsfEjIXKxKvBsLKeAZ9uaR1RkCKhFLRJQ2xbl\nCId9bBv82bPJNrTQY9RhNdVh+C4FbLppJD9Mq1KXRpj7hMkTyheJ53qozrTh5gqlEWT7MMmRiOwf\ntahFAmjbohzh8HbBZ5o4s+fQ3vwhisVa2t77XyqLFfidLjV0EqM4LPeqTSwsfNJ2jHXx2bTEQnjN\nzVieB8kkflPTMJRCZPxQi1okgLZflGMn4RAd4ckUjBhFzyJOHxVbH9Ma6vasA6SxyVKNGQbXCOFY\nFmYsAuEw5rapyUSkbBTUIgFkWaUJQ3Ylu66DqhqD7vrpGJ5DlmqyRLdO7Tn0qklTQS9WJkd1z2Zy\ndpwtXgLf3/rhQkEtUlYKapEAsqzShCHeB4ZzF3MuZj7LhKkRWuqzhH2XohEjS5ThiEcHgzwV9FBP\nXzRBfupMqust6EvT0WGUPlxoQJlIWSmoRQKqsdEnFPIpFg0KBSgWDSwcamt8LFya5tQRqbWprixi\nYmAy9I9pWfh4lLrYXd/Cj8awEjXYfo5C2sWxIwpqkTLTYDKRgDIMaGrycV0X192af67JlnYDPI+G\nljDtiQoK2UqqUhYuBtYQ36U2ABMDD4+8WcGWbDWNRhRncxoz7lGsbRy2yVdExgu1qEUCzrJg2yPK\nlm0RqrDw80XCuRSTZobwKqtJVzSSoxKHoRtQtq0lnaKSVMVEjMY63JCF4xoUDZueLgh3J/WIlkiZ\nKahFRgvfx2hvpzFewE53wubNmI011DSYRKpMilTgMXTd39uOmyWGE6vG9vLYxRyFcCVOPIFnR6Dg\nYCSTQ1QCkfFJXd8io4SRTJZW07Jtamc34scOxF/zF3z/w1QafaTyDYR68oTo2dpFPTQS9JKyHFwi\nuEUPp7YOK1FHbZWPh0kony/NT6571SJloRa1yGjguqVnlLct9myahGIRzAOm4dfWYlRWEKmNkY3W\nkyU6JPN/u0AfNWSsakxMctUJKqcmaG4xqKvzMU0D09QjWiLlpqAWGQ1cd8cJUCwLz7YxE3Gs2hr8\nyhghiliGg0t066IZZTz91j8+Jm4ogmsZFBITaJlTjen50NmBbXt/ff5brWmRslFQi4wGlrXTBCh+\nfT1+KETVjFqcCRPJtswg3TCNXqOWXmpxyzj/d2nVrAhZo4ZIfRXVDVHC8ShePo/V2UE4n6E+4YDn\n4UX0iJZIOeketchoYFl4kUjpHvW27m/DwK+rw6+rI9zdR6YtTSQSxc/1EEpm8DC2toLBYnCzlvlA\njgq6K1rwqicRnjaZVMtsDgx30uAksdp9nM3gTZmC39g46OqKyF+pRS0ySviNjfihEEaxCIUCRrGI\nHwrhT5xI8YgjSTXNwLRNIjVRqK4gQw15LIpAnoHdsy61pNm6Slc9abOKLfYEkr0V0NNHId5Ab+Uk\nvJpa/ESi9PD3bqY+FZGBUYtaZLQwDPymptKI6m2jqrd2MXsTJpCfeSAZL0vErqRorCMfTpPvXE+I\nDDFS5HGJUcRk77/4pXWnwcUmR4QiEVLhRtLhWojGCNVUUFHtUel04bo19HgVVNk2pkZ8i5Sdglpk\ntNkuoPs3GR4VH5pE3svTUTEVsztMzsiSC1dDcgtVXg9ReoFOwhRxKXWHf/A/gP7pQYEslaX70lSR\npYreWBOb6g7Cjk9kmruZerODUL6CYk2CbKyWGCaWv92HCBEpCwW1yBhg4WJHDXpr6jEqm3Edl5zd\nTt63MXIRnHSYmBPBwsMiRyUZvK3LeGzrEncp3cf2gAIVuNj42BQjNaTsJtzJEwjVNxKxPQpVDcQa\nC7jxWtzKajwrgodFyPAU0iJlpqAWGQssi9pGiy3roqR68viJyRgNHtVdG8CxKbpVeIUQbiFMlAxh\nNhEhTQ4fc+swM58QLgYZKsmaNYTCJr2hWjqjEzCsEJnEVMJhn3ixi7gFhhfG9zw828aNN2DiaMS3\nyBBQUIuMBZaFGwqTKLZT2d2G5eTxazwy9dW02xPJGmmiPR00pt6FvE8f1YTwMPBwKLWmM1TRR4JQ\nKEQuUknINukxG/ENi56a6WCaGLaNG6sl01CBa3VRbGjGra7FNguYEVsjvkWGgIJaZIwI9XSC61GM\n14Pl4hdd3GyICXVp1juVRKw0nlVDqKOPnFtBFg+fImFcPAwsM0xlZQjP9WitnUxFpUWbMYNo2MGq\nr6PeTpGPVxCvcukpVJKrnkB21uFEYz6JZhM/pJa0yFBQUIuMBYUCoXSK8OQG8iYQ8vBNE8etJLT+\nHRpmRHH6QrA5RyZkQm83oUInZiFLwQLbKOBbYRzfJGx61FgZIuEIk531+HY1bjhKZY2FVxfFrI9S\nUT+BykOqSUzysGJa2FJkKCmoRcaCQgE8j7qET2+3QTZrYXgGVlOC3OYKLC9HT7iasBkhXg2RxgRs\nsXAzfeSMCAUvR8TPgx3D9wuEfY+iFcWL1WI4eWw/i+9Vko5PwPrQJOJ1IWKVObDVihYZagpqkbHA\ntkv3kA2YPsOnowPyeejtMQh/eCZeJETD5iRe/EOENxbozVnUTYzR1xMnlTWJuWlShk8xlqDRb8Ox\nqkgXo/jFHEZ1hFjIwTF8EhPD1DSaWKaDEYvga+CYyJArW1A/+eST/OpXv6K5uRmAQw89lM9//vPl\nOryI7Ilt41VXY+ZyGKEQDQ0+xaIPjo85oYrcjA+z5fV2sjUpvFgFXevSUG/htW8k3Wtgdm0kFAvh\nFTyKkQS2ZVJ0XCKRPMVpU3AKeQqVdZiZHHbXZmrmtGjgmMgwKWuL+phjjuH0008v5yFFZB95s2fD\nm29i9vWB52G4JlTEcT90EKZv4Da1kLF8+uqn0V14m2izT8gIUTTSuNUm2bpGSGfJxV3cbJGol6K+\n3qQ4sQbXtilMnEK+biI1k1xoqNdUoSLDRF3fImOFaeJ9+MN4hQIUCriWTWFLjLDpYwHRqIdlm7hE\nSU4+lKYPVdPbeABeWxKnsIE6Ow2dHVTFukn7cSJeDK+5Er+mASNiYyWqMSoiuH4OS7OPiQybsgb1\nmjVruP7663Fdl6985SvMmDGjnIcXkX1h22DbWEAk4uE4BqYJ9fU+4LFpk0kkYpIq2BSj1dTPrcDs\nsUn2FLBafOriHRTfayNS2IzrG/h2CKe+Ebe2DsMxsEJab1pkOBn+DqvR75ulS5eybNmyHbYdd9xx\nTJgwgblz5/LWW2/x05/+lJtvvrlsBRWR/ef70N4O2Wzpa8Mo5Xg8Dm+8AZlM6X2m4VNTaGfWpCw9\n3T5tmz28zk6oqiEcj1HfZOH7ELY8GlvCsHUsiogMvQEF9b4455xz+OlPf4pp7nklzU2bNg3F6QOh\npaVF9RvFxlL9drHgFi0tLbz33iYKBbAtF9tyaW83cN1Si7mjy4L2DshlCVnQ1OxR12JDU2Pg70+P\npWu3K6rf6NbS0rJf7y9b1/fDDz9MfX09H/vYx1i3bh3xeHyvIS0iw2MXC24BYId9It1JzHwex/Hx\nkxahyghuXSMNjeDW1eMVXXzHo2Yqem5aZASULag/9rGPcfvtt/PEE0/geR4LFy4s16FFpNxcFwoF\njNZWDN/HD4fxAM8yCTkOdCZx65u2BrxFoWDh4qGYFhl+ZQvq+vp6rr766nIdTkSGgu9jJEstaPJ5\nrPXr8Ssr8evqME0DwzDANDHzOdztRnYbhqHxYyIjRH3TIuOIkUxiOA5+OAyhEFgWhuNgdHZiWWDb\nHp4Hhg+Gt3W9aq80elxBLTIyFNQi44XrllrS28aOWFZ/C5p8HlyX+nofy/IpugZ5x6JYNAiFfBob\nh2TMqYjsA014IjJeuC47PORhWXi2jeG6pfvUnodhWTTUubiNYYoJA8ty1ZIWGWFqUYuMF9ta0Nvx\n6+tLC2u4LjgORrGIHwphNDVi25rXRCQI1KIWGS8sCy8SwXCcv3Z/GwZ+XR1eYyN+IrH757hEZMSo\nRS0yjviNjaUWc7EI+Xx/C9pvakJNaJFgUotaZDwxDPymptKjV83NuFtHfotIcKlFLTIelZ7FUkiL\njAIKahERkQBTUIuIiASYglpERCTAFNQiIiIBpqAWEREJMAW1iIhIgCmoRUREAkxBLSIiEmAKahER\nkQBTUIuIiASYglpERCTAFNQiIiIBpqAWEREJMAW1iIhIgCmoRUREAkxBLSIiEmAKahERkQBTUIuI\niASYglpERCTAFNQiIiIBNuCgXr16NV/96ld58cUX+7e99957LFq0iKuuuoq77rqrLAUUEREZzwYU\n1K2trTz66KPMnj17h+0///nPOfPMM7nuuuvIZDK8/PLLZSmkiIjIeDWgoE4kElxyySVUVFT0b3Mc\nh/b2dmbNmgXAkUceyWuvvVaeUoqIiIxToYHsFIlEdtrW29tLZWVl//c1NTV0dXUNvGQiIiKy96Be\nunQpy5Yt22HbggULOPzww/e4n+/7+1SAlpaWfXrfaKX6jW6q3+g1lusGqt94stegnj9/PvPnz9/r\ngeLxOH19ff3fd3Z2kkgk9rrfpk2b9vqe0aqlpUX1G8VUv9FrLNcNVL/Rbn8/hJTt8axQKMSkSZN4\n4403AHj++ef32uoWERGRPRvQPeqXXnqJRx55hI0bN7J27Voee+wxFi1axJlnnsnPfvYzfN9n1qxZ\nHHrooeUur4iIyLgyoKCeO3cuc+fO3Wn75MmT+e53vzvoQomIiEiJZiYTEREJMAW1iIhIgCmoRURE\nAkxBLSIiEmAKahERkQBTUIuIiASYglpERCTAFNQiIiIBpqAWEREJMAW1iIhIgCmoRUREAkxBLSIi\nEmAKahERkQBTUIuIiASYglpERCTAFNQiIiIBpqAWEREJMAW1iIhIgCmoRUREAkxBLSIiEmAKahER\nkQBTUIuIiASYglpERCTAFNQiIiIBpqAWEREJMAW1iIhIgCmoRUREAkxBLSIiEmAKahERkQALDXTH\n1atXc8stt7Bw4UKOPPJIAK655hry+TyRSASA008/nZkzZ5anpCIiIuPQgIK6tbWVRx99lNmzZ+/0\n2sKFC5k6deqgCyYiIiID7PpOJBJccsklVFRUlLs8IiIisp0Btai3dW3vyoMPPkhfXx+TJk3izDPP\nxLbtPR6rpaVlIEUYNVS/0U31G73Gct1A9RtP9hrUS5cuZdmyZTtsW7BgAYcffvhO7/30pz/N1KlT\nmTBhAnfddRePP/44p5566h6Pv2nTpv0s8ujR0tKi+o1iqt/oNZbrBqrfaLe/H0L2GtTz589n/vz5\n+3Swo446qv/rI488khUrVuxXYURERGRHZXs8y/d9rrvuOtLpNFAaFa5BZSIiIoMzoHvUL730Eo88\n8ggbN25k7dq1PPbYYyxatIj58+fz3e9+l2g0SiKRYMGCBeUur4iIyLgyoKCeO3cuc+fO3Wn7scce\ny7HHHjvoQomIiEiJZiYTEREJMAW1iIhIgCmoRUREAkxBLSIiEmAKahERkQBTUIuIiASYglpERCTA\nFNQiIiIBpqAWEREJMAW1iIhIgCmoRUREAkxBLSIiEmAKahERkQBTUIuIiASYglpERCTAFNQiIiIB\npqAWEREJMAW1iIhIgCmoRUREAkxBLSIiEmAKahERkQBTUIuIiASYglpERCTAFNQiIiIBpqAWEREJ\nMAW1iIhIgCmoRUREAkxBLSIiEmChgezkui533nknbW1teJ7HV77yFQ466CDee+897r77bgzDYOrU\nqZxzzjnlLq+IiMi4MqAW9dNPP000GuW6667j61//Oj//+c8B+PnPf86ZZ57JddddRyaT4eWXXy5r\nYUVERMabAQX13/7t33L66acDEI/HSaVSOI5De3s7s2bNAuDII4/ktddeK19JRURExqEBdX2HQn/d\n7Q9/+APHHXccvb29VFZW9m+vqamhq6tr8CUUEREZx/Ya1EuXLmXZsmU7bFuwYAGHH344jz/+OO++\n+y6XXXYZvb29O7zH9/19KkBLS8t+FHf0Uf1GN9Vv9BrLdQPVbzzZa1DPnz+f+fPn77R92bJlvPji\ni1x66aWEQiHi8Th9fX39r3d2dpJIJMpbWhERkXFmQPeo29raeOKJJ7jkkkuwbRsodYdPmjSJN954\nA4Dnn3+eww8/vHwlFRERGYcMf1/7qLfzwAMPsGLFChoaGvq3LVq0iNbWVn72s5/h+z6zZs3ijDPO\nKGthRURExpsBBbWIiIgMD81MJiIiEmAKahERkQAb0HPUgzEeph9dvXo1t9xyCwsXLuTII48E4Jpr\nriGfzxOJRAA4/fTTmTlz5kgWc8B2Vb+xdP0AnnzySX71q1/R3NwMwKGHHsrnP//5ES5Vedx77728\n/fbbGIbBmWee2T9J0ViwatUqbrnlFqZMmQLA1KlTOeuss0a4VIO3bt06brzxRk455RROPvlktmzZ\nwu23347nedTW1nL++ecTDodHupgD9sH63XHHHaxdu5bq6moATj31VObOnTvCpRy4+++/nzVr1uB5\nHp/97Gc54IAD9uv6DXtQbz/96Pr16/nxj3/M97///f7pR2fNmsXixYt5+eWXOeKII4a7eIPW2trK\no48+yuzZs3d6beHChUydOnUESlU+u6vfWLl+2zvmmGP6Z+AbK1avXk1rayvXX389GzZs4M477+T6\n668f6WKV1Zw5c/jWt7410sUom1wuxz333MPBBx/cv+3BBx/kk5/8JMcccwwPPPAAy5cv56STThrB\nUg7cruoH8OUvf7m/ITCavf7666xfv57rr7+evr4+vv3tb3PIIYfs1/Ub9q7vsT79aCKR4JJLLqGi\nomKkizIkdlW/sXT9xrrXXnuNj3zkIwBMnjyZdDpNJpMZ4VLJnoTDYa644ood5qVYtWoV8+bNA2De\nvHm8+uqrI1W8QdtV/caSOXPmcPHFFwNQWVlJPp/f7+s37C3qsT796Lau7V158MEH6evrY9KkSZx5\n5pn9z6CPJruq31i6fttbs2YN119/Pa7r8pWvfIUZM2aMdJEGrbu7e4dbLvF4nO7u7jH1wXLDhg38\n8Ic/JJVKsWDBAg499NCRLtKgWJaFZVk7bMvn8/1dpduu4Wi1q/oBPP744zz66KPU1NRw1llnEY/H\nR6B0g2eaJtFoFChNFHbEEUfw5z//eb+u35AG9VBPPzrS9lS/D/r0pz/N1KlTmTBhAnfddRePP/44\np5566nAVdUD2p37bGy3Xb5td1fO4445jwYIFzJ07l7feeovbb7+dm2++eYRKOHRG27Xam4kTJ7Jg\nwQKOOeYY2trauPbaa1myZMkODQQJvo9//ONUV1czffp0fv/73/PQQw9x9tlnj3SxBuWFF15g2bJl\nLFq0iAsuuGC/9h3Sf71jffrR3dVvV4466qj+r4888khWrFgxVMUqm32t32i9ftvsrZ4HHnggvb29\neJ6HaY7uByUSicQOn967urpG1bXam7q6Oo499lgAJkyYQG1tLZ2dnTQ1NY1wycorGo1SKBSwbXvU\n/b7ti0MOOaT/63nz5nHXXXeNYGkG75VXXuG3v/0t3/nOd6ioqNjv6zfs/+uMx+lHfd/nuuuuI51O\nA6UBPaN9UNn2xuL1e/jhh3n22WeB0ojUeDw+6kMa4LDDDmPlypUArF27lkQiQSwWG+FSlc8zzzzD\nI488ApS6+Xt6eqirqxvhUpXfIYcc0n8dV65cOep/3z7opptuoq2tDSjdj982in80ymQy3H///Vx+\n+eVUVVUB+3/9hn1msrE+/ehLL73EI488wsaNG4nH4yQSCRYtWsSKFSt4+OGHiUajJBIJFi5cuMf7\n2UG1u/pt2LBhTFy/bTo6Ovofn/A8jzPOOGPMPMb0n//5n6xZswbDMDj77LOZPn36SBepbLLZLIsX\nLyaTyeA4Dl/84hdH9WM9UPpAdd9995FMJrEsi7q6Oi644ALuuOMOisUiDQ0NnHvuuaO2e39X9Tv5\n5JN5+OGHsW2baDTKueeeS01NzUgXdUD++Mc/8tBDDzFx4sT+bd/4xjf4yU9+ss/XT1OIioiIBNjo\n78sTEREZwxTUIiIiAaagFhERCTAFtYiISIApqEVERAJMQS0iIhJgCmoREZEAU1CLiIgE2P8PCTol\nbfjZR3MAAAAASUVORK5CYII=\n",
            "text/plain": [
              "<matplotlib.figure.Figure at 0x7fc422f9dd68>"
            ]
          },
          "metadata": {
            "tags": []
          }
        }
      ]
    },
    {
      "metadata": {
        "id": "YdtvGLvmg341",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "## Variational Inference\n",
        "\n",
        "MAP can be used to find the mode (or one of the modes) of the posterior distribution, but does not provide any other insights about it. We next use variational inference, where the posterior distribtion $p(\\mathbf{W}, \\mathbf{Z} \\mid \\mathbf{X})$ is approximated using a variational distribution $q(\\mathbf{W}, \\mathbf{Z})$ parametrised by $\\boldsymbol{\\lambda}$. The aim is to find the variational parameters $\\boldsymbol{\\lambda}$ that _minimize_ the KL divergence between q and the posterior, $\\mathrm{KL}(q(\\mathbf{W}, \\mathbf{Z}) \\mid\\mid p(\\mathbf{W}, \\mathbf{Z} \\mid \\mathbf{X}))$, or equivalently, that _maximize_ the evidence lower bound, $\\mathbb{E}_{q(\\mathbf{W},\\mathbf{Z};\\boldsymbol{\\lambda})}\\left[ \\log p(\\mathbf{W},\\mathbf{Z},\\mathbf{X}) - \\log q(\\mathbf{W},\\mathbf{Z}; \\boldsymbol{\\lambda}) \\right]$.\n"
      ]
    },
    {
      "metadata": {
        "id": "BfXGOtI9g7bG",
        "colab_type": "code",
        "colab": {}
      },
      "cell_type": "code",
      "source": [
        "tf.reset_default_graph()\n",
        "\n",
        "def variational_model(qw_mean, qw_stddv, qz_mean, qz_stddv):\n",
        "  qw = ed.Normal(loc=qw_mean, scale=qw_stddv, name=\"qw\")\n",
        "  qz = ed.Normal(loc=qz_mean, scale=qz_stddv, name=\"qz\")\n",
        "  return qw, qz\n",
        "\n",
        "\n",
        "log_q = ed.make_log_joint_fn(variational_model)\n",
        "\n",
        "def target_q(qw, qz):\n",
        "  return log_q(qw_mean=qw_mean, qw_stddv=qw_stddv,\n",
        "               qz_mean=qz_mean, qz_stddv=qz_stddv,\n",
        "               qw=qw, qz=qz)\n",
        "\n",
        "\n",
        "qw_mean = tf.Variable(np.ones([data_dim, latent_dim]), dtype=tf.float32)\n",
        "qz_mean = tf.Variable(np.ones([latent_dim, num_datapoints]), dtype=tf.float32)\n",
        "qw_stddv = tf.nn.softplus(tf.Variable(-4 * np.ones([data_dim, latent_dim]), dtype=tf.float32))\n",
        "qz_stddv = tf.nn.softplus(tf.Variable(-4 * np.ones([latent_dim, num_datapoints]), dtype=tf.float32))\n",
        "\n",
        "qw, qz = variational_model(qw_mean=qw_mean, qw_stddv=qw_stddv,\n",
        "                           qz_mean=qz_mean, qz_stddv=qz_stddv)\n",
        "\n",
        "energy = target(qw, qz)\n",
        "entropy = -target_q(qw, qz)\n",
        "\n",
        "elbo = energy + entropy\n",
        "\n",
        "\n",
        "optimizer = tf.train.AdamOptimizer(learning_rate = 0.05)\n",
        "train = optimizer.minimize(-elbo)\n",
        "\n",
        "init = tf.global_variables_initializer()\n",
        "\n",
        "t = []\n",
        "\n",
        "num_epochs = 100\n",
        "\n",
        "with tf.Session() as sess:\n",
        "  sess.run(init)\n",
        "\n",
        "  for i in range(num_epochs):\n",
        "    sess.run(train)\n",
        "    if i % 5 == 0:\n",
        "      t.append(sess.run([elbo]))\n",
        "\n",
        "  w_mean_inferred = sess.run(qw_mean)\n",
        "  w_stddv_inferred = sess.run(qw_stddv)\n",
        "  z_mean_inferred = sess.run(qz_mean)\n",
        "  z_stddv_inferred = sess.run(qz_stddv)"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "metadata": {
        "id": "CpSsruVIqAv5",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 783
        },
        "outputId": "a47e3bb0-25ac-4544-87b4-13ae65d6096c"
      },
      "cell_type": "code",
      "source": [
        "print(\"Inferred axes:\")\n",
        "print(w_mean_inferred)\n",
        "print(\"Standard Deviation:\")\n",
        "print(w_stddv_inferred)\n",
        "\n",
        "plt.plot(range(1, num_epochs, 5), t)\n",
        "plt.show()\n",
        "\n",
        "with ed.interception(replace_latents(w_mean_inferred, z_mean_inferred)):\n",
        "  generate = probabilistic_pca(\n",
        "      data_dim=data_dim, latent_dim=latent_dim,\n",
        "      num_datapoints=num_datapoints, stddv_datapoints=stddv_datapoints)\n",
        "\n",
        "with tf.Session() as sess:\n",
        "  x_generated, _ = sess.run(generate)\n",
        "\n",
        "plt.scatter(x_train[0, :], x_train[1, :], color='blue', alpha=0.1, label='Actual data')\n",
        "plt.scatter(x_generated[0, :], x_generated[1, :], color='red', alpha=0.1, label='Simulated data (VI)')\n",
        "plt.legend()\n",
        "plt.axis([-20, 20, -20, 20])\n",
        "plt.show()"
      ],
      "execution_count": 11,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Inferred axes:\n",
            "[[-0.09666733]\n",
            " [ 1.4165041 ]]\n",
            "Standard Deviation:\n",
            "[[0.01517736]\n",
            " [0.01317804]]\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfsAAAFKCAYAAAAaHfYCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3X94FOW99/H3PdlASJYNiYCYIKIG\nVMCQgAWp2lONz3MsVtujUWuP0lThVKlCFW1pbcGacqz1OacNAirSIqiHqq1VzqHVg1AVEcESNAmI\niIoYEAIkK4HlR5K5nz8WFoJRSLLJ7I/P67q4Nrl3ZvzO18AnM7P3jLHWWkRERCRhOV4XICIiIh1L\nYS8iIpLgFPYiIiIJTmEvIiKS4BT2IiIiCU5hLyIikuAU9iIiIgnO53UB0bB169Y2r5uTk9Ou9aVl\n6mv0qacdQ32NPvU0+o7taU5OTqvW15G9iIhIglPYi4iIJDiFvYiISIJT2IuIiCQ4hb2IiEiCU9iL\niIgkOIW9iIhIglPYi4iIJDiFvYiISIJT2IuIiCS4hLhdroiItJ21FhobobEBrAXHgEk59OpEXo0x\nHfffb2qEJvfQaxO4TUe+jnx/6LWxEVz30J8msIe+tvaocRfrukfeO/rrZmPN1wHAAMaEvzCEe2AO\nvfFl482+PrKsyToJM6igQ3p3ohT2IiLHEQmjw4F47GtTIzQc89rYgI283wTYQxuLbJRmA/bor49Z\nFtvy8o0N0HDw0GvDke8bGrGNB48a+4LlGhuPjJ8IY8Bxjnp1Wvx+q89Hk2uPvGdM8zA/HNyH/1i3\n1f9P4okFnP98AtM907MaFPYikpRsw0HYvoXQxircdRXwaTV22xbYH2o50OOZMZCaCr5USO0Sfs3w\nH/W978i4Mc2PeI89aj7u9xaT4sDBg+Gxpqbwa0oKdO0Sfk1JASflyNcpvvAvBim+yJhxjnrv2HWO\nfj38i8bhP8Y55vujfxFJaTZmWloHwvUe/uXL2vCre+j1qPHwW+5Rv5wdtYx1I8uazCxPgx4U9iKS\n4GxoL3z6CXZbdfj10/ArO2vAuuw6euEuXcHfPfzaLSMckCm+I6++VPD5MIdeaek1stzRYw5HTgEf\n5fD3kXFz1JdHxo68HLP8sQF++PvDY6k+8IUDtqNOwbfklCR56l3ndbT9FPYiEvestfBZbfjo/NNP\njrxu2xIeP5Y/AHlnY045lcyzBrM7PQCnnApZJ4WP9kQSjMJeROKKtRa2bsZWrYatnxwJ9X17P7/w\nSb1hyDBMn1PhlL6YU06FPn0x3QORRbrn5FCfBEehktwU9iISF+zWzdh/vI79x/LwafjDUnzQ+xQ4\nZyimT9+jQj0X0zXNu4JFYojCXkRilt1WfSTgt3wcHkztAsO+ihk2CnPamdCzD8anf8pEvoz+hohI\nTLHbtx4K+NehelN40JcKBedjzrsAM/QrmLR0T2sUiTcKexHxnK3Ziv3H8nDAf/JReNDng6EjMOdd\niBk6AtNNAS/SVgp7EfGE3bHtSMBv/iA8mOKDc88LB3zBCEy639siRRKEwl5EOo3dVXMk4De9Hx5M\nSQl/Yv68izAFIzEZCniRaFPYi0iHsrvrsG++Gg74jzaEBx0HBhWGr8EXno/xB758IyLSLgp7EekQ\n1lrsG0uwT/8+PAfeOOHpceddiCkc1Wyuu4h0LIW9iESdrd2J+8RMqFoNad0w19yEOf/rmEAPr0sT\nSUoKexGJms8dzQ8qwBlzO+akXl6XJpLUFPYiEhW2blf4aL7yH9C1G+bG8ZiL/rlTH8AiIi1T2ItI\nu3zuaP6coTjfux1zUm+vSxORQxT2ItJmOpoXiQ8KexFptfDR/FLs03N0NC8SBxT2ItIqOpoXiT8K\nexE5ITqaF4lfCnsROa7PHc3fMB7zNR3Ni8QLhb2IfCFrLXbFUuwfdTQvEs8U9iLSIh3NiyQOhb2I\nNKOjeZHEo7AXkQgb3IU7X0fzIolGYS8iANh1b+M++hsI7QkfzY+5DdPzZK/LEpEoUNiLCO4rf8Mu\neBQcB/Ovt2D+6Rs6mhdJIAp7kSRmm5qwz/4Bu+S/wR/A+eHPMHmDvC5LRKJMYS+SpOy+EO7sB8PP\nnD/lVJzbf4Hp1cfrskSkAyjsRZKQ3bENd8avYOtmGDIMZ9zdmPQMr8sSkQ6isBdJMnbju7iz/h3q\nP8MUXYG55iZMSorXZYlIB1LYiyQR981XsPOmg+tivnsLzsWjvS5JRDqBwl4kCVjXxS78L+yiZ6Bb\nBs4tP8YMKvS6LBHpJAp7kQRnDxzAnftbWP0G9OoT/iDeKad6XZaIdCKFvUgCs8Ha8AfxPt4IAwbh\n3PozTPeA12WJSCdT2IskKLv5A9wZ06BuJ+arReFb36amel2WiHigTWHf1NTEww8/zPbt23Fdlxtv\nvJGzzz6bTZs2MWfOHIwx9OvXj3HjxgGwcOFCVqxYgTGG4uJihg0bRigUoqysjFAoRFpaGhMnTsTv\n91NRUcGCBQtwHIfCwkKKi4ujusMiycCueRN3zn9Aw0HM1d/D/PNVuiOeSBJz2rLSa6+9RlpaGqWl\npdxyyy3MmzcPgHnz5lFSUkJpaSmhUIg1a9ZQU1PD8uXLKS0tZfLkycyfPx/XdVm0aBGDBw+mtLSU\nkSNH8vzzzwMwd+5cJk2aRGlpKRUVFVRXV0dvb0USnLUW98U/4z58PwDOrT/FuexqBb1IkmtT2F90\n0UWMGTMGgEAgwJ49e2hsbKSmpoa8vDwAhg8fTmVlJVVVVRQWFuLz+QgEAvTq1Yvq6mqqqqoYMWJE\ns2W3b9+O3++nZ8+ekSP7ysrKKO2qSGKzjQ3YedOxf54Hmdk4P/k1pvB8r8sSkRjQptP4Pt+R1RYt\nWsQFF1zA7t27ycg4cgeuzMxM6urq8Pv9BAJHPhAUCASoq6sjGAxGxjMzMwkGg83GDo9v27atLSWK\nJBVbvxv3kfthw1o4LQ/ntnswPU7yuiwRiRHHDfslS5awdOnSZmPXXHMNBQUFvPjii3z00Uf85Cc/\nYffu3c2Wsda2uL2WxluzbEtycnJOaLmOWl9apr5GX0s9bfhkEzsfnIz7aTXdLigi+85f4qSleVBd\n/NLPavSpp9HXnp4eN+yLioooKir63PjSpUtZvXo1d999d+QUfX19feT92tpasrKyyM7OZuvWrZHx\nuro6srKyyMrKIhgMkp6eHln28NjR28jOzj7uThy9/dbKyclp1/rSMvU1+lrqqV33Nu4jD8C+vZjL\nr+XAld9lW22tRxXGJ/2sRp96Gn3H9rS1wd+ma/bbt29n8eLF3HXXXXTp0gUIn9rPzc1l/fr1AKxa\ntYqCggKGDBlCeXk5jY2N1NbWUltbS9++fcnPz2fFihUArFy5koKCAnr37s2+ffuoqamhqamJ8vJy\n8vPz21KiSMJzX/kbbtm90HAAc/MdON++AeO06a+0iCS4Nl2zX7JkCfX19dx///2RsZ///OeUlJQw\ne/ZsrLXk5eVFgrqoqIipU6cCMHbsWBzHYfTo0UyfPp0pU6aQnp7OhAkTIu+XlZUBMGrUKJ0KEjmG\ntRb753nYl56D7pk443+qZ9CLyJcy9kQvjMcwncaPPepr9OXk5LCl+hPsE7Owry+GPrk4E6bqGfTt\npJ/V6FNPo6+9p/F1Bz2ROGEbDuLOfjB8j/vT8nAmTsV0z/S6LBGJAwp7kThg9+9jxy/vgDUrYeAQ\nnNt+jumW7nVZIhInFPYiMc7urcedfh8HPnwPho7A+be7MV26el2WiMQRhb1IDLPBWtzfTYUtH5N+\n8TfYf+04jE9/bUWkdTRPRyRG2R3bcH8zGbZ8jLnkm2Tf+UsFvYi0if7lEIlBtnoT7u/uhc9qMVd8\nB3PF9ZpDLyJtprAXiTH2g/W40++D0B7MdWNxLr3S65JEJM4p7EViiF23Bnfmv0NjA+b7P8L56iVe\nlyQiCUBhLxIj7OrluI/9BxiDc+tkTIEeTysi0aGwF4kB7rL/xT4xC7p0DT+e9mw9E0JEokdhL+Ix\n96W/YP80F/zdcSbei+k/wOuSRCTBKOxFPGKtxf7lCezf/gQ9TsK58z7MKad6XZaIJCCFvYgHrNuE\nfepR7GsvQu+ccNCf1NvrskQkQSnsRTqZbWzA/uF32LeWwamn4/zoXkwgy+uyRCSBKexFOpE9cAD3\nkfuhqhzyBuHc/nNMut/rskQkwSnsRTqJDe3BfagUNr4LQ4bj3DIZ01UPtBGRjqewF+kE9rO68O1v\nqz/CfOUizE0/wvhSvS5LRJKEwl6kg9md23F/OwVqPsX802WY7/4A46R4XZaIJBGFvUgHsrvrcP/f\nPbCrBjP6Gsy3b8AY43VZIpJkFPYiHcQeOIA7Y1o46K/4Ds6V3/W6JBFJUnpmpkgHsK6L+4ffwkcb\nMKMuwVxxvdcliUgSU9iLdAD73HwofwMGDsGM+aFO3YuIpxT2IlHmvvYS9qXn4ORcnPE/1afuRcRz\nCnuRKLLr1mCfejj8UJsJUzAZ3b0uSUREYS8SLXbLZtxHHgDHwRl/D6b3KV6XJCICKOxFosLursN9\n6D7YF8KUTMQMGOR1SSIiEQp7kXZqNsXuW9/FGflPXpckItKMwl6kHZpPsbsYc/l1XpckIvI5CnuR\ndrB/eeLIFLsbb9MUOxGJSQp7kTZyX3sJ++Kfj0yxS9UUOxGJTQp7kTaw694+aordLzTFTkRimsJe\npJXCU+x+fdQUuxyvSxIR+VIKe5FW0BQ7EYlHCnuRE2QPHjXF7kpNsROR+KGwFzkBzabYnX8x5pua\nYici8UNhL3IC7F+egNVvwMDBmDGaYici8UVhL3Ic7rL/DU+x652Dc6um2IlI/FHYi3yJyBS7jENP\nsfMHvC5JRKTVFPYiX8BuPfQUO2Nwxv8Mc7Km2IlIfFLYi7TA7q7DnX4f7NuL+d4EzMDBXpckItJm\nCnuRYzSbYnfF9Tjnf93rkkRE2kVhL3KU5lPsvo654jtelyQi0m4Ke5Gj2OePnmJ3u6bYiUhCUNiL\nHOIufgH7N02xE5HEo7AXAdwl/4195vfQIxtnoqbYiUhiUdhL0nP/vgj7x8cgMwtn0q/0FDsRSTgK\ne0lq7qsvYv/rUQj0wJk0DdOnr9cliYhEna8tK3322WfMnDmThoYGGhsbGTNmDAMGDGDTpk3MmTMH\nYwz9+vVj3LhxACxcuJAVK1ZgjKG4uJhhw4YRCoUoKysjFAqRlpbGxIkT8fv9VFRUsGDBAhzHobCw\nkOLi4qjusMhh7rL/xT45C7pn4tz5K8wpCnoRSUxtOrJftmwZX/va15g6dSrXX389Tz/9NADz5s2j\npKSE0tJSQqEQa9asoaamhuXLl1NaWsrkyZOZP38+ruuyaNEiBg8eTGlpKSNHjuT5558HYO7cuUya\nNInS0lIqKiqorq6O3t6KHOIuX4J9Yib4u+PcWYrJ7ed1SSIiHaZNYf/Nb36TCy+8EICdO3eSnZ1N\nY2MjNTU15OXlATB8+HAqKyupqqqisLAQn89HIBCgV69eVFdXU1VVxYgRI5otu337dvx+Pz179owc\n2VdWVkZpV0XC3Df/jp03HdL94SP6vv29LklEpEO16TQ+QDAY5IEHHmDfvn1MmTKF3bt3k5GREXk/\nMzOTuro6/H4/gcCRTzYHAgHq6uoIBoOR8czMTILBYLOxw+Pbtm1ra4kin+OufBX7hzLolo5zx32Y\nU0/3uiQRkQ533LBfsmQJS5cubTZ2zTXXUFBQwP333095eTmzZs1i/PjxzZax1ra4vZbGW7NsS3Jy\n2vfp6fauLy2Ltb6Gli1m1x9+i0lPp/e0WXQZMMjrklot1nqaKNTX6FNPo689PT1u2BcVFVFUVNRs\nbN26dezZswe/38+wYcOYOXMmgUCA+vr6yDK1tbVkZWWRnZ3N1q1bI+N1dXVkZWWRlZVFMBgkPT09\nsuzhsaO3kZ2dfdydOHr7rZWTk9Ou9aVlsdZXu/oN3Nm/gS5dMROmsjOjB8RQfSci1nqaKNTX6FNP\no+/YnrY2+Nt0zX7lypW8+uqrAGzevJmTTjoJn89Hbm4u69evB2DVqlUUFBQwZMgQysvLaWxspLa2\nltraWvr27Ut+fj4rVqyIbK+goIDevXuzb98+ampqaGpqory8nPz8/LaUKBJh334T97EHIbUrzsR7\nMWec5XVJIiKdqk3X7K+++mpmzpzJqlWraGhoYOzYsQCUlJQwe/ZsrLXk5eVFgrqoqIipU6cCMHbs\nWBzHYfTo0UyfPp0pU6aQnp7OhAkTIu+XlZUBMGrUKJ0Kknax76zCfeQ34EvFmTgVk3eO1yWJiHQ6\nY0/0wngM02n82BMLfbWV/8Cd9e/gODgT7sWcNcTTetorFnqaiNTX6FNPo8+T0/gisc5WlePOuh+M\ng3PbL+I+6EVE2kNhLwnHrns7fEQPOLfdgzlnqMcViYh4S2EvCcWur8Cd+SuwLs4Pf4YZVOh1SSIi\nnlPYS8KwG6pwHyqFJjf8PPohw70uSUQkJijsJSHYjetwp98HTY04t/wEk/8Vr0sSEYkZCnuJe/aD\n9bi/+yU0HMT5tx9jCkZ6XZKISExR2Etcsx9twC27FxoO4Iy7CzNslNcliYjEHIW9xC378Ubc306F\n/fsxN9+JOe9Cr0sSEYlJbX7qnYiX7Ifv4Zb9EvaHMDfdgTPia16XJCISsxT2Eles24R96S/YF54C\n18WUTMA5/+telyUiEtMU9hI3bO0O3N//FjZUQWY2zk0TNY9eROQEKOwlLrhvvY59ciaE9kLB+Thj\nbsN0D3hdlohIXFDYS0yz+0PY/5qNXbE0/Cz6G3+Iuej/YozxujQRkbihsJeYZT9Yj/v7/4Qd2+C0\nPJyxkzB9cr0uS0Qk7ijsJebYpibsX5/F/s8fwVrMN4oxV16P8aV6XZqISFxS2EtMsTu2hY/mP1gP\n2T1xbrpTj6cVEWknhb3EBGstduUr2Kcegf37MOddiLlhPCbD73VpIiJxT2EvnrOhPdinHsGueg26\ndsN8/0eYURfrQ3giIlGisBdP2Q1rw6fta3fAmWfj3Hwnplcfr8sSEUkoCnvxhG1sxP73Auzf/gyA\nueI7mMuvw6SkeFyZiEjiUdhLp7Pbt+LO+Q/Y9D6c1Ds8pS7vHK/LEhFJWAp76TTWWuzri7FPz4ED\n+8PX5a//AaZbuteliYgkNIW9dAq7ZzfuEzOhfAV0y8CMu0tPqhMR6SQKe+lw9t13cP/wOwjugoGD\nw3PnT+rldVkiIklDYS8dxu4PYf88D/vK3yAlBfMvN2Iuuwrj6EN4IiKdSWEvHWL/mpW4v/0l7KqB\nnH4435+I6T/A67JERJKSwl6iyu4LYf80lx2vvQSOgxl9Leab12FSdV97ERGvKOwlauzaNbjzH4La\nnaT2z6PphvGY0/K8LktEJOkp7KXdbGgv9tk/YF9fHL42f8V3OPnmiXy6Y4fXpYmICAp7aSdbuTo8\npa5uJ/Q9PXxtvt8ZOm0vIhJDFPbSJnbvHuwzv8e+sQRSfJgrvxt+7rxPP1IiIrFG/zJLq9l33sJ9\nciYEa6HfmTjfn4Dpe7rXZYmIyBdQ2MsJs3vrsX+cg33z7+Gj+W/fgPnnq3Q0LyIS4/SvtJwQ+/ZK\n3CdnwWd1cFpe+Np87mlelyUiIidAYS9fyu7ZjV3wGHbVq+DzYa4ag/m//6JH0YqIxBGFvXwhW/4G\n7pMPQ/1ncPpAnJIJmJx+XpclIiKtpLCXz7H1n2EXzMa+tQx8qZjiEsyl39LRvIhInFLYSzN29XLc\npx4JH82fcRZOyUTMKX29LktERNpBYS8AWNfFPjsX+/ILkNoFc81NmEuv0BPqREQSgMJesA0HcX//\nn7D6DTjlVJzxP8P0yfW6LBERiRKFfZKze+txZ0yDjetg4GCc8fdgMvxelyUiIlGksE9idsc23Om/\nhG1bMOddiLnpR5jULl6XJSIiUaawT1L244240++D3cHwvPmrv4dxHK/LEhGRDqCwT0K28h+4j/4G\nDh7AfOffcIq+6XVJIiLSgRT2ScZ97SXsUw9Dig/nlsmYYaO8LklERDqYwj5JWGuxLzyFXfQM+Lvj\n3PYLzJlne12WiIh0AoV9ErCNDdj5M7Ar/g69+uBMmKqpdSIiSURhn+DsvhDuw/fDu+9A/wE4t/8C\nE+jhdVkiItKJ2hX2wWCQO+64g7vuuovBgwezadMm5syZgzGGfv36MW7cOAAWLlzIihUrMMZQXFzM\nsGHDCIVClJWVEQqFSEtLY+LEifj9fioqKliwYAGO41BYWEhxcXFUdjQZ2bpd4al11Ztg6AiccXdh\nuqZ5XZaIiHSyds21evLJJ+ndu3fk+3nz5lFSUkJpaSmhUIg1a9ZQU1PD8uXLKS0tZfLkycyfPx/X\ndVm0aBGDBw+mtLSUkSNH8vzzzwMwd+5cJk2aRGlpKRUVFVRXV7dvD5OU3fIx7v13Q/UmzNe/gXPr\nTxX0IiJJqs1hX1VVRVpaGv36hR952tjYSE1NDXl5eQAMHz6cyspKqqqqKCwsxOfzEQgE6NWrF9XV\n1VRVVTFixIhmy27fvh2/30/Pnj0jR/aVlZVR2M3kYtdX4D4wGep2hp8//91b9MQ6EZEk1qbT+I2N\njTz77LP8+Mc/5vHHHwdg9+7dZGRkRJbJzMykrq4Ov99PIBCIjAcCAerq6ggGg5HxzMxMgsFgs7HD\n49u2bTtuPTk5OW3ZjaitH0v2vvIitWX3ApB9VykZF3/Ds1oSqa+xQj3tGOpr9Kmn0deenh437Jcs\nWcLSpUubjRUUFFBUVNQs3I9lrT3h8dYs25KtW7ee0HItycnJadf6scJai33xz9jn5kO3dJxbf8pn\nZw3lM4/2LVH6GkvU046hvkafehp9x/a0tcF/3LAvKiqiqKio2dgvfvELXNflpZdeYtu2bWzcuJGJ\nEydSX18fWaa2tpasrCyys7ObFVhXV0dWVhZZWVkEg0HS09Mjyx4eO3ob2dnZrdqhZGTdJuyC2dhX\n/gZZPXEmTsXknuZ1WSIiEiPadM2+tLSUadOmMW3aNIYNG8bYsWPp378/ubm5rF+/HoBVq1ZRUFDA\nkCFDKC8vp7GxkdraWmpra+nbty/5+fmsWLECgJUrV1JQUEDv3r3Zt28fNTU1NDU1UV5eTn5+fvT2\nNgHZAwdwZ90fDvq+/XF++qCCXkREmonqPPuSkhJmz56NtZa8vLxIUBcVFTF16lQAxo4di+M4jB49\nmunTpzNlyhTS09OZMGFC5P2ysjIARo0apes+X8LuDuLO+BV8tAHOGRq+/W36F19aERGR5GTsiV4Y\nj2HJeM3e1mzF/d29sGMbZtTFmDG3YXypXpcVEa99jWXqacdQX6NPPY2+Dr9mL7HHuk24jzwQDvrL\nr8V8618xxnhdloiIxCiFfRyyr70En3yEGXUxzrdv8LocERGJce26g550PrtnN/b5p6BbOubqEq/L\nERGROKCwjzP2L0/C3nrMFddjMrO8LkdEROKAwj6O2I8/wC57CU45FXPx5V6XIyIicUJhHyes6+Iu\neBSsxbn+3zA+fdxCREROjMI+Ttg3X4EP1mOGX4A5Z6jX5YiISBxR2McBG9qL/fPj0KUL5pqbvC5H\nRETijMI+Dtj/+SPsDmK+cQ3mpF5elyMiInFGYR/j7NbN2KX/A736YP75X7wuR0RE4pDCPoZZa3H/\n+Bg0NeFcNxaT2sXrkkREJA4p7GNZ+Qp49x049zzI/4rX1YiISJxS2Mcoe+AA7jO/B58vfFSve9+L\niEgbKexjlH3xT1C7A/N/vo05WY/5FRGRtlPYxyC7Yxv2xeegx0mY0dd4XY6IiMQ5hX0Mcp+eA40N\nmGtvwqR187ocERGJcwr7GGMrV8M7q+CsczHnXeh1OSIikgAU9jHENjSEp9o5Ds53xulDeSIiEhUK\n+xhiX34BarZiLr4c07e/1+WIiEiCUNjHCFu3C7voGeieibnyeq/LERGRBKKwjxH2T3PhwH7MVWMw\n6X6vyxERkQSisI8B9r0q7KrX4PSBmK8WeV2OiIgkGIW9x2xTE+6CR8EYnOt/gHH0v0RERKJLyeIx\n+8rfYMvHmAsuxZw+wOtyREQkASnsPWR3B7ELn4L0DMxVY7wuR0REEpTC3kP2L09AaC/myn/FdM/0\nuhwREUlQCnuP2I82YJe/DLmnYb7+Da/LERGRBKaw94B1Xdz/ehSsDX8oLyXF65JERCSBKew9YN9Y\nApvex3zlIsxZQ7wuR0REEpzCvpPZ0B7sc/Ohaxqm+PtelyMiIklAYd/J7MIFUP8Z5vJrMdk9vS5H\nRESSgMK+E9nqTdi/L4LeOZhLv+V1OSIikiQU9p3EWou7YDa4bvjxtampXpckIiJJQmHfSew/XocN\nVTB0BObc4V6XIyIiSURh3wns/n3YZ/4AvlSc68Z6XY6IiCQZhX0nsH99FoK7MJddhenVx+tyREQk\nySjsO5jdF8IufgGye2IuK/a6HBERSUIK+462vgIaGzBfvRTTtavX1YiISBJS2HcwW1UOgBkyzONK\nREQkWSnsO5C1Fru2HNIzoL+eVS8iIt5Q2Hek7VtgVw2cM1QPuxEREc8o7DuQXbsGADNYp/BFRMQ7\nCvsOpLAXEZFYoLDvILbhILxXATn99MAbERHxlMK+o7y/Dg4exAwu9LoSERFJcgr7DqJT+CIiEisU\n9h3Eri2HLl1g4GCvSxERkSTna8tKr7zyCk8//TQnn3wyAPn5+Vx11VVs2rSJOXPmYIyhX79+jBs3\nDoCFCxeyYsUKjDEUFxczbNgwQqEQZWVlhEIh0tLSmDhxIn6/n4qKChYsWIDjOBQWFlJcHH+3mLV1\nu2DLxzBkGCa1i9fliIhIkmtT2AOMGjWKMWPGNBubN28eJSUl5OXlUVZWxpo1a8jNzWX58uVMmzaN\nUCjElClTKCgoYNGiRQwePJgrr7ySl19+meeff54bbriBuXPncs8995Cdnc29997L+eefT9++fdu9\no53Jrj101zxdrxcRkRgQtdP4jY2N1NTUkJeXB8Dw4cOprKykqqqKwsJCfD4fgUCAXr16UV1dTVVV\nFSNGjGi27Pbt2/H7/fTs2TNyZF9ZWRmtEjtP5Hq9nlsvIiLea/OR/bvvvsu0adNoamrixhtvJDMz\nk4yMjMj7mZmZ1NXV4ff7CQSRFtMLAAAOQUlEQVQCkfFAIEBdXR3BYDAynpmZSTAYbDZ2eHzbtm1t\nLdET1m3CrnsbsntBn1yvyxERETl+2C9ZsoSlS5c2G7vgggu45pprGDZsGBs2bGDGjBncc889zZax\n1ra4vZbGW7NsS3Jyck5ouY5a/2gH1ldSE9pDxkX/h+zc5A77aPZVwtTTjqG+Rp96Gn3t6elxw76o\nqIiioqIvfH/gwIHs3r2b7t27U19fHxmvra0lKyuL7Oxstm7dGhmvq6sjKyuLrKwsgsEg6enpkWUP\njx29jezs7OPuxNHbb62cnJx2rX8s95X/BWDfGWdFdbvxJtp9FfW0o6iv0aeeRt+xPW1t8Lfpmv0L\nL7zA66+/DsDmzZsJBAKkpqaSm5vL+vXrAVi1ahUFBQUMGTKE8vJyGhsbqa2tpba2lr59+5Kfn8+K\nFSsAWLlyJQUFBfTu3Zt9+/ZRU1NDU1MT5eXl5Ofnt6VEz9h1a8Bx4OyhXpciIiICtPGa/YUXXsiM\nGTNYvHgxruty6623AlBSUsLs2bOx1pKXlxcJ6qKiIqZOnQrA2LFjcRyH0aNHM336dKZMmUJ6ejoT\nJkyIvF9WVgaEP/EfT6eC7N498OEGOPNsTHrG8VcQERHpBMae6IXxGBYrp/HtP17HffQ3mG/9K843\nr4vKNuOVTuNFn3raMdTX6FNPo8+T0/jSMt0iV0REYpHCPkqstdiqcvAH4LQzvS5HREQkQmEfLVs/\ngeAuzKACjKO2iohI7FAqRYlduzr8hW6RKyIiMUZhHyWR6/WDFPYiIhJbFPZRYA8cgA1roe/pmB7H\nvwmQiIhIZ1LYR8OGKmhs0FPuREQkJinsoyDySNshmnInIiKxR2EfBXbtGuiaBnnneF2KiIjI5yjs\n28nuqoFt1XB2PsaX6nU5IiIin6Owb6fIKXxdrxcRkRilsG+nI7fIVdiLiEhsUti3g21shHffgV59\nML3j5+l8IiKSXBT27fHRBtgX0oNvREQkpins28FW6Xq9iIjEPoV9O9i15ZDig7PP9boUERGRL6Sw\nbyNb/xls/gDyzsGkpXtdjoiIyBdS2LeRXfc2WKtT+CIiEvMU9m0VmV+vD+eJiEhsU9i3gXXd8Pz6\nQA/o29/rckRERL6Uwr4tqjfB7iBmcCHGUQtFRCS2Kana4PBd89ApfBERiQMK+zawa8vBGMygAq9L\nEREROS6FfSvZ/ftg47vQ70xM90yvyxERETkuhX1rvVcJTY36FL6IiMQNhX0r6Ra5IiISbxT2rWTX\nlkO3dDjjLK9LEREROSEK+1awNVthxzY4Ox/j83ldjoiIyAlR2LfC4Sl3ul4vIiLxRGHfCkfCXtfr\nRUQkfijsT5BtbID1FdAnF9PzZK/LEREROWEK+xO18V04sF+n8EVEJO4o7E+QrteLiEi8UtifILu2\nHHypMHCI16WIiIi0isL+BNjP6uCTj2DgYEzXrl6XIyIi0ioK+xOgT+GLiEg8U9ifiLWHb5Gr6/Ui\nIhJ/FPbHYV0Xu+5t6HES5PTzuhwREZFWU9gfz+YPYM9uzOBCjDFeVyMiItJqCvvjiDzlbohO4YuI\nSHxS2B+HXbsGjAPnFHhdioiISJso7L+EDe2FD9fD6QMwGX6vyxEREWkThf2XWf8OuK4+hS8iInFN\nYf8lNL9eREQSgcL+C1hrw2Gf7ofTB3hdjoiISJsp7L/Iti2wqwYzqADjpHhdjYiISJsp7L+APXTX\nPHQKX0RE4pyvrSsuXLiQZcuW4fP5uPnmm8nLy2PTpk3MmTMHYwz9+vVj3LhxkWVXrFiBMYbi4mKG\nDRtGKBSirKyMUChEWloaEydOxO/3U1FRwYIFC3Ach8LCQoqLi6O2s61xOOzNIIW9iIjEtzYd2X/y\nySe88cYb/PrXv2bcuHGUl4eDcd68eZSUlFBaWkooFGLNmjXU1NSwfPlySktLmTx5MvPnz8d1XRYt\nWsTgwYMpLS1l5MiRPP/88wDMnTuXSZMmUVpaSkVFBdXV1dHb2xNkGw7ChirIPQ2T3bPT//siIiLR\n1KawX716NaNGjSIlJYUzzjiDa6+9lsbGRmpqasjLywNg+PDhVFZWUlVVRWFhIT6fj0AgQK9evaiu\nrqaqqooRI0Y0W3b79u34/X569uwZObKvrKyM3t6eqPfXwsGD+hS+iIgkhDadxt+xYweO4zBt2jSa\nmpoYM2YMgUCAjIyMyDKZmZnU1dXh9/sJBAKR8UAgQF1dHcFgMDKemZlJMBhsNnZ4fNu2bW3dtzaL\n3CJX8+tFRCQBHDfslyxZwtKlS5uNBYNBCgoK+NnPfsZ7773Ho48+yt13391sGWtti9trabw1y7Yk\nJyfnhJY70fU/3VBJU9eu5HytCNOla7u2ncza+/9FPk897Rjqa/Spp9HXnp4eN+yLioooKipqNvbM\nM8+Qk5ODMYazzz6bmpoaAoEA9fX1kWVqa2vJysoiOzubrVu3Rsbr6urIysoiKyuLYDBIenp6ZNnD\nY0dvIzs7+7g7cfT2WysnJ6fZ+rZ2J+7HH8KQ4Xy6c1ebt5vsju2rtJ962jHU1+hTT6Pv2J62Nvjb\ndM2+oKCAd955B4AtW7bQs2dPfD4fubm5rF+/HoBVq1ZRUFDAkCFDKC8vp7GxkdraWmpra+nbty/5\n+fmsWLECgJUrV1JQUEDv3r3Zt28fNTU1NDU1UV5eTn5+fltKbLPIp/D1lDsREUkQbbpmP3DgQN5+\n+23uueceAG6++WYASkpKmD17NtZa8vLyIkFdVFTE1KlTARg7diyO4zB69GimT5/OlClTSE9PZ8KE\nCZH3y8rKABg1alTnnwrSLXJFRCTBGHuiF8ZjWLRO49umJtw7b4Ru6Tj3P4YxJlolJh2dxos+9bRj\nqK/Rp55Gnyen8RPWpvchtAczeJiCXkREEobC/mjWhZQUzIiveV2JiIhI1LT5drmJyOQNwpn5J0yK\nHnwjIiKJQ0f2x1DQi4hIolHYi4iIJDiFvYiISIJT2IuIiCQ4hb2IiEiCU9iLiIgkOIW9iIhIglPY\ni4iIJDiFvYiISIJT2IuIiCQ4hb2IiEiCU9iLiIgkuIR4nr2IiIh8MR3Zi4iIJDiFvYiISIJT2IuI\niCQ4hb2IiEiCU9iLiIgkOIW9iIhIgvN5XYCXHn/8cd5//32MMZSUlJCXl+d1SXHrySef5N1338V1\nXb797W9z5plnMmPGDFzXpUePHtx+++2kpqZ6XWbcOXjwIJMmTeLqq69myJAh6mk7LVu2jIULF+I4\nDtdddx39+vVTT9tp//79zJgxg71799LQ0EBxcTE9evRgzpw5GGPo168f48aN87rMuLF582YefPBB\nLr/8ci677DJ27tzZ4s/osmXL+Otf/4oxhksvvZRLLrnkS7ebcu+9997bObsQW9atW8fq1auZOnUq\nAwcO5LHHHqOoqMjrsuJSVVUVb731FlOnTmXkyJE8+OCD7Ny5kwsvvJAxY8bw0UcfUVNTw5lnnul1\nqXHnmWeeoa6ujry8PF555RX1tB3q6+uZMWMG06ZN4/zzz2fx4sVUVVWpp+308ssv4/P5uP322yko\nKKCsrIz33nuPMWPGUFxczOuvv07Xrl055ZRTvC415h3+xemMM86gR48e5OXl8fjjj3/uZzQ3N5fp\n06dz3333cckll/Dwww9zwQUX0KVLly/cdtKexq+srOQrX/kKAH379mXv3r2EQiGPq4pPgwYN4o47\n7gAgIyODAwcOsHbtWs477zwAzjvvPCoqKrwsMS5t2bKF6upqCgsLAdTTdqqsrOTcc8+lW7duZGVl\n8YMf/EA9jYLu3btTX18PwN69e/H7/dTU1ETOlA4fPpzKykovS4wbqamp/PSnPyUrKysy1tLP6MaN\nGznzzDNJT0+nS5cunHXWWaxfv/5Lt520YR8MBgkEApHvA4EAwWDQw4ril+M4pKWlAbB06VIKCws5\ncOBA5HSoets28+fP53vf+17ke/W0fWpqajhw4AAPPPAAU6ZMobKyUj2NggsuuICdO3dy++23M3Xq\nVG688UYyMjIi72dmZlJXV+dhhfEjJSXlc0fnLf2MtiW/kjbsj6W7BrffW2+9xdKlS7n55pu9LiXu\nvfrqqwwcOJDevXt7XUpCqa+v56677mL8+PHMmjVLf++j4LXXXqNnz5489NBDTJkyhYceeqjZ++px\nbEjaD+hlZWU1+02orq6u2akTaZ23336b5557jnvuuYf09HTS0tI4ePAgXbp0oba2Vr1tpfLycmpq\naigvL2fXrl2kpqaqp+2UmZnJWWedRUpKCn369KFbt26kpKSop+303nvvMXToUAD69+/PwYMHaWpq\niryvvrZPS3/vj82v2tpaBgwY8KXbSdoj+6FDh/Lmm28C8OGHH5KVlUW3bt08rio+hUIhnnzySSZP\nnozf7wfg3HPPjfT3zTffpKCgwMsS484dd9zB/fffz7Rp07jkkku4+uqr1dN2Gjp0KFVVVbiuS319\nPfv371dPo6BPnz5s3LgRgB07dtCtWzdyc3Mj15BXrVqlvrZDSz+jAwYM4IMPPmDv3r3s37+f9957\nj3POOedLt5PUT7176qmnePfddzHGcPPNN9O/f3+vS4pLL7/8Ms8++2yzT9v+8Ic/5JFHHqGhoYGe\nPXsyfvx4fL6kPZHULs888wy9e/dm6NChzJgxQz1th8WLF7N06VIArr766sgUUfW07fbv38+sWbP4\n7LPPcF2X6667jh49ejB79mysteTl5TX77Il8sQ8//JD58+ezY8cOUlJSyM7OZsKECcycOfNzP6Nv\nvvkmCxcuxBjDZZddxkUXXfSl207qsBcREUkGSXsaX0REJFko7EVERBKcwl5ERCTBKexFREQSnMJe\nREQkwSnsRUREEpzCXkREJMEp7EVERBLc/wfl7L/3En7U/wAAAABJRU5ErkJggg==\n",
            "text/plain": [
              "<matplotlib.figure.Figure at 0x7fc422f814e0>"
            ]
          },
          "metadata": {
            "tags": []
          }
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeoAAAFOCAYAAABJ4F7CAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3XmUXGWd//H3vbfWXqq7urs6nc6+\nSJAhEJKILArjCYMoZxB14nCcI0SQOQQFRFEB4xHkhwvIHkRFRRjGBWZUMih4JGFTjCKLQBIWyZ5O\nd1d6ra79Lr8/ijRZSdKp7r7V/Xl5cpLcqrr1fVK0n3qe+9znMTzP8xARERFfMke7ABEREdk/BbWI\niIiPKahFRER8TEEtIiLiYwpqERERH1NQi4iI+FjgcF58//33s27dOlzX5eyzz2bWrFksX74c13Wp\nr6/nkksuIRgMlqtWERGRcccY6n3Ur7zyCv/3f//HVVddRSqV4stf/jJz587luOOO48QTT+RnP/sZ\nTU1NnH766eWuWUREZNwY8tD3UUcdxeWXXw5AdXU1+XyeNWvWsHDhQgAWLlzISy+9VJ4qRURExqkh\nB7VpmkQiEQBWrVrFcccdRz6fHxzqjsVi9Pb2lqdKERGRceqwJ5M9++yzrFq1igsuuKAc9YiIiMgu\nDmsy2YsvvsivfvUrvvrVr1JVVUUkEqFQKBAKheju7iYejx/wHG1tbYdTgq+1traqfRVM7atcY7lt\noPZVutbW1kN6/pB71JlMhvvvv58rr7ySmpoaAObOncvq1asBWL16NfPmzRvq6UVERITD6FE/88wz\npFIpbrnllsFjn/3sZ/n+97/PY489RlNTE6eeempZihQRERmvhhzUp512Gqeddtpex7/2ta8dVkEi\nIiLyNq1MJiIi4mMKahERER9TUIuIiPiYglpERMTHFNQiIuPYypUrOe200+jr63vH57355pts2bLl\nkM9/zjnnkM1m9/v4XXfdxaOPPrrfx//+97/T09NzyO87liioRUQqhONAoVD6vVxWrlxJa2srTz75\n5Ds+7+mnn2br1q3le+OD9Mgjj4z75agPa2UyEREZfp4HyaRBPm+yc8PD3l6DRMLDMIZ+3v7+fl59\n9VW+9KUv8Ytf/IKzzjoLgDfeeINbb70VwzA4+uijOf3001mxYgX19fXU19dz7bXXcs899xCNRrnr\nrruYMWMG73//+/l//+//kcvlyOVyXHrppbz73e/e5/v+4Q9/4Oc//zmJRIJwOMyMGTNIp9ODr3dd\nl4suuoh0Os0f//hHNm7cyLXXXsuTTz7Jk08+ieu6nHDCCZx33nlDb3wFUY9aRMTnkkkD2zYIBj1C\nIQiFwLYNksnDSGngiSee4IQTTuD4449n27ZtJJNJAO644w6+8IUvsHz5crq7u6mqquL444/nwgsv\n3G/4dnd3c+aZZ3LLLbdw4YUX8vOf/3yfz/M8jx/96EfcdNNNXH/99Wzbtm2v13/hC1/g5z//OQsX\nLmT27Nl85StfYcKECQDcfvvtfO973+PRRx8lnU4fVvsrhXrUIiI+5jiQz5sEg95ux02zdNxxHCxr\naOdeuXIl5557LpZlccopp/D444/ziU98gi1btjBr1iwArr766oM6Vzwe57/+67/45S9/SbFYHNxd\ncU/9/f1Eo9HBvSCOPvrovV4PYO2jUeFwmM9//vNYlkVfXx+pVIrq6upDbnelUVCLiPiY4zA43L0n\nz/NwHIYU1MlkknXr1vG9730PwzDI5/PU1NTwiU98AuMA4+m7Pm7bNgD/+7//S1NTE1dffTWvvfYa\nd911135rNs23B3Nd193r9V1dXVx33XW7va69vZ0HH3yQu+++m2g0yqc//elDb3SF0tC3iIiPWRb7\nDU7DMA6rN3322Wfz4x//mB/96Efcd9999Pf3s23bNqZPn87atWsBuOGGG9i0aROGYeC8NYutqqqK\nrq4uHMcZfF5fX9/grlBPP/30YIDvKRaLkU6nGRgYwLZtXnnllb1e/9hjjw2+3jRLowZ9fX3U19cT\njUZ5/fXX6ejooFgsDq3xFca65pprrhnNAlKp1Gi+/bCqra1V+yqY2le5xlLbSkPc4DjG4MSxmpoa\n+vsHCIU8amuHdt7bb7+dT3/604ND0IZhkE6n2bRpEx//+Me54447eOSRR5g8eTKnnXYaqVSKH/7w\nh8ycOZPJkydz++238/zzzzNhwgRaWlpYsGABP/jBD3jmmWd4//vfz1NPPUVVVRV///vfOfvsswkG\ng4PvU11dzc0338zf/vY3mpqamDRp0m6vP+OMM/j9739PVVUVtbW1/PjHP+Zf/uVfePHFF1mxYgWZ\nTIYjjjiC3/3ud3zwgx8sxz/ziKo9xA/N8PY3pjJCxvqeo2pf5VL7KtdYa9ues75bWlro7d1+2LO+\n/WqsfX57OtT9qHWNWkTE5wwDmps9HMfBcWDKFAiFRrWPJSNIQS0iUiEs6+1fMn5oMpmIiIiPKahF\nRER8TEEtIiLiYwpqERERH1NQi4iMU7/+9a+5+OKLueyyy1i6dCnPPfccAMuXL2f79u1DPu/nP/95\nNmzYsN/H//SnPx30YiUbNmzg85///Ds+5yMf+cg7Pn6gncH29M1vfpMXX3yRCy64YLfjnudxzjnn\n0NPTM7h951133cUf//jHQzr/oVJQi4hUijLuc9ne3s5vf/tbbr/9dm677Ta++tWvct999wHwuc99\njokTJx72e+zPgw8+uN+Vy8qtWCzy4IMPHvTz//znPxMOh5k3bx6BQIBNmzYNPvbyyy8zZcqUwUVi\nAC644ALuuececrlcWevelW7PEhHxO8/DSCYx8/nBdb+N3l68RIKhrngyMDBAoVCgWCwSCASYPHky\nt912G1DqEV922WU8+eST9PX1sW3bNrZv387555/PI488Qnt7O9/+9rfp7Ozk17/+Nddeey1Q6tk+\n9NBDg++RTCb55je/CZTWBL/yyitZs2YNa9eu5Stf+Qo33XQTDz/8MCtXrsQ0Td73vvfxiU98gvb2\ndj772c8SDAYHNwfZleM4XHfddSSTSebMmTN4/LnnnuMnP/kJgUCA2tpavv71r3PnnXeyfv16brnl\nFv7zP//zgFtx/s///A9Lly4FYNGiRTz++OMsWbIEKO02tmjRot2eHwqFOOmkk1i5ciVnnnnmkD6L\nA1GPWkTE54xkEsO28YJBdu5zadg2xlvbUg7F7NmzOfLII/nkJz/Jt7/9bR5//PHBtbx31d/fzw03\n3MCpp57K73//e2644Qb++Z//mWeeeeaA79HV1cW5557LLbfcwoc+9CEeeughTj/9dBoaGvjOd77D\njh07eOqpp7jjjju47bbbeOqpp+jo6OC+++7jAx/4ALfeeitNTU17nffZZ5/FcRzuvPNOTjvtNPr7\n+4HSktTLli3jtttuo6qqimeffZZzzjmHKVOmcPnllx9wK07bttmwYcPgl4MPfOADg8Pmruvyl7/8\nhVNOOWWveo455hheeOGFA/+jD5F61CIifuY4pZ70W2tlDzJNzHy+FK5DXAHl6quvZtOmTTz77LP8\n4he/YMWKFdx88827PWdnj7OxsXFwc5B4PD4Yju+koaGBO+64g5/+9KekUimOOOKI3R5ft24dW7du\n5fLLLwcgk8nQ3t7Om2++yb//+78DcOyxx/KXv/xlt9dt2rSJf/qnfwLgqKOOIhwOA1BfX8+NN96I\n4zhs376d+fPn7/a6A23F2d/fTywWG2xnIpGgvr6e9evX09/fz7ve9S6qqqr2amcikRjcy3s4KKhF\nRPzMcd5xm8uh7nPpeR7FYpFp06Yxbdo0PvrRj3LeeefR0dGx2/N23Rd61z/vq6Y9e+T33HMP73nP\nezjrrLN48skn+fOf/7zb48FgkBNOOIEvfvGLux3/9a9/PRiW+3qf/W2VecMNN/Ctb32LadOmDQ7j\n7+pgt+Lc1aJFi3jiiScYGBjYa9h7pGjoW0TEzyzrHbe5HGpv+ne/+x3f/e53B4MwnU7jed5uE6UO\npLq6mq6uLgDefPNNMpnMbo/v3LrS87zdZnrv3LryiCOO4MUXXySXy+F5HnfccQf5fJ4ZM2bw2muv\nAexzSHnKlCmDj7/yyiuD502n0zQ3NzMwMMALL7xAsVjcbXvOA23FGYvF6O/v3+3Lwamnnspzzz3H\nSy+9xHvf+959/jvs2LFjn0P05aIetYiIn1kWbjiMYdulPS93cl3ccHjIQX3GGWewefNmLr74YqLR\nKLZt87nPfW5wGPlgzJo1i2g0yuc+9zmOPvpoWlpadnv8X//1X7n99ttpaWnhox/9KDfffDPPPvss\nxx57LJdccgm33norH//4x7nssssGJ5OFw2HOPfdcLr74Yp5++ul9TiZ773vfyyOPPMJll13GrFmz\nBkPyIx/5CJdccgmTJ0/mnHPO4d577+W9730vtm1zzTXXcM455/Ctb32LJ598krPPPptVq1bxyCOP\n8KEPfQiAQCDA9OnTWb9+/eD71tbWEo/HicVihEKhff47vPTSSxx33HEH/e92qLTN5TAaD1u1qX2V\nayy3b8y1bY9Z3y0tLWw/zFnffjaan9+f/vQn/vrXvw5eNz+QQqHA0qVLWb58OdFo9KBec6jbXGro\nW0TE7wwDr7kZp7UVt6UFpkzBa24ekyE92k4++WSy2Sxr1649qOf/5Cc/4bzzzjvokB6Kwxr63rx5\nMzfeeCNnnnkmZ5xxxuD9arW1tQCcddZZe826ExGRIdI+lyPi6quvPujnXnTRRcNYScmQgzqXy3HP\nPfdw9NFH73b8k5/8JAsWLDjswkREROQwhr6DwSBXXXXVIc0QFBERkUNz2JPJHnjgAWKx2ODQd29v\nL7ZtU1dXx/nnn08sFitXrSIiIuNOWW/POuWUU6itrWX69On85je/4cEHH9xr95E9jamZmXsYczNP\n96D2Vbax3L6x3DZQ+yrdqM76njt3LtOnTwdg4cKFbN68uZynFxERGXfKGtTf/e53B5efW7NmDVOm\nTCnn6UVERMadIQ99r1+/nvvuu49kMollWaxevZozzjiDW2+9lVAoRCQS4eKLLy5nrSIiIuPOkIN6\n5syZXHPNNXsdP+GEEw6nHhEREdmFViYTERHxMQW1iIiIjymoRUREfExBLSIi4mMKahERER9TUIuI\niPiYglpERMTHFNQiIiI+pqAWERHxMQW1iIiIjymoRUREfExBLSIi4mMKahERER9TUIuIiPiYglpE\nRMTHFNQiIiI+pqAWERHxMQW1iIiIjymoRUREfExBLSIi4mMKahERER9TUIuIiPiYglpERMTHFNQi\nIiI+pqAWERHxMQW1iIiIjymoRUREfExBLSIi4mOBw3nx5s2bufHGGznzzDM544wz2LFjB8uXL8d1\nXerr67nkkksIBoPlqlVERGTcGXKPOpfLcc8993D00UcPHnvggQf44Ac/yDe+8Q1aWlp4/PHHy1Kk\niIjIeDXkoA4Gg1x11VXE4/HBY2vWrGHhwoUALFy4kJdeeunwKxQRERnHhjz0bVkWlmXtdiyfzw8O\ndcdiMXp7ew+vOhERkXHusK5Rl0Nra+tolzCs1L7KpvZVrrHcNlD7xpOyBnUkEqFQKBAKheju7t5t\nWHx/2traylmCr7S2tqp9FUztq1xjuW2g9lW6Q/0SUtbbs+bOncvq1asBWL16NfPmzSvn6UVERMad\nIfeo169fz3333UcymcSyLFavXs2ll17KnXfeyWOPPUZTUxOnnnpqOWsVEREZd4Yc1DNnzuSaa67Z\n6/jXvva1w6lHREREdqGVyURERHxMQS0iIuJjCmoREREfU1CLiIj4mIJaRETExxTUIiIiPqagFhER\n8TEFtYiIiI8pqEVERHxMQS0iIuJjCmoREREfU1CLiIj4mIJaRETExxTUIiIiPqagFhER8TEFtYiI\niI8pqEVERHxMQS0iIuJjCmoREREfU1CLiIj4mIJaZJxyHCgUSr+LiH8FRrsAERlZngfJpEGhAO3t\nJoZhEA67JBIehjHa1YnIntSjFhlnkkkD2zYIhSAUgmDQw7YNkkmltIgfKahFxhHHgXzexNzjJ980\nS8c1DC7iPwpqkXHEccDzvH0+5nmeglrEhxTUIuOIZYGxnwvRhmFgWSNckIgckIJaZByxLAiHXVx3\n9+OuWzquoBbxHwW1yHiwy71YiYRHIOBRKJQOFYsGgYBHIrHvIXERGV1lvT1rzZo13HzzzUyZMgWA\nqVOncv7555fzLUTkUHgeRjKJmc/jeR6GYeCGwzQnEkxoASj1otWTFvGvst9HfdRRR/HFL36x3KcV\nkSEwkkkM28YLBgHwAMO2IZnEmjSJUGh06xORA9PQt8hY5TiY+Tz7uhfLzOe1JJlIhSh7j3rr1q18\n5zvfYWBggMWLF3PMMceU+y1Exi3HKf06qOFqx3nHW7EU1CKVwfD295M8BN3d3bz66quceOKJdHR0\ncO2113LHHXcQCGilUpHD4XnQ2QnZbOnPhgHRKDQ3s99lP52Cg7NpC1YktHeoFwowZcpg2u+cWBYK\nQcg6lG8DIjLcypqgDQ0NnHTSSQC0tLRQX19Pd3c3zc3N+31NW1tbOUvwldbWVrWvgvmpfZ2dpWU/\ndx3Fdl3Yvt2juXn379o71/LO502sHf2YjkMoAo2Nb63l7bp4gQATp09n69Y2XnvNJJUycR2PSH+S\nWCjDzBkuGAZ2MIzRnMAKVNbyon767IaD2lfZWltbD+n5Zb1G/fTTT7NixQoAent76evro6GhoZxv\nITLuHOyynzvvwGpvL4V60LSxGuuwQuDminRvL2IUi3iBAF4iAcBrr5nkciahENRkkoRNm6wd4flX\norR3RUhuc+l4pZvOToPyjb2JyKEoa4964cKF3Hbbbfztb3/Dtm0+85nPaNhb5DAdaNlP24aurlIP\nOp/36Gg3aHA6aa7JYgCeAUUrSCraSCRuYIUsLBcG+hyS22xc18AtegS7srjBKLlcafGTiRNdQhET\nilnyeZdk0tyr9y4iw6+sKRqNRrnyyivLeUqRce/tZT/fDknHKQ19u65Bd7eB6xr090M6bZDZ2EUR\nl1RNhKlTXXp6DOyCy46+FG1dCcIhl8zGJMFChva/vUEk4OCZFvGYTTFSR2exgUDQYtYsl1AIDA8s\nHPL5AI7j6LK1yAhTd1fE53Yu+2nbBoYB3d0G+XwppC3LwzAsAgGPfN6gr9sl15UlH7DYshnefDNA\naysMDFh42RxGXY4NL/dBZxcTjQwDm3swYxaeAZm0S3BymOpcNymvkb4+g6oqD88Az7TwHG9wjpmI\njBwFtUgFSCQ8kkloa7MoFsGyDCIRl1jMo63NoK3NxHUg2NFJcWM7qUwQy/LIGAE2DDRS4/RBaoCu\n7UV4cz2OZ9EbDpJKWQQHeqkJ5XALaWyniFUI4k2ppVAwyWVcrGgELAvD1aYdIqNBQS1SAQyjNGs7\nm/UwTQ/TLPVsHQf6+kz6+w3oSJLfUSDc65DvGiBEjqBXILblFZzmZvLVjWxtt2jY2E3AyVMT3MG7\nswN4dh6iUaL5FEYgQzQcJmzm6IpMgWNaYEKCYNKguVmbdoiMBgW1SIVwHDCMUkjvuvtVXx+8sdZh\n2oa1OH1pJqReoTaXJGqnyHkRooE8mzdOod2eSq+dYEr6TSLkqTO3U00fDgHqsylCZpHIwDb6rYl0\n9ENbbyuNAIaF57n7K0tEhpmCWqRCmGZpdrfnmXieh+cZ9PR45HIGiY41xLLtWPkUDYUO4n3rCTt5\ncoQpEKLKNplBkndhY+IRo58WdwtB8oQoEs3bZEL1OLkiVjBNfaKAkf0HLa6H19SAFbIoFExNJhMZ\nBQpqkQrR3W0QCIDjePT1GWSzJm+8YbBtg80/DWwm4mWpym4nnmkjUkgT8jJU4ZAlTJR2DGCAOmws\nJrKNOH3YWEQoYmCSL0Qp9g5gFAxyVjfVzeCsDxM8ag6EqvA8TSYTGQ0KapEKsHPRk6Ymjzdfdymm\nXXp7ghhYTMhtJNHzD8L5FDM6/kTC3kI1/ZTy1MHFwsPBxsTAJEM1IfLUkMfBwMDDwsCiSCGTJmfE\nyaQboadIqKeGvvU7iB09FcPQZDKR0aCgFqkAjgOe62K0t9PUXwrYoB2gIZ9iW7afXMajvnc79XYH\nMXqJkAfAAFzA2uX3GJnBJQkdwMbCI0CIIrgpCsUothfAK3qQL0BHF5lpk6iOmQpqkVGgoBapAJbp\nEV7/GoFsmkC/SQCId2fo39JNzbYUicwmJqdepYoBrLcWRgmw7zWCd83aABDAwcZ567Ec2Vw18d6N\n7HCm8saGIFmjQHVrkQlTo7iuS0uLt9+NQESk/BTUIhUg0NVJtDhAMRiGoEVoYAdOspPAth00ZFKQ\nzWKSw6KIiUuQUm/6YFjAzg0vTTwipOnpG6DHMSl0peh1O9ncapFxTOrrobfX48gjXYW1yAhRUIv4\nneNgptPEG6AnDbg2VrKDTNolm3FxBgoMFKuxCREgTxD7oEN6pxBQpDQ87uIS8FLMSL1M6h8dBJv7\n6X++lf7YXMKBEGDR2ekxYYLW/RYZCQpqEb9zHDzDwDINGuIu9bkOMv/opDsbxu7L4A5ksIp5IvQS\npnjIIb1TgFIvvJoBIECOPJF8Aa/HwnzVpNj7MuFFxxGc1Uo6nsBp8nTNWmQEKKhF/M6yMAIB3FAI\nM5kk4DlgmbheAMcziZp5JvM6UQoY7H4N+lAYvB3W0Ac45IkQHkiRjxRJtecI/WM7sRqHUNzGmTJB\nQS0yAhTUIn5nWbjhMEZdHUZHB65n0N1nYfT3ELIcOkOTqE9vogVnyL3p3d4OCONhkSKIQ4AcTneA\nHdk+Ms+mCQe3U2O0E5w0D1onoIvVIsNLQS1SAbxEArZvx6muoX3A4o0c9HRsobvTo3ZgOw10EiY9\n5N70nix2zhjPYGITdfOQNejI1JJNRqjt6CbQ0YYbNPGam8v0riKyLwpqkUrguri1MdbviNG2PcS6\npEfU7qUx/xrVdBAhh4dRlh71TgalSWYmBQz6qXF76HZyFL0AhaIByR2YDQ04Wq5MZFgpqEX8zPMw\nOjsx29roX9dB9OWNTGrvp2NDPU5vminF1wjST5QsUbJl61HvZFDqWQdwqKGbatLYrkHBrMHO2liF\nAlpXVGR47Ws9BBHxCSOZxOzowOvuo5jKY9cl6M+Eqc7uYFrfy0zlH0xnC7X0ER2uGoAAHvX0UJfe\nhtvdi110sG1wPUMhLTLM1KMW8SvHwcxmMfN5nHwe0yni2lDMuzRmNjODddTSR4wUFkOf7X0gJqV7\nrD087Gg1RjZL4Y0tZBomkh6oI9AVIJHQamUiw0U9ahG/chw82y6t8e05BAyPSCFFONVDa3oTNaSJ\nUiBIKaSHMydLS40WcfqzGAWbGiuLa5p4TU3YtkEyqZQWGS7qUYv41Vv3TxumiWGZBHN9hPq6qMu2\n0ehto4oUIQoj8kNsAgUC5IohuotVpJlGdboVttqE64PU1hraq1pkmKhHLeJXloUbjZbuoc7lqHb6\n8Xr6iA10ECFPiOKwDXfvyQMCFIjSjxWyyIfryOUMMhkT2zbo6irNKROR8lNQi/iYl0jgNjXh5fP0\n5cJYToGY2YdHCBdjxILaBCJ4TGETie2vULdtDYVkH86OfjzXo1jU0LfIcNHQt4ifGQZufZzu6mls\nmXQEmW0NzPQ2UEUWkyDht/adHgkhoJoUXjoLPQMEujrJbE1iNnoEmrToichwUVCL+NyOHQambWBu\nbaMu20Eo10eEAULkh3UC2Z5K71XAxaFjR4CMW6A3UCTR0Ebj7EZdnxYZJhr6FvExx4GcGyG4o51Q\nTzs1DOAaBg4WLiYjfVk4DExiK82FrczIrCPsZRhIFsj25RXUIsNEPWoRH3Mc8GwbNxQikukh1LaV\nom3hYQImAdwRrccAYqRozm0Gs4bWjhco1M8mGNACZSLDRT1qER+zLDDyeayubqLtWwh0d1AkSJ4w\nBqXZ2CMtTJY43YQCBSYXN1Nt99G9oVezvkWGiYJaxKc8D7q6DLJr15Pd0kU27WB7QaqCNkCZNrU8\nNBYQIU+cXsx8noBXwAuFiXZt1/1ZIsOk7EPfP/3pT3njjTcwDIMlS5Ywe/bscr+FyLiQTBrYmSIt\n9nay1TUUCWHmgQKE8AAPh5H/tl0gjINFwCkQ69tEgX7q+j3MtkaYOVn7U4uUWVl/xteuXUt7ezvX\nX389F110Effcc085Ty8ybjgO5PMmgWKWQDFPvnU66WgjgRBEyBImhznC16ehdI3awKWIRcaoxvJM\nJlT1U2elsdwiRjI54jWJjHVlDeqXX36Z97znPQBMnjyZdDpNJpMp51uIjAuOA57n4ZkWjmdiY5Fu\nmEIm2kCh4OECNuERn/Vd2kmriIOFGY1SqKknRIFwPgXBIPZAHqegIXCRcirr0Hdvby8zZ84c/Hss\nFqO3t5eqqqr9vqa1tbWcJfiO2lfZRqt9Oy/3hiwH1r+J80YSoymKtbWGYrCOYj6PiT3ide3swxcI\nUjRDmAMeqWwDgeYY/dsbqGsK4BUmEIqGaG4e3VFw/bdZ2cZ6+w7FsN6e5XkHnpPa1tY2nCWMqtbW\nVrWvgo12+3p7DWzbwGxowk69Qah3A9gp7ECAXhqI59NUj3BN5lu/LBzcTJ4iJl0buuhzIFK/hnR7\nhHA0iNsfYPt2j+bm0ZiXPvqf3XBT+yrboX4JKevQdzwep7e3d/DvPT09xOPxcr6FyLiRSHgEAh52\nfYL8xMlgO3i1ddBYhxUOkic6KkPfHqWZ31V2H5loAguHNFV09EbJFwPQ3Y1plq6xayK4yOEra1Af\ne+yxrF69GoD169cTj8eJRqPlfAuRccMwoLnZY+IUgxnTHJpnRDEnJChE67FyaarJ4ozSmkUGRUJu\nlmihn0iuG6+mlv58mHRVAiOXK+2l7XkKapEyKOtP+Zw5c5g5cybLli3DMAwuuOCCcp5eZFyyLAhE\nTAITIgR6k2TjObI7qil2BTEIU8AmNIL1mEAImxp66c6nwTQYaJlFzqrHcW1MA1zXwTACWqlMpAzK\n/nX8P/7jP8p9SpHxzXFwEwm8+nqi4c3YhQwhK0cRiwJhCqRHNKhLM78LWIQJ2HnM/ADVyU30N8eo\nqgLTgiIW4bCroBYpA61MJuJqPL93AAAXY0lEQVR3loURDGLPnIlbU4PZFCNbHacQa6KXRoqYI7qU\naBAI4RAlT8LZTijVhZsaoNHupD7mkDejBMImicToTCQTGWu0KYeI31kWbjiM0dgIsRiYXeRqmiAf\noK9/gHpiuPQykp1XF3AIAhZRJ8vs4Eaa58UpTqzHmNiAFVBIi5SLglrE5xwHnLoEwbyN1VCPFwlB\nNEBf2CKLi4VNARipaZsO4GJgE8CsCWM21NDSYuA1NxCY0AgBLSEqUk4KahGf8rzSet/5vInneRhu\nKzXTF9C8pY2OnjS9bpAwSfJUYVAgSmFE6jIAhwC2GcGzgriuiVnI4Pb2gqmraSLlpp8qEZ9KJksL\nngSDHqEQBCMmhWgdfXWTMKdPoy+SwDBNcgQJ4I3IderSkHcpqItGCNsKQTSEE09AKITR3T0CVYiM\nLwpqER/auSnHXh3URIKcWY1lFIkXOmgJJAliYsGIbHpZWmPcJEMNthUmZdWTjbXghoJ4ra2Y+by2\nuxQpMwW1iA/t3JRjTx4GdsskesITKBCh4IWxCWJgDfteWi5QBLKEsY0ohVCEtFVD57Rj8Y6cA65b\nqllBLVJWCmoRH7IsMPaxo4VhQCobJFwbhGiYXjOOh/nWztTDy6QU1j00sMNrpC0wlTcbFjDQOIPe\nPhPPcUs16+ZpkbJSUIv4kGVBOOzi7tFNLjgWxZo6GhpcmupsnECAMFnAG5HdqS1Kq5KBR9ArMDnS\nQ1N9Adcz6Oo2cMNhBbVImSmoRXxq56YcxaJBoQDFooEZNIlNrMJrbCQ+MUx1rYFdVU+R6LAHtQN4\nmAQoDi5cWj+whbCZxQya5IhiNySGuQqR8Ue3Z4n41M5NORzHwXHe7qi2ORNxUp3UNJrkq0yKhRrs\nTBUF0gRwh+3bd+kadYQCUXJmFflQLdmsTTGaIN00h5p6C8d1R3ThFZHxQD1qEZ+zLAiFSr9bFoQj\nHvkZR2LPmE38XXEKNQ30VrXST5wMDMvVapfS9pYFLLppJlXbQjhRS29sKgOxCdgE6ekxNeotMgwU\n1CIVJpHwsMIWAzOOhuZ6oq3VhGIhXDNGgRgFSrOzyxXY3lvnyxBkLUexvuYY0nWt9NRMJheK4Vql\nNdH2NUtdRA6fhr5FKoxhQHOLgbG+HW9iPYH+fsJ2lL5CjP5+i4JtEGGAKM5hD0OXetHQR4Q8MfqC\nk6iJByAaIm+E6aiexqxoFZblUV/v7TZELyLloaAWqUBGZyeBUBBn7rtxuvO4A31EalNkcgFSdiN5\nwpj0YFDEZOg/6M5bv2yidNFEwQjj2dAXm0l+4jQCs49gynQIBkuT3hTSIuWnoBapNI6DmU7jmSZm\nopHC+9+PV19HIRcgne+imB8g5cRIEWcymwiTw2Ro17kMII/BP5jF85xEf3gqtfFaIlVhrEQzs09O\nAB6ui/afFhkmCmqRSuM4eIaBaRh4QLw5RN+0aWS3ZfB61pHJRMilPaqcPrbj0EQ7FgMEObRlRh0g\nD/SRIMkkcjVNNLZaTHp3iFDUIxtJUcy5uK5JOOxq/2mRYaKgFqk0loURCOCGQhiOg2GaVE+Ok6pv\noLYliuca0JEjM1BPvhgmTJYwBVwKhDm4sC6t6W1QIEyWCCEcwpZDIGxQIEBDg4lLnsmBNiZNnqie\ntMgw0qxvkUpjWbjhMF48jmdZGMUieC65ugSFWXOItNZRH/eIRKEhkMIhTJIWBqgiD2+tK1ay8/e3\nd8V6e01vFw8XiwgeCbqYwhYak6/hJbvp29RHphDEy9t0dbpowrfI8FFQi1QgL5HACwahrg4nHset\nj5OZM5d800TMhnqCVSFqIzZmJIhLGBMPlxAFgtiUgnhntjq8fTtXEShgUiRCjig7aMTFohCowauq\nxjOD1BeShKMGE2qzVNUYOAWXZHIk9u4SGZ809C1SiQwDr7kZx3He2q3KotYCZ0czbdNqoL+GaHod\nlhVgwA7TmMuTI4yLh0eOEO5b//MI4L411B3EwGaAOgaoJUqBLKVNP2zDAsfGiFZjFbPkInU0TIhi\nGi5GwCKfN3EcR0PgIsNAQS1Syd5arswCqgJpMkB4UoLeYA123qPQ2Y+XSxLN9RHEwyZAAQOPAg4B\nSleaC+SJUCBCEJcMMYqEKQZqyVZPoCfSSirSAJNacawwuVAGr6qKeJ2LGygtmeY5uodaZLgoqEXG\niPrmIP2uQSprghWC1okYgSgWDtliE2YuSNBOUyBEgBwuASIMMEA9W5hBgSjNJOkLNlBVF8QORNhS\n826yVg1BilRbRQJOhpidI+B1YmWqKTSWNuEwDN1DLTJcFNQiY4QbCBFtqcHuKGJWQ8CtpQ/wepPY\n1TV0m1GaMlso2hYGUQyggwn8zTiJXHUDvVacWdYmBupaiEZsImHoduJE7TS1ZoZiJIGNixELEW6a\nSqHOwOrvoRhv1j3UIsNIQS0yRlgWZKe/G7fvNcL5fizTJVAfIjDvXWwxm2jo2Uh3JkqmNwe2R9Bx\n6Dbr8OqbCEVCNFSFCFW30FjMEwraFCdOJdbWjVUVoDoYIlidxY3UEHjXZNJuiGxNPeF0nkCDTSKh\neakiw0VBLTJGWBZU10L3tH+iaOQxi3ky2/uxOjtocLdR/WYnXh+Y0XoGUg6eF2XANgnWmGxMHENj\nVY4NxnTqrX4iIZfYpBqy0YkEIhbJhioaJofxDItw3KDOMkk0O0RMB6PRBiM02s0XGbMU1CJjSHOz\nR3e3QzodxAiGMKbU0Ok10JDN0546Fi+aIurlidcH2eFMYyCZZiAQJ1DVQMDtoiZm4YZaSJr1uDUu\nxZCJbReprnVxg1EcxyUU8YjFIBKBgGvgaMxbZFgpqEXGEMOAI490aW83yGZNenqgaVKYmlCCmphF\nOjsJDI+OboN4ooG+jT1s72siHW0iZJtkIlEaEwaBlEGyFzwPmmodquoC4DpMmAANDR6W5WEZLm44\nrKneIsOsbEH9xBNP8Mtf/pIJEyYAcMwxx/Cxj32sXKcXkYNkGDBxokehULq/OhwGIxAl1WVAMIhh\nQH8WGk2D2kSEWtOiEApTZ4VI5QwiEY9o1GPHDoPGRpemOpf+2glU53rw8hkCjktjHLxAGC+RGO3m\niox5Ze1Rn3jiiZx77rnlPKWIHIZAwMOywJnQgruhHzOdxfUMDMfEMYPYtfWEnRBTW4PU5UNUDdg0\nN3tkMibxuMPM6TaEgtgNJobRiFOIE2st4IYs9aRFRoiGvkXGKMsq3d8MHlbAwJkzh4LZiZvKkq2y\nsKvrcLI5SCSIRj1qWpuoak8yoWGAgYBJtNojWB3GaUgQemuFUM8zcayQMlpkBJU1qNetW8f111+P\n4zh86lOfYsaMGeU8vYgcAssq7RFt2wamCY1NgDGB7Vs9chnoiibIBfqY2uhgmpDJQKi5CXNqPcVO\nl7oJBoWQhbXLMt5a2ERk5Bmed+j73qxcuZJVq1btduzkk0+mpaWF+fPn8/rrr/ODH/yAm266qWyF\nisih8zzo7IRstvRnw4BQCGIxePXVUjgDmCbU1MCcOdDTA9u2QbG4M+whkSi9PhiEt6ahiMgIGVJQ\nH4wLL7yQH/zgB5jmOy+E0NbWNhxv7wutra1qXwUbS+3buXeHtcul5dbWVjZubKNQKIV3KASdnQa2\nbWAY0NVlUCiYOI5HMAitrQ6JhIdRARtljaXPbl/UvsrW2tp6SM8v29D3Qw89RGNjI+973/vYvHkz\nsVjsgCEtIiPD2s/cr50BDaUgz+dNgsHSd/emJg/HcXBdcF2DxsbKCGmRsaZsQf2+972P5cuX84c/\n/AHXdVm6dGm5Ti0iI8BxYM8Btp0BXyhodyyR0VK2oG5sbOTrX/96uU4nIiNs11nie9IkMpHRo7Fp\nEQHeniXuursfd120O5bIKFJQi8igRMIjEPAoFg0KBSgWDQIBj0RiWOacishB0IInIjLIMEobeziO\ns9cscREZHQpqEdmLAlrEPzT0LSIi4mMKahERER9TUIuIiPiYglpERMTHFNQiIiI+pqAWERHxMQW1\niIiIjymoRUREfExBLSIi4mMKahERER9TUIuIiPiYglpERMTHFNQiIiI+pqAWERHxMQW1iIiIjymo\nRUREfExBLSIi4mMKahERER9TUIuIiPiYglpERMTHFNQiIiI+pqAWERHxMQW1iIiIjymoRUREfGzI\nQb127Vo+85nP8Nxzzw0e27hxI8uWLeNrX/sad999d1kKFBERGc+GFNTt7e08/PDDzJkzZ7fj9957\nL0uWLOG6664jk8nwwgsvlKVIERGR8WpIQR2Px7niiiuoqqoaPGbbNp2dncyePRuABQsW8PLLL5en\nShERkXEqMJQXhcPhvY719/dTXV09+Pe6ujp6enqGXpmIiIgcOKhXrlzJqlWrdju2ePFi5s2b946v\n8zzvoApobW09qOdVKrWvsql9lWsstw3UvvHkgEG9aNEiFi1adMATxWIxUqnU4N+7u7uJx+MHfF1b\nW9sBn1OpWltb1b4KpvZVrrHcNlD7Kt2hfgkp2+1ZgUCASZMm8eqrrwLw17/+9YC9bhEREXlnQ7pG\n/fzzz7NixQq2bdvG+vXreeSRR1i2bBlLlizhhz/8IZ7nMXv2bI455phy1ysiIjKuDCmo58+fz/z5\n8/c6PnnyZL7xjW8cdlEiIiJSopXJREREfExBLSIi4mMKahERER9TUIuIiPiYglpERMTHFNQiIiI+\npqAWERHxMQW1iIiIjymoRUREfExBLSIi4mMKahERER9TUIuIiPiYglpERMTHFNQiIiI+pqAWERHx\nMQW1iIiIjymoRUREfExBLSIi4mMKahERER9TUIuIiPiYglpERMTHFNQiIiI+pqAWERHxMQW1iIiI\njymoRUREfExBLSIi4mMKahERER9TUIuIiPiYglpERMTHAkN94dq1a7n55ptZunQpCxYsAOCaa64h\nn88TDocBOPfcc5k5c2Z5KhURERmHhhTU7e3tPPzww8yZM2evx5YuXcrUqVMPuzAREREZ4tB3PB7n\niiuuoKqqqtz1iIiIyC6G1KPeObS9Lw888ACpVIpJkyaxZMkSQqHQO56rtbV1KCVUDLWvsql9lWss\ntw3UvvHkgEG9cuVKVq1atduxxYsXM2/evL2e++EPf5ipU6fS0tLC3XffzaOPPspZZ531judva2s7\nxJIrR2trq9pXwdS+yjWW2wZqX6U71C8hBwzqRYsWsWjRooM62fHHHz/45wULFvDMM88cUjEiIiKy\nu7LdnuV5Htdddx3pdBoozQrXpDIREZHDM6Rr1M8//zwrVqxg27ZtrF+/nkceeYRly5axaNEivvGN\nbxCJRIjH4yxevLjc9YqIiIwrQwrq+fPnM3/+/L2On3TSSZx00kmHXZSIiIiUaGUyERERH1NQi4iI\n+JiCWkRExMcU1CIiIj6moBYREfExBbWIiIiPKahFRER8TEEtIiLiYwpqERERH1NQi4iI+JiCWkRE\nxMcU1CIiIj6moBYREfExBbWIiIiPKahFRER8TEEtIiLiYwpqERERH1NQi4iI+JiCWkRExMcU1CIi\nIj6moBYREfExBbWIiIiPKahFRER8TEEtIiLiYwpqERERH1NQi4iI+JiCWkRExMcU1CIiIj4WGMqL\nHMfhrrvuoqOjA9d1+dSnPsWRRx7Jxo0b+dGPfoRhGEydOpULL7yw3PWKiIiMK0PqUT/11FNEIhGu\nu+46LrroIu69914A7r33XpYsWcJ1111HJpPhhRdeKGuxIiIi482Qgvr9738/5557LgCxWIyBgQFs\n26azs5PZs2cDsGDBAl5++eXyVSoiIjIODWnoOxB4+2W//e1vOfnkk+nv76e6unrweF1dHT09PYdf\noYiIyDh2wKBeuXIlq1at2u3Y4sWLmTdvHo8++igbNmzgK1/5Cv39/bs9x/O8gyqgtbX1EMqtPGpf\nZVP7KtdYbhuofePJAYN60aJFLFq0aK/jq1at4rnnnuNLX/oSgUCAWCxGKpUafLy7u5t4PF7eakVE\nRMaZIV2j7ujo4A9/+ANXXHEFoVAIKA2HT5o0iVdffRWAv/71r8ybN698lYqIiIxDhnewY9S7+NnP\nfsYzzzxDU1PT4LFly5bR3t7OD3/4QzzPY/bs2Zx33nllLVZERGS8GVJQi4iIyMjQymQiIiI+pqAW\nERHxsSHdR304xsPyo2vXruXmm29m6dKlLFiwAIBrrrmGfD5POBwG4Nxzz2XmzJmjWeaQ7at9Y+nz\nA3jiiSf45S9/yYQJEwA45phj+NjHPjbKVZXHT3/6U9544w0Mw2DJkiWDixSNBWvWrOHmm29mypQp\nAEydOpXzzz9/lKs6fJs3b+bGG2/kzDPP5IwzzmDHjh0sX74c13Wpr6/nkksuIRgMjnaZQ7Zn++68\n807Wr19PbW0tAGeddRbz588f5SqH7v7772fdunW4rsvZZ5/NrFmzDunzG/Gg3nX50S1btvC9732P\nb33rW4PLj86ePZvbbruNF154geOOO26kyzts7e3tPPzww8yZM2evx5YuXcrUqVNHoary2V/7xsrn\nt6sTTzxxcAW+sWLt2rW0t7dz/fXXs3XrVu666y6uv/760S6rrI466ii++MUvjnYZZZPL5bjnnns4\n+uijB4898MADfPCDH+TEE0/kZz/7GY8//jinn376KFY5dPtqH8AnP/nJwY5AJXvllVfYsmUL119/\nPalUii9/+cvMnTv3kD6/ER/6HuvLj8bjca644gqqqqpGu5Rhsa/2jaXPb6x7+eWXec973gPA5MmT\nSafTZDKZUa5K3kkwGOSqq67abV2KNWvWsHDhQgAWLlzISy+9NFrlHbZ9tW8sOeqoo7j88ssBqK6u\nJp/PH/LnN+I96rG+/OjOoe19eeCBB0ilUkyaNIklS5YM3oNeSfbVvrH0+e1q3bp1XH/99TiOw6c+\n9SlmzJgx2iUdtt7e3t0uucRiMXp7e8fUF8utW7fyne98h4GBARYvXswxxxwz2iUdFsuysCxrt2P5\nfH5wqHTnZ1ip9tU+gEcffZSHH36Yuro6zj//fGKx2ChUd/hM0yQSiQClhcKOO+44/v73vx/S5zes\nQT3cy4+Otndq354+/OEPM3XqVFpaWrj77rt59NFHOeuss0aq1CE5lPbtqlI+v5321c6TTz6ZxYsX\nM3/+fF5//XWWL1/OTTfdNEoVDp9K+6wOZOLEiSxevJgTTzyRjo4Orr32Wu64447dOgjif6eccgq1\ntbVMnz6d3/zmNzz44INccMEFo13WYXn22WdZtWoVy5Yt49JLLz2k1w7rf71jffnR/bVvX44//vjB\nPy9YsIBnnnlmuMoqm4NtX6V+fjsdqJ1HHHEE/f39uK6LaVb2jRLxeHy3b+89PT0V9VkdSENDAyed\ndBIALS0t1NfX093dTXNz8yhXVl6RSIRCoUAoFKq4n7eDMXfu3ME/L1y4kLvvvnsUqzl8L774Ir/6\n1a/46le/SlVV1SF/fiP+/zrjcflRz/O47rrrSKfTQGlCT6VPKtvVWPz8HnroIf74xz8CpRmpsVis\n4kMa4Nhjj2X16tUArF+/nng8TjQaHeWqyufpp59mxYoVQGmYv6+vj4aGhlGuqvzmzp07+DmuXr26\n4n/e9vTd736Xjo4OoHQ9fucs/kqUyWS4//77ufLKK6mpqQEO/fMb8ZXJxvryo88//zwrVqxg27Zt\nxGIx4vE4y5Yt45lnnuGhhx4iEokQj8dZunTpO17P9qv9tW/r1q1j4vPbqaura/D2Cdd1Oe+888bM\nbUz//d//zbp16zAMgwsuuIDp06ePdkllk81mue2228hkMti2zb/9279V9G09UPpCdd9995FMJrEs\ni4aGBi699FLuvPNOisUiTU1NXHzxxRU7vL+v9p1xxhk89NBDhEIhIpEIF198MXV1daNd6pA89thj\nPPjgg0ycOHHw2Gc/+1m+//3vH/TnpyVERUREfKzyx/JERETGMAW1iIiIjymoRUREfExBLSIi4mMK\nahERER9TUIuIiPiYglpERMTHFNQiIiI+9v8BlEb5bDTXr0YAAAAASUVORK5CYII=\n",
            "text/plain": [
              "<matplotlib.figure.Figure at 0x7fc422f9de48>"
            ]
          },
          "metadata": {
            "tags": []
          }
        }
      ]
    },
    {
      "metadata": {
        "id": "rWrL1mCNIcNy",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "## Acknowledgements\n",
        "\n",
        "This tutorial was originally written in Edward 1.0 ([source](https://github.com/blei-lab/edward/blob/master/notebooks/probabilistic_pca.ipynb)). We thank all contributors to writing and revising that version."
      ]
    },
    {
      "metadata": {
        "id": "ZJd8RbQxvyp6",
        "colab_type": "text"
      },
      "cell_type": "markdown",
      "source": [
        "#### References\n",
        "\n",
        "<a name='1'></a>[1]: Michael E. Tipping and Christopher M. Bishop. Probabilistic principal component analysis. _Journal of the Royal Statistical Society: Series B (Statistical Methodology)_, 61(3): 611-622, 1999."
      ]
    }
  ]
}