nf-tutorial-jax.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Real-NVP\n",
    "\n",
    "Real-NVP implementation in less than 100 lines of [JAX](https://github.com/google/jax)\n",
    "\n",
    "Author: Eric Jang"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# !pip install --upgrade jax jaxlib scikit-learn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import jax.numpy as np\n",
    "from jax import random"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data\n",
    "\n",
    "$p(x)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-2, 2)"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from sklearn import cluster, datasets, mixture\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "import matplotlib.pyplot as plt\n",
    "n_samples = 2000\n",
    "noisy_moons = datasets.make_moons(n_samples=n_samples, noise=.05)\n",
    "X, y = noisy_moons\n",
    "X = StandardScaler().fit_transform(X)\n",
    "xlim, ylim = [-2, 2], [-2, 2]\n",
    "plt.scatter(X[:, 0], X[:, 1], s=10, color='red')\n",
    "plt.xlim(xlim)\n",
    "plt.ylim(ylim)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/usr/local/Cellar/python/3.7.3/Frameworks/Python.framework/Versions/3.7/lib/python3.7/site-packages/jax/lib/xla_bridge.py:130: UserWarning: No GPU/TPU found, falling back to CPU.\n",
      "  warnings.warn('No GPU/TPU found, falling back to CPU.')\n"
     ]
    }
   ],
   "source": [
    "# Initialize parameters, not committing to a batch shape\n",
    "rng = random.PRNGKey(0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# N(1, 1)\n",
    "# X = random.normal(rng, (n_samples, 2)) + 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Base Distribution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def sample_n01(N):\n",
    "  D = 2\n",
    "  return random.normal(rng, (N, D))\n",
    "\n",
    "def log_prob_n01(x):\n",
    "  return np.sum(-np.square(x)/2 - np.log(np.sqrt(2*np.pi)),axis=-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([ 11., 110., 465., 462., 152.,  98.,  89., 109., 157., 347.]),\n",
       " array([-3.796753 , -3.608996 , -3.421239 , -3.233482 , -3.0457249,\n",
       "        -2.8579679, -2.6702106, -2.4824536, -2.2946966, -2.1069396,\n",
       "        -1.9191825], dtype=float32),\n",
       " <a list of 10 Patch objects>)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD8CAYAAAB5Pm/hAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAD6xJREFUeJzt3X2MnNV1x/HvCQbSkgQbbCHHdrMkcRshVQXqIlKiKgU15aXCrkRSaFpcYslSBSqIVokpUkvVRoJWDWnUCmLFCEOjJJQ0wg1E4ABRVSmQLgmY97BQI3AhOLwGIdqSnP4xd5thu+uZ3Z3ZGY6+H2m097nPnXnOXs/89pk7L47MRJJU19tGXYAkabgMekkqzqCXpOIMekkqzqCXpOIMekkqzqCXpOIMekkqzqCXpOKWjboAgJUrV+bExMSoy5Ckt5R77rnnh5m5qte4sQj6iYkJJicnR12GJL2lRMST/Yxz6UaSijPoJak4g16SijPoJak4g16SijPoJak4g16SijPoJak4g16SihuLT8Zqfia23TyyY++9/IyRHVvSwhj0mpdR/ZHxD4y0cC7dSFJxBr0kFWfQS1JxBr0kFWfQS1JxBr0kFWfQS1JxBr0kFWfQS1JxBr0kFWfQS1JxBr0kFWfQS1JxBr0kFWfQS1JxBr0kFWfQS1JxBr0kFWfQS1JxBr0kFWfQS1JxBr0kFWfQS1JxfQd9RBwUEd+LiK+37aMj4u6ImIqIr0TEIa3/0LY91fZPDKd0SVI/5nNGfyHwcNf2FcCVmfl+4EVgS+vfArzY+q9s4yRJI9JX0EfEWuAM4AttO4CTgRvbkJ3Aptbe2LZp+09p4yVJI9DvGf1ngU8CP2nbRwIvZeYbbftpYE1rrwGeAmj7X27jJUkj0DPoI+K3gOcy855BHjgitkbEZERM7t+/f5A3LUnq0s8Z/UnAmRGxF/gynSWbvwOWR8SyNmYtsK+19wHrANr+w4HnZ95oZm7PzA2ZuWHVqlWL+iUkSXPrGfSZeUlmrs3MCeBs4I7M/DhwJ3BWG7YZuKm1d7Vt2v47MjMHWrUkqW+LeR/9p4CLI2KKzhr8jta/Aziy9V8MbFtciZKkxVjWe8hPZea3gG+19hPACbOMeR346ABqkyQNgJ+MlaTiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKm5e/5WgJFU0se3mkR177+VnDP0YntFLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnE9gz4i3h4R34mI+yLiwYj4i9Z/dETcHRFTEfGViDik9R/atqfa/onh/gqSpAPp54z+v4CTM/OXgGOBUyPiROAK4MrMfD/wIrCljd8CvNj6r2zjJEkj0jPos+PVtnlwuyRwMnBj698JbGrtjW2btv+UiIiBVSxJmpe+1ugj4qCIuBd4DtgNPA68lJlvtCFPA2taew3wFEDb/zJw5Cy3uTUiJiNicv/+/Yv7LSRJc+or6DPzx5l5LLAWOAH4wGIPnJnbM3NDZm5YtWrVYm9OkjSHeb3rJjNfAu4EPggsj4hlbddaYF9r7wPWAbT9hwPPD6RaSdK89fOum1URsby1fwb4DeBhOoF/Vhu2GbiptXe1bdr+OzIzB1m0JKl/y3oPYTWwMyIOovOH4YbM/HpEPAR8OSL+CvgesKON3wFcHxFTwAvA2UOoW5LUp55Bn5l7gONm6X+Cznr9zP7XgY8OpDpJ0qL5yVhJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKs6gl6TiDHpJKq5n0EfEuoi4MyIeiogHI+LC1n9EROyOiMfazxWtPyLicxExFRF7IuL4Yf8SkqS59XNG/wbwx5l5DHAicH5EHANsA27PzPXA7W0b4DRgfbtsBa4aeNWSpL71DPrMfCYzv9vaPwIeBtYAG4GdbdhOYFNrbwSuy467gOURsXrglUuS+jKvNfqImACOA+4GjsrMZ9quZ4GjWnsN8FTX1Z5ufZKkEeg76CPiHcBXgYsy85XufZmZQM7nwBGxNSImI2Jy//7987mqJGke+gr6iDiYTsh/MTP/uXX/YHpJpv18rvXvA9Z1XX1t63uTzNyemRsyc8OqVasWWr8kqYd+3nUTwA7g4cz8TNeuXcDm1t4M3NTVf257982JwMtdSzySpCW2rI8xJwG/D9wfEfe2vj8FLgduiIgtwJPAx9q+W4DTgSngNeC8gVYsSZqXnkGfmf8GxBy7T5llfALnL7IuSdKA+MlYSSrOoJek4gx6SSrOoJek4gx6SSrOoJek4gx6SSrOoJek4gx6SSrOoJek4gx6SSrOoJek4gx6SSrOoJek4vr5PnpJWhIT224edQkleUYvScV5Rq+3hFGd6e29/IyRHFcaJM/oJak4g16SijPoJak4g16SijPoJak4g16SijPoJak4g16SijPoJak4g16SijPoJak4g16SijPoJak4g16SijPoJak4g16SijPoJak4g16SivO/ElwE/yNjVeT9up6eZ/QRcU1EPBcRD3T1HRERuyPisfZzReuPiPhcRExFxJ6IOH6YxUuSeutn6eZa4NQZfduA2zNzPXB72wY4DVjfLluBqwZTpiRpoXoGfWb+K/DCjO6NwM7W3gls6uq/LjvuApZHxOpBFStJmr+Fvhh7VGY+09rPAke19hrgqa5xT7e+/ycitkbEZERM7t+/f4FlSJJ6WfSLsZmZEZELuN52YDvAhg0b5n19aSmM8oXJvZefMbJjq5aFntH/YHpJpv18rvXvA9Z1jVvb+iRJI7LQoN8FbG7tzcBNXf3ntnffnAi83LXEI0kagZ5LNxHxJeDDwMqIeBr4c+By4IaI2AI8CXysDb8FOB2YAl4DzhtCzZKkeegZ9Jl5zhy7TpllbALnL7YoSX5wSYPjVyBIUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVZ9BLUnEGvSQVt2zUBSzWxLabR12CJI21oZzRR8SpEfFoRExFxLZhHEOS1J+BB31EHAT8A3AacAxwTkQcM+jjSJL6M4wz+hOAqcx8IjP/G/gysHEIx5Ek9WEYQb8GeKpr++nWJ0kagZG9GBsRW4GtbfPViHh0SIdaCfxwSLc9CONeH4x/jeNeH1jjIIx7fbCAGuOKRR3vPf0MGkbQ7wPWdW2vbX1vkpnbge1DOP6bRMRkZm4Y9nEWatzrg/GvcdzrA2schHGvD8a3xmEs3fw7sD4ijo6IQ4CzgV1DOI4kqQ8DP6PPzDci4gLgVuAg4JrMfHDQx5Ek9Wcoa/SZeQtwyzBuewGGvjy0SONeH4x/jeNeH1jjIIx7fTCmNUZmjroGSdIQ+V03klRciaCPiL+MiD0RcW9E3BYR755lzK+3/dOX1yNiU9t3bUT8R9e+Y5e6vjbux1017OrqPzoi7m5fKfGV9iL3QPU5h8dGxLcj4sE29ne69o3LHG6OiMfaZXNX/y9HxP1tDj8XETHI+tox/iYiHml1fi0ils8y5hdm3A9fiYiL2r7LImJf177Tl7q+Nm5vm6t7I2Kyq/+IiNjd5nZ3RKwYZH391hgR6yLizoh4qN0XL+zaNy5zOOvXwCzFY3lWmfmWvwDv6mr/EXB1j/FHAC8AP9u2rwXOGnV9wKtz9N8AnN3aVwN/OIoagZ8H1rf2u4FngOXjMoft3/WJ9nNFa69o+74DnAgE8A3gtCHU+BFgWWtfAVzRY/xBwLPAe9r2ZcCfDHEO+6oP2AusnKX/r4Ftrb2t1+83rBqB1cDxrf1O4PvAMeMyh+3f9XHgvcAhwH1d9Q39sTzbpcQZfWa+0rV5GNDrhYezgG9k5mvDq+qnFlDf/2lnnicDN7auncCmwVXX0U+Nmfn9zHystf8TeA5YNehaFlof8JvA7sx8ITNfBHYDp0bEajp/KO7KziPsOoYzh7dl5htt8y46nyE5kFOAxzPzyUHXMpsF1DfTRjr3Pxje/bBnjZn5TGZ+t7V/BDzMEn36vs85nPVrYJbqsTybEkEPEBGfjoingI8Df9Zj+NnAl2b0fbo9HbsyIg4dUX1vj4jJiLhrelkJOBJ4qevONbSvlJjPHEbECXTOVh7v6h71HM719RtrWntm/zB9gs4zhwOZ7X54QZvDa4axNNLlQPUlcFtE3BOdT7BPOyozn2ntZ4Gjhlgf9DGHETEBHAfc3dU96jmc6364ZI/lmd4yQR8R34yIB2a5bATIzEszcx3wReCCA9zOauAX6bzPf9olwAeAX6HztP9TI6rvPdn5VN3vAp+NiPfNt44lqHF6Dq8HzsvMn7TucZnDoepVYxtzKfBGq3Ou2zkEOBP4p67uq4D3AcfSWRb72xHV96HMPJ7ON9CeHxG/NnNAe2a0oLfsDXAO3wF8Fbio6xnfuMzheFmK9aGlvAA/BzxwgP0XAtsPsP/DwNdHVV/XuGvpLDEFne/OmF4X/CBw66jmEHgX8F0OsB4/qjkEzgE+37X9+da3GnhkrnEDru0PgG/TXv85wLiNwG0H2D/Rz/1kWPV1jb+MtuYNPAqsbu3VwKOjmkPgYDonaxeP2xzOfIzSOQm6ZBSP5enLW+aM/kAiYn3X5kbgkQMMP4cZT5fbGer0evgm4IGlri8iVkwvd0TESuAk4KHs3CPupBP6AJuBmwZZ3zxqPAT4GnBdZt44Y9/I55DOA/8jbS5X0Hnh7NbsLDe8EhEntvrOZThzeCrwSeDM7P36z5z3w+a3Gfwc9qwvIg6LiHdOt+nM4XQdu+jc/2B498N+agxgB/BwZn5mxr6RzyFzfA3MUj2WZ7UUf02GfaHz9O0BYA/wL8Ca1r8B+MKMv/D7gLfNuP4dwP3tNv4ReMdS1wf8aqvhvvZzS9f130vnXSNTdJ7qHzqKOQR+D/gf4N6uy7HjModt+xNtnqboLC3RNe4BOq8p/D3tw4IDrnGKztrs9Nxc3frfDdzSNe4w4Hng8BnXv77N4R46obp6qetr97X72uVB4NKu6x8J3A48BnwTOGIUcwh8iM6y0Z6ucaePyxy27dPpvBvo8RlzOPTH8mwXPxkrScWVWLqRJM3NoJek4gx6SSrOoJek4gx6SSrOoJek4gx6SSrOoJek4v4Xlb5SBxs5FPEAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.hist(log_prob_n01(X))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Real-NVP\n",
    "\n",
    "Affine coupling layer\n",
    "\n",
    "$x_2=x_1$\n",
    "\n",
    "$y_2=x_1*e^{s(x_1)}+t(x_1)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "from jax.experimental import stax # neural network library\n",
    "from jax.experimental.stax import Dense, Relu # neural network layers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "def nvp_forward(net_params, shift_and_log_scale_fn, x, flip=False):\n",
    "  d = x.shape[-1]//2\n",
    "  x1, x2 = x[:, :d], x[:, d:]\n",
    "  if flip:\n",
    "    x2, x1 = x1, x2\n",
    "  shift, log_scale = shift_and_log_scale_fn(net_params, x1)\n",
    "  y2 = x2*np.exp(log_scale) + shift\n",
    "  if flip:\n",
    "    x1, y2 = y2, x1\n",
    "  y = np.concatenate([x1, y2], axis=-1)\n",
    "  return y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def nvp_inverse(net_params, shift_and_log_scale_fn, y, flip=False):\n",
    "  d = y.shape[-1]//2\n",
    "  y1, y2 = y[:, :d], y[:, d:]\n",
    "  if flip:\n",
    "    y1, y2 = y2, y1\n",
    "  shift, log_scale = shift_and_log_scale_fn(net_params, y1)\n",
    "  x2 = (y2-shift)*np.exp(-log_scale)\n",
    "  if flip:\n",
    "    y1, x2 = x2, y1\n",
    "  x = np.concatenate([y1, x2], axis=-1)\n",
    "  return x, log_scale"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "$\\log p(x) = \\log p(f^{-1}(x)) + |\\log det J(f^{-1})(x)| $"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "def init_nvp():\n",
    "  D = 2\n",
    "  net_init, net_apply = stax.serial(\n",
    "    Dense(512), Relu, Dense(512), Relu, Dense(D))\n",
    "  in_shape = (-1, D//2)\n",
    "  out_shape, net_params = net_init(rng, in_shape)\n",
    "  def shift_and_log_scale_fn(net_params, x1):\n",
    "    s = net_apply(net_params, x1)\n",
    "    return np.split(s, 2, axis=1)\n",
    "  return net_params, shift_and_log_scale_fn\n",
    "\n",
    "def sample_nvp(net_params, shift_log_scale_fn, base_sample_fn, N, flip=False):\n",
    "  x = base_sample_fn(N)\n",
    "  return nvp_forward(net_params, shift_log_scale_fn, x, flip=flip)\n",
    "\n",
    "def log_prob_nvp(net_params, shift_log_scale_fn, base_log_prob_fn, y, flip=False):\n",
    "  x, log_scale = nvp_inverse(net_params, shift_log_scale_fn, y, flip=flip)\n",
    "  ildj = -np.sum(log_scale, axis=-1)\n",
    "  return base_log_prob_fn(x) + ildj"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "def init_nvp_chain(n=2):\n",
    "  flip = False\n",
    "  ps, configs = [], []\n",
    "  for i in range(n):\n",
    "    p, f = init_nvp()\n",
    "    ps.append(p), configs.append((f, flip))\n",
    "    flip = not flip\n",
    "  return ps, configs\n",
    "\n",
    "def sample_nvp_chain(ps, configs, base_sample_fn, N):\n",
    "  x = base_sample_fn(N)\n",
    "  for p, config in zip(ps, configs):\n",
    "    shift_log_scale_fn, flip = config\n",
    "    x = nvp_forward(p, shift_log_scale_fn, x, flip=flip)\n",
    "  return x\n",
    "\n",
    "def make_log_prob_fn(p, log_prob_fn, config):\n",
    "  shift_log_scale_fn, flip = config\n",
    "  return lambda x: log_prob_nvp(p, shift_log_scale_fn, log_prob_fn, x, flip=flip)\n",
    "\n",
    "def log_prob_nvp_chain(ps, configs, base_log_prob_fn, y):\n",
    "  log_prob_fn = base_log_prob_fn\n",
    "  for p, config in zip(ps, configs):\n",
    "    log_prob_fn = make_log_prob_fn(p, log_prob_fn, config)\n",
    "  return log_prob_fn(y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "ps, cs = init_nvp_chain(4)\n",
    "y = sample_nvp_chain(ps, cs, sample_n01, 1000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([ 14.,  62.,  77., 169., 371., 546., 208., 116., 195., 242.]),\n",
       " array([-4.383059 , -4.1515055, -3.9199514, -3.688398 , -3.456844 ,\n",
       "        -3.2252903, -2.9937365, -2.7621827, -2.5306292, -2.2990754,\n",
       "        -2.0675216], dtype=float32),\n",
       " <a list of 10 Patch objects>)"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD8CAYAAAB0IB+mAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAADtNJREFUeJzt3W2MXFd9x/HvrwmBClKc4MUNtoVTYYlSVTx0lVLBC0oKTQiq8wIQCBVDLVmVoKUCCQyRilpayVGlpkWtqCyCairaENFGsUhK4+ZBCFWBrCFPJLRZUke2leCFhFCEoAr998Uctxuzm5lZz3jsPd+PNJpzzz0z978nN7/c3Jl7J1WFJGn9+5lZFyBJOj0MfEnqhIEvSZ0w8CWpEwa+JHXCwJekThj4ktQJA1+SOmHgS1Inzp11AQAbN26sbdu2zboMSTqrHDp06DtVNTfq+DMi8Ldt28bCwsKsy5Cks0qSR8YZ7ykdSeqEgS9JnTDwJakTBr4kdcLAl6ROGPiS1AkDX5I6YeBLUicMfEnqxBlxpa10ptq256aZbfvw3itmtm2tTx7hS1InDHxJ6oSBL0mdMPAlqRMGviR1wsCXpE4Y+JLUCQNfkjph4EtSJwx8SeqEgS9JnTDwJakTBr4kdWKkwE9yOMl9Se5OstD6LkxyMMlD7fmC1p8kn0iymOTeJK+a5h8gSRrNOEf4v15Vr6iq+ba8B7i1qrYDt7ZlgMuB7e2xG/jkpIqVJK3dqZzS2QHsb+39wJXL+j9TA3cCG5JcdArbkSRNwKiBX8AtSQ4l2d36NlXVo639GLCptTcDR5a99mjre5oku5MsJFlYWlpaQ+mSpHGM+otXr62qY0leCBxM8s3lK6uqktQ4G66qfcA+gPn5+bFeK0ka30hH+FV1rD0fB24ALgG+feJUTXs+3oYfA7Yue/mW1idJmqGhgZ/kuUnOP9EG3gjcDxwAdrZhO4EbW/sA8K72bZ1XA08uO/UjSZqRUU7pbAJuSHJi/N9X1ReT3AVcn2QX8Ajwtjb+ZuBNwCLwQ+A9E69akjS2oYFfVQ8DL1+h/7vApSv0F/DeiVQnSZoYr7SVpE4Y+JLUCQNfkjph4EtSJwx8SeqEgS9JnTDwJakTBr4kdcLAl6ROGPiS1AkDX5I6YeBLUicMfEnqhIEvSZ0w8CWpEwa+JHXCwJekThj4ktQJA1+SOmHgS1InDHxJ6oSBL0mdMPAlqRMGviR1wsCXpE4Y+JLUCQNfkjph4EtSJ0YO/CTnJPl6ki+05YuTfCXJYpLPJTmv9T+7LS+29dumU7okaRzjHOG/H3hw2fLVwDVV9RLgCWBX698FPNH6r2njJEkzNlLgJ9kCXAF8qi0HeD3w+TZkP3Bla+9oy7T1l7bxkqQZGvUI/y+ADwH/05ZfAHyvqp5qy0eBza29GTgC0NY/2cZLkmZoaOAneTNwvKoOTXLDSXYnWUiysLS0NMm3liStYJQj/NcAv5XkMHAdg1M5fwlsSHJuG7MFONbax4CtAG3984HvnvymVbWvquaran5ubu6U/ghJ0nBDA7+qPlJVW6pqG/B24LaqeidwO/CWNmwncGNrH2jLtPW3VVVNtGpJ0thO5Xv4HwY+kGSRwTn6a1v/tcALWv8HgD2nVqIkaRLOHT7k/1XVHcAdrf0wcMkKY34EvHUCtUmSJsgrbSWpEwa+JHXCwJekThj4ktQJA1+SOmHgS1InxvpapjQr2/bcNOsSpLOeR/iS1AkDX5I6YeBLUicMfEnqhIEvSZ0w8CWpEwa+JHXCwJekThj4ktQJA1+SOmHgS1InDHxJ6oSBL0mdMPAlqRMGviR1wsCXpE4Y+JLUCQNfkjph4EtSJwx8SeqEgS9JnTDwJakTQwM/yXOSfDXJPUm+keSPWv/FSb6SZDHJ55Kc1/qf3ZYX2/pt0/0TJEmjGOUI/8fA66vq5cArgMuSvBq4Grimql4CPAHsauN3AU+0/mvaOEnSjA0N/Br4QVt8VnsU8Hrg861/P3Bla+9oy7T1lybJxCqWJK3JSOfwk5yT5G7gOHAQ+Bbwvap6qg05Cmxu7c3AEYC2/kngBZMsWpI0vpECv6p+UlWvALYAlwAvPdUNJ9mdZCHJwtLS0qm+nSRpiLG+pVNV3wNuB34N2JDk3LZqC3CstY8BWwHa+ucD313hvfZV1XxVzc/Nza2xfEnSqEb5ls5ckg2t/bPAG4AHGQT/W9qwncCNrX2gLdPW31ZVNcmiJUnjO3f4EC4C9ic5h8F/IK6vqi8keQC4LsmfAF8Hrm3jrwX+Lski8Djw9inULUka09DAr6p7gVeu0P8wg/P5J/f/CHjrRKqTJE2MV9pKUicMfEnqhIEvSZ0w8CWpEwa+JHXCwJekThj4ktQJA1+SOmHgS1InDHxJ6oSBL0mdMPAlqRMGviR1wsCXpE4Y+JLUCQNfkjph4EtSJ0b5iUNJWte27blpZts+vPeK07Ytj/AlqRMGviR1wsCXpE4Y+JLUCQNfkjph4EtSJwx8SeqEgS9JnTDwJakTBr4kdcLAl6RODA38JFuT3J7kgSTfSPL+1n9hkoNJHmrPF7T+JPlEksUk9yZ51bT/CEnScKPcPO0p4INV9bUk5wOHkhwE3g3cWlV7k+wB9gAfBi4HtrfHrwKfbM+SxjCrG3qdzpt56fQaeoRfVY9W1dda+7+AB4HNwA5gfxu2H7iytXcAn6mBO4ENSS6aeOWSpLGMdQ4/yTbglcBXgE1V9Whb9RiwqbU3A0eWvexo65MkzdDIgZ/kecA/An9QVd9fvq6qCqhxNpxkd5KFJAtLS0vjvFSStAYjBX6SZzEI+89W1T+17m+fOFXTno+3/mPA1mUv39L6nqaq9lXVfFXNz83NrbV+SdKIhn5omyTAtcCDVfXny1YdAHYCe9vzjcv635fkOgYf1j657NSPJK1qlr881YNRvqXzGuC3gfuS3N36Psog6K9Psgt4BHhbW3cz8CZgEfgh8J6JVixJWpOhgV9VXwayyupLVxhfwHtPsS5J0oR5pa0kdcLAl6ROGPiS1AkDX5I6YeBLUicMfEnqhIEvSZ0Y5cIr6f94JaR09vIIX5I6YeBLUicMfEnqhIEvSZ0w8CWpEwa+JHXCwJekThj4ktQJA1+SOmHgS1InDHxJ6oSBL0mdMPAlqRMGviR1wsCXpE4Y+JLUCQNfkjph4EtSJwx8SeqEv2kr6Wn83eL1yyN8SerE0MBP8ukkx5Pcv6zvwiQHkzzUni9o/UnyiSSLSe5N8qppFi9JGt0oR/h/C1x2Ut8e4Naq2g7c2pYBLge2t8du4JOTKVOSdKqGBn5VfQl4/KTuHcD+1t4PXLms/zM1cCewIclFkypWkrR2az2Hv6mqHm3tx4BNrb0ZOLJs3NHWJ0masVP+0LaqCqhxX5dkd5KFJAtLS0unWoYkaYi1Bv63T5yqac/HW/8xYOuycVta30+pqn1VNV9V83Nzc2ssQ5I0qrUG/gFgZ2vvBG5c1v+u9m2dVwNPLjv1I0maoaEXXiX5B+B1wMYkR4GPAXuB65PsAh4B3taG3wy8CVgEfgi8Zwo1S5LWYGjgV9U7Vll16QpjC3jvqRYlSZo8r7SVpE4Y+JLUCQNfkjph4EtSJ7w98lnI29dKWguP8CWpEwa+JHXCwJekThj4ktQJP7Q9BX54Kuls4hG+JHXCwJekThj4ktQJA1+SOmHgS1InDHxJ6oSBL0mdMPAlqRMGviR1wsCXpE4Y+JLUCQNfkjph4EtSJwx8SeqEgS9JnTDwJakTBr4kdeKs/8Urf3VKkkbjEb4kdWIqgZ/ksiT/nmQxyZ5pbEOSNJ6JB36Sc4C/Bi4HXga8I8nLJr0dSdJ4pnGEfwmwWFUPV9V/A9cBO6awHUnSGKYR+JuBI8uWj7Y+SdIMzexbOkl2A7vb4o+T3D+rWs4gG4HvzLqIGXMOBpyHgXU/D7l6pGGrzcOLx9nWNAL/GLB12fKW1vc0VbUP2AeQZKGq5qdQy1nFeXAOTnAeBpyHgUnNwzRO6dwFbE9ycZLzgLcDB6awHUnSGCZ+hF9VTyV5H/AvwDnAp6vqG5PejiRpPFM5h19VNwM3j/GSfdOo4yzkPDgHJzgPA87DwETmIVU1ifeRJJ3hvLWCJHVipoGf5INJKsnGVdb/JMnd7bFuP/gdYR52JnmoPXae7vqmKcnHk9zb/hnfkuRFq4xb1/vCGPOwbvcFgCR/luSbbS5uSLJhlXGHk9zX5mvhdNc5bWPMw1i3sZnZKZ0kW4FPAS8FfqWqfuo7pkl+UFXPO+3FnUbD5iHJhcACMA8UcKiNe+J01zoNSX6uqr7f2r8PvKyqfneFcet6XxhlHtb7vgCQ5I3Abe3LH1cDVNWHVxh3GJhfKTfWg1Hmod3G5j+ANzC4wPUu4B1V9cBq7zvLI/xrgA8x2HF7NmwefhM4WFWPt3+xDwKXna7ipu1EyDXPpdP9YcR5WNf7AkBV3VJVT7XFOxlcx9OdEedh7NvYzCTwk+wAjlXVPUOGPifJQpI7k1x5Omo7nUach3V/q4okf5rkCPBO4A9XGbau9wUYaR7W/b5wkt8B/nmVdQXckuRQu2p/PVttHsbeH6Z2a4Uk/wr8/AqrrgI+CrxxhLd5cVUdS/ILwG1J7quqb02yzmmb0Dyc1Z5pDqrqxqq6CrgqyUeA9wEfW2Hsut4XxpiHs96weWhjrgKeAj67ytu8tu0PLwQOJvlmVX1pOhVPx4TmYSxTC/yq+o2V+pP8MnAxcE8SGPyvyteSXFJVj530Hsfa88NJ7gBeCZxV/5JPYB6OAa9btrwFuGMqxU7JanOwgs8yuH7jp4JuPe8LK1htHs76fQGGz0OSdwNvBi6tVT5kXLY/HE9yA4PTG2dV4E9gHka6jc3JG53pAzgMbFyh/wLg2a29EXiIwQdZM6/5NM/DhcB/tvm4oLUvnHW9E/y7ty9r/x7w+R73hRHnYV3vC+1vvAx4AJh7hjHPBc5f1v434LJZ1z6DeTgXeJjBgeN5wD3ALz3T+55R38NPMp/kU23xF4GFJPcAtwN76xk+fV5Pls9DVT0OfJzBJ/B3AX/c+taLvUnuT3Ivg9Nb74cu94Wh89DBvgDwV8D5DE7T3J3kbwCSvCjJiav3NwFfbvvDV4GbquqLsyl3aobOQw0+1D1xG5sHgetryG1svNJWkjpxRh3hS5Kmx8CXpE4Y+JLUCQNfkjph4EtSJwx8SeqEgS9JnTDwJakT/wv4YB6vML9KZQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.hist(log_prob_nvp_chain(ps, cs, log_prob_n01, X))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Train It"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "from jax.experimental import optimizers\n",
    "from jax import jit, grad\n",
    "import numpy as onp\n",
    "ps, cs = init_nvp_chain(4)\n",
    "\n",
    "def loss(params, batch):\n",
    "  return -np.mean(log_prob_nvp_chain(params, cs, log_prob_n01, batch))\n",
    "opt_init, opt_update, get_params = optimizers.adam(step_size=1e-4)\n",
    "\n",
    "@jit\n",
    "def step(i, opt_state, batch):\n",
    "  params = get_params(opt_state)\n",
    "  g = grad(loss)(params, batch)\n",
    "  return opt_update(i, g, opt_state)\n",
    "\n",
    "iters = int(1e4)\n",
    "data_generator = (X[onp.random.choice(X.shape[0], 100)] for _ in range(iters))\n",
    "opt_state = opt_init(ps)\n",
    "for i in range(iters):\n",
    "  opt_state = step(i, opt_state, next(data_generator))\n",
    "ps = get_params(opt_state)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "y = sample_nvp_chain(ps, cs, sample_n01, 1000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "DeviceArray(2.8329544, dtype=float32)"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "y.max()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJztnX+MZUd157/n9fN73t4h7XhnGH5l2maIgKQBmx4xZkC7VuxtrOeoweHZ2lGSDdqs3NH2SBkJy4ZuyZEsDWCvQBYZNh4w8ZIIzcY29nZv/GsYBwuhxLPuIWM8YJtgOwZbMEwW0WC11E131/5Rr3Lr1quqW/fdH+++985Huno/7q+6de89derUOadICAGGYRhmeKj1uwAMwzBMvrBgZxiGGTJYsDMMwwwZLNgZhmGGDBbsDMMwQwYLdoZhmCGDBTvDMMyQwYKdYRhmyGDBzjAMM2TUsx6AiC4E8E0Azc7x7hdC/Jlvn507d4pLLrkk66kZhmFGitOnT/+LEGJX0naZBTuAdQC/I4R4jYguAPAtInpECPGka4dLLrkEKysrOZyaYRhmdCCil0O2yyzYhUw281rn5wWdhRPQMAzD9IlcbOxENEZEZwD8FMDXhRCnLNvcSEQrRLRy/vz5PE7LMAzDWMhFsAshtoQQlwF4C4D3EdGUZZsvCiH2CSH27dqVaCJiGIZheiRXrxghxM8BfAPANXkel2EYhgkns2Anol1EdFHn+78B8B8BPJf1uAzDMExv5OEV80YAXyGiMciG4l4hxN/mcFyGYRimB/LwivkOgMtzKAvDMAyTAxx5yjB5srwMHDokPxmmT7BgZ5i8WF4GDh4EvvAF+cnCnekTLNgZJi9OnADW1uT3tTX5m2H6AAt2hsmLmRlgfFx+Hx+XvxmmD+ThFcMwDADMzgLHj0tNfWZG/maYPsCCnWHyZHaWBTrTd9gUwzAMM2SwYGcYhhkyWLAzDMMMGSzYGYZhhgwW7AzDMEMGC3aGYZghgwU7wzDMkMGCnWEYZshgwc4wWeGMjkzFYMHOMFngjI5MBWHBzjBZ4IyOTAVhwc4wWeCMjkwF4SRgDJMFzujIVBAW7AyTFc7oyFQMNsUwDMMMGSzYGYZhhgwW7AzDMEMGC3aGYZghgwU7wzDMkMGCnWHSwOkDmAGABTvDhMLpA5gBgQU7w4TC6QOYAYEFO8OEwukDmAGBI08ZJhROH8AMCCzYGSYNnD6AGQDYFMMwDDNksGBnGIYZMjILdiL6DSL6BhF9j4i+S0R/mkfBGIZhmN7Iw8a+CeDjQohvE9HrAJwmoq8LIb6Xw7EZhmGYlGTW2IUQPxZCfLvz/ZcAngXw5qzHZRiGYXojVxs7EV0C4HIAp/I8LsNUAk4nwAwIuQl2ItoB4GsADgshfmFZfyMRrRDRyvnz5/M6LcOUA6cTYAaIXAQ7EV0AKdS/KoR4wLaNEOKLQoh9Qoh9u3btyuO0DFMenE6AGSDy8IohAF8G8KwQ4nPZi8QwFYTTCTADRB5eMR8A8IcAniGiM53/FoQQD+dwbIapBpxOgBkgMgt2IcS3AFAOZWGYasPpBJgBgSNPGYZhhgwW7AzDMEMGC3aGYZghgwU7M3zkHUjEgUnMgMGCncmXfgvBvAOJODCJGUBYsDP50asQzLMxyDuQiAOTmAGEBfsok7d23YsQ1BuDG24Arr3WXp7QsuYdSMSBScwgIoQofZmenhZMn1laEmJ8XAhAfi4t+bedn/dvk/aYivl5ub2+mPumPW5oeUPJ+3gM0yMAVkSAjGWNfVQJ1a6Xl6UmrTRqn8asojPn54HDh+UxFxfjmrapeU9MdB/HLE/ansDsLHD0qPyeR49EHY+Dk5gBgSezHiWWl6OQ+JkZ4J57pKC0mRiWl4Fjx4CnngLW1+V/6+vyP0AeZ2ICWF2Nh9irz4MHI2EMyHMdPgzceaf8/557ZCOwutpdTrM8SWV1Xasqw113AbfcAhw5ElZP+jE4hQAziISo9XkvbIrpAzZzhsvEsLQkRLPZbSIBhJiejo6jlmZTiFYrOo7NvAIIMTUV/z05KcTevUKMjcnfjUb8OHr50ppDzDLU6+lMKb2YlRimYMCmGCaGzZxhmhiUmeTYsUhL12k0gN2745o4ILd9+GHg+uvlMfQBR513vCP+/8svAy+8AGxtAbWaLMell0Zl0T1sTnXmbjl1KnwQdWws+r25mc6jxawv1VNhmAGATTHDhst8YNqyzd+66aJel4J2eztaX6sBN90E7N8PPP64XfBvbEgB+NBD0uzy6U9LfVmxtibNL4uLwNmz8X23t4EHH5RC/p57gCuvjAvWO+6QwlmhTDmzs9E1m6ahW26J9kvr0WLWz4kT8jxskmEGABbsg44uyIFIOOuCD+i2ZZ85E/+ta6i6AFVsb8tz7d8PXHWV1NBtfPe7crszZ+JCHQAefRSYm5O27htuiDcOtZoU6kBUjvFx+b1W6y6TPojqsuevrgI339w9DhCCWV9K42fBzgwCIfaavBe2seeEaQduteJ25fn5+LaNhtvmrB9LX5T9W3dFXFiQ+9vs6GqbvXvt61qt6HzT00IQReVRdn19DKDVsp+rVpP7m9esr89iHzfrq9lkOzvTd8A29hHAtAMD7mCa2VngXe+Kfm9udtuNr7wSmJ6WtnR1jFtuAfbsibZZW5Pa7HveE993x474Ni+9ZC+z0uhnZ4Erroi0+s1Nab9vtaKehrK5u3oQp0/LOlDlNdersvQSLTo7C9x3nyxPqwXcey9r68zAwKaYQWZmBvjSl6Rtu9GQZo65ObeL3u7d9uPo9vXxcWlLX12VduYzZ4BXXolvb/M9371bClNlo7cJY0AOmB48KIW37sYIAD/8IfCTn8hr0K9R38Zkc1MKXgB45JFu80+WaFGeWIMZUFhjH3SUZqo+Z2elIFODfTpzc0CzKb83m5EAPXYsrvk/9hjw5JPA7bdLW7o+iApIYW8Ofr7yirRrz89Lu7bec2i3AdIm2dK9co4fB173umidGoBVzM7K49YdOki9Djz/PPDNb8aF+thYXPtnmBGCNfZBQx8sPXYs0ow3N4Hf/33gmmukMF5bA+6+Ww50XnZZNIB4771xjV4dT+f0aff5iYBvfavbK2Z9XR7ryBF53P37IwG9thYXuvW61PoPHQLOnQN++Uv/Na+udvcAajVppnnpJekyaZbxuuukKYVhRpEQQ3zeCw+e9og5WDo97R7AtC31ejwASAj34KM5gLp3bzTQ6VvUYOXCQjzwSA2Mjo0J0W5H12E7pi1gSm2vB0O5AqHyCiriHDFMxUDg4Clr7IOEOVj685+n239zU2rzTzwhTRSA9En3UatJ7fe557rt1wDwhjcAO3dGphkVzHPiROS+uLEhB2WvuCIyE6nrMI/ZbsdNJ6pHodwXZ2ZkkNLiogx4ctnzdXNPL+jjDqbrKMNUHBbsVcYMNtIHEhsNt+dJEkrovfSSPdDowAHgoovk98sui/K72Dh/HvjgB+M293PnuoXt2bPArbdGwlE3I+nn1c0n5qDu8eNSqH/qU9Ex2215HefPA295i1y/tZU9xa4t8pQFOzMg8OBpVbFNWqFnT7z66u5BTZOa4/Y2GtLGffKkff1FF8no0YceklqyS6gDUog+9VT8v927u1MKrK9HtvxTp+xa9j/8g9TEAXm9i4tx4bq4CHz1q/F9nnpKnm/XLjk2sLUltfjDh7MJ4pmZuBvlo48Cb3ubO188w1QIFuxVxaYxHjokfx89KjVpl+BWvP719v+vvloK7I2N6D/9WI8+GglY07Wx3ZZmFZWHpdkEXn01Wq+8bQ4flsJWRx3LJRiFAD7zGZlz5qMf7fa8OXsW+NGP4v/96EfSvHT6dNT72NzsjqzVCZm0Y3ZW1pNie1sO0uo5cRimorBgz0raWYh6mQmo0ZC2cKW9Ly5K80iSxn7uXLeb4Pi4FLzmzECXXx5ts70tc6wsL3eH1u/eDaysAA88IHsOU1Nx7fuqq+TnnXdK84iOOpZPk97eBu6/P37Miy+Orze3t3HypHsmptDp++bm7G6WGxs8RR5TbUJGWPNehsYrJu0sRK1Wd9h80vHn57s9V8z0t75lejoKyR8bk94q5vGVF4u5r1rnukYzvW+j4fZWMfddWJDXsbDg9+5RZfalMHAtekoFhVm2Vsvv+WKrl1qNPWWYvoBArxgW7FkwhYRNkAjhzsPi2t63f60mXQ+VoKvX48JbX4i6habrnOa16MLL5fZnE5Jmec1c7a7rUw1EvR6/FpXTRneftAlaW4NgO6dZNpUPRm94zOvVz61y1CwssCskUzos2MsgRGNfWrJr2D6NXWn3SiAuLbmTainf9HZbiB07utfrPuQ+IWtq3yETUyRp8z7BZ67Xr9nVGJlJwVQjZ6sTvWfiOrctaZp5TQsLUZmmp+OJwUJ7XgyTEyzYy8InwGyaujlLkG0fMwtjL2YIU5NutezaqY4vO2Qv1+/bR5VFmW/0utIbI1tZ1TkXFrobo6TegX4Mm2nM1nPx9QzS9LwYJiMs2KuAKSSmpuxCStfOfdGUocvkZLcQDzEblTUdnNmAKDu37T/XtHi2npAyBSVhNiLmdHxqXUiDqhomhikBFuz9wpynUzeDmAJgYSEeUt9syv9CNPZ63W6GcM1nGiq0XUI0T3uyTbDbtHhXuc1eTdqGKKmRUwO77bbbrq8Wc0CaYQokVLBz5GmemGHohw/L1x+Qn6dORYmxxseBr30tWg9IP+zVVRl9qba77DLpk/3889J98bXX5P+bm8CvfhU//9RUlIQLiLsVquAmV0pffTszpD/v0Pq5uWh6PT3LpMoAqWeCtM3V+tJLcR/8Cy8E3vnO8PPrEbxmhOrychRp+/3vJ7uUbm3JKQABWfcMUwVCpH/ey9Bp7KFuiS4brb4+ZFYjpdXqmnxRJoFQz5+0mL0A13lsGrsrcVmaOnD1QkKSooXcO4YpALAppiR87nMLC9G6pC49ILv+Oi57+9SU3ZxR9PUVYXdXglr3OEnysLH5lufR8JieQUTpBq6LugcM0yFUsOdiiiGivwTwuwB+KoSYyuOYA4NuKlhfl5M7XHqpDJ9fXZXmmDNnpBnl6afdMwsBwIMPynD13bulecA2e9D4eNTlf+KJyJygzzqUJ6EmnF5YXpbXq8wq9bqsv7m5bjOS/tuMhlU0GtkTf+lJ0YSQqRZ27eqOomWYCpOXjf1/AjgK4K9yOt7gYNprlYBVdulmUwqIjQ35fedOmW73wguBX/wibsPd2pLh9IC0sd98cyRUVUOhC9eiBK5JUVPEnTgRt5VvbspGUZ8AxHZtZpbLd71LNoZmg5AWW0O6sQH87Gf27Ws1uWxuys9z56IUBWXcF4ZxEaLWhywALgFwNmTboTLFCJFsL/bZZaen3bZ3M3py2CZ+ML1blOdQaOBXEXWhTEOqXC4T2t69UUSsfv/q9XgUrfKY0VMoMEyPoGwb+0gLdhNfiLzNJuzLhaKEV9o8M4OC6ccvRHEDtmnQXR7N9Ah6OoGkwVZbdOz09PDcP6ZUQgV7ae6ORHQjgBsBYM+ePWWdtn+o2YO2t6OUuLVa3PRSr0fd9f37Zff9H/9R5iUXQpp2JiYis44i6+xAVcJm5vG5I5aB7vL44ovx2ZuA6H7cdZe8Tzr1uvxPv//mnKynT8tj8KxMTEGUlrZXCPFFIcQ+IcS+XWae7mFDnx1oezv+XUEkbei6z/nMjBxoFULmO1cCxZzooh/Crkz0CUX6IfxM3/nVVZkDf3Y2vm5rK35P9+yRsQm33BLlq3ehGmeGKQDOx94vhOj27jCFxuqq1NiVkGg0pNfIKGh6s7ORME2b8z4rZq56vRGdmbHnaB8fB/78z2V5jxyRwl0PtFKo3ps6btnXxowEebk7HgdwJYCdRPQKgD8TQnw5j2MPJHNzUkj7XBuB+IxCyvNlfDwyQUxMSJOAmu7tpptGL7qxH5NK+1w8Z2dlT+v22+V9aTbl5CKXXRZp4LOzUc/L5Ior5KQmplmHJ8xmciQXwS6EOJjHcQYa0z3vPe+RtlQfn/0s8Nhjcsq39XUpzHV7rq7Bb266/beHGVtKgTKEn8/F88iRaExkZkamirjjDnmPlIB28eSTwJveJGfBajb7c23M0MO5Ymz4fKhd2yvN6+67o+nhklhfjwt/3Z6r6OcgYhnodQ1013uagdS09y0LSvAvL0faOxAJ6Lk5OXesbSo/FasASHONGigfxvvL9IcQ15m8l0q7OyZlZLRh81s33RcvusjvFudyYxw233UdM0Wu8gdPO2mHeays7qBp6ty89/oEJQcOJN9z09+dYTyAc8X0SEgOFvPFb7eTX949e+K/f+3X4mlqQyeIGCZ8gVxpfdfz8n1P20CYeWtUvp+lpfR5Zoa1AWdyI1Sws1dMEmfPxj0W9Fnuf+/3gA98QLq4+Wg2gT/4g8gjApBmmJtuki59990HPPTQ6NlXZ2Zk3ZiMjaU3S/g8WdJgs+n7OHMm/nv37ug45uD59DSwsNDtVdNsAidPymfq4EH2kGEyw4LdZG5OuhUqfvjD+MtmuiT+/d/bvR90pqbkYJuOyr2uXPpGFbPu6nXpKpi2TvLyfU/TQCwuAo88Ev1WAWfqODa3yNVV4CMfkc9Euy3Le9VVUc6ctTUZB8EukEwWQtT6vJdKm2KEsE+7pucHT8qr7upmu2yxo4pZJzt29N/WbEtx4NrOzCNDFN8nyUSnTD2mOUeZcIYpdQSTC2BTTA8sLwPXXis1ptnZSHNrNOSsPcvL8v/LL/cfR2X9A6IZgnRNcGwsHnU6aqignImJeO/otdekC2i/NFVlZnv4YZkS2ceJE5EnjEKIuOlGmWVcKFOP6caqTDisvTM9wu6OCjM3+OOPAx//uLShPv549LIfPw7ceivw0Y+6A5CuuSYKUpqYiFzwykqzW2V019DxcZlyV3f5XF+PhGPZdZXGZ35mBviLv4i7M9Zq3VGqd98dz/Gu02jI5+PMGakArK/L/4ii72oKQeVGmzU1MTMahKj1eS+VNMXYPDSmp4W4+GK7SSZkeraiZx8aJFzTB+opcgFp3tAzKpZZb714xChzjGtS66UlISYn7c9KvR5de70ut1tYiMxBtv1G/TkaccDujikxc4Pb7Oj6S+VyZ9PdI6uQfrYK+KYPVPnX9bz0pu26zHoLtbHr26fxsSeyC3n9GWu3k1M9MyMJC/Ze0F9qU1u6+OLkiY/NAVHW2CVmA2fz2bYNLvej3nz3LEuwmNpXnwdXb+R8y+Skez5YZqRgwZ4FpUHqL5fqIusvti4EzK64/iKPeuBJ6IxI+jZJ9VZURK6rl1VUZKs5iYqpqauo1DS9CCY/Khb5HSrYSW5bLvv27RMrKyulnzcIfXCvVpMDfO9+N3DRRdFAlp6sS81Fqs9JCsQHCDlrnzuPy/Ky9PwAZIZEsy5t9WYOwOZZv65jHzokA4gU8/PxnD55nFcNtqtJPsbGgOuukwP3ZnnKzIszqhT5nPUIEZ0WQuxL3DBE+ue9VFpjT5rqzLTBK+1S1+bMY7BN1I5tzlOzLm2aUp7pA2zamO1/m8aeh2nGdhz9t+1a2cSXD0n3r4JjZGBTTA/oCcDSLGYwU6vFL54Q6V8cV12GCNleypb2GKYw7rUM+r6NRtwMs2dPt0nPPE8FBc7AkdY86PJ6KhkW7L3gS0rl82KwaZkVs82Vjpm50RRWrZYcx9BtyqbGXqu5B1Gz1m9W4Zi0v698Ic+Zqi+blw5r7NkJvf+6S2totlednOUAC/ZQTC0sTbqAyUl795npfnHGxqI61ntF9boU8LrgWliwu/vlqZlmFY5J3jOuRs1cr2vs5rOl+/jn3bCNOmYaB5c2HpLt1UUBDTAL9iRMbwRV8aY3jGshqkTXrLLYGsn5ebu2agps2zZFaKYuW3qowHRt62rUXPvanjtXw8YCPR9C03Ob6bbTCPYCTGYs2G3YfIltL47Skoj8Gjx3g90sLcWFU60WuYzqGqoeqavvq/t6Z3XzCxWGeWlY5rXbhLKtTAsLcozBply4TH5MbyQJdv1ZMJ/VhQXZo5qe9t+DtC68AbBgNzG7xz4hbXsBXdGAPHBlx6d1214Mn2dIFtII6zw1LN02awpl3fzimjXLbNiU0OdnL4yk50dXMGy2c/NZmJqKnl2zZ5Uk3E1lMkOjzILdxNY9Vi+ZTSPUb4gZTMJRgMnYNB6lGdkGmovSRNMI69BypOkBuNwWk8rk8sDRBYo5NpGmbMNMmvtoC/zyvedm4xrawOakNLBgNzFvdrstb5IrcZP5Irk0esaN7QWx+fgX6b6XttEI0fR6aYRcDV2Ip0WSF40yERTZQFYdXyOqtG3bPj4lw2YGNDV2IGysLad7M7qC3fdi+rpFIdoVd33DcdmTk16msgZJe8XXCIU0CjY//SR0s45P6x81/3bX+2wbQ7M9V7b6CqlDcwxEeTBlVRoCGE3BHiIgXC9XLw8GY0e/D7bADtsDnvah71evyfWMpen+p+1BqO1dmSEbjSifjM92r4436L1N1RN05ddR12cm8jOFtO5Wa1MyGg23ucvW+yohiGk0BXtI0IhNYNvMA/pNVznCB/2FKANbw5lmGsAQwbO0FJl3bF41RWMrY1pbfuiz5DPD7NgRDazqwkgJItOGXGbPqKgGxCdUTQFtpqvQy9JuRw2lEsh6D6DV6g6eM4W7zd5uc23NkdEU7OaDa2bFc9nelpaim6iEkGuwlTV3N66XztXI9poqIEvQSFFksb2H2vRNzVR5FbnyyZhupUXWm6npJvUcesXlbWW6Eprb7dkTrbPZyc1euy1nlO0ZtqVdLtAENpqCXYh4q+tLMKUPjNii0Hwv1LDbLnvFfJmURuSLzNTXhWq9VRTsQvRmTvI1BvqzrI6rXEV1RcQ245RNAOrmi7x7O2m9fnrFbEBctm1X79ymaRMlz+zlGuReWopPElOw4je8gj305XFFOKruqe6psWtXfDs1qYZvsJXpxtZjConMVC9+Gju1zwe5TLKYHMx60N8LX12YQkiZEsw4AJsQD8k82ss16QO8RWrsacpnE+Lqfdb/a7ft9W2asvIoU0aGU7Cn6e6a3STdzmu+UL5p8GxaE+PG9YDrL4mvoUzz0vb7fmS19esmQL23KIS/8bPllnEJ6KQEYspcqQYJ9XuTNBhoU3xUbpwq3B8h7IFiKuDQdHeuSpk9DKdgT+vOtbAQCW3XJNM+NzKbt0zFb3xf8Ql1U3vMqy77+TL2YhIyy2vmItF7izazoa0nmrZR0QWy2UiYSo5r4Ftv1Gw5gaqAaU51TZIeMs5Rkfd+OAV7qE1S/e/rdqptDxywC/V63d7lZXOMHd+9cZnFbMcIfYGUNpp3lz9NGdIKdlsd2QbyTFOAbja0+bTb3qeQ63DZxU23Stu9Mq9dCXdbNst+YZbRZpYJkSkVMsMOp2AXwq8Vmtp1SLIpm8Zeq8nWPclXlolICtxJk/grqXdkGxjL4564XmLzmbNpvL1EkKryLixITd22ztxHmU10Ldnmihc6VmE+181mPC+Sa4DSFJp792ZTfPLWik1TrGt2riQrQMWCvoZXsLuwuTImaVOuafDe8IbooVCaOw+g+gnRfHwDUWncS10ub1nviU2Imr0CW5SjTegnNUihnkJmPIX67jN/hAojXfuv1bqTsvl6RHpjrUxFvQrAIrRi271U5zLvF2vsjoMA1wB4HsAPAHwiaftCBLtN40u6YbZsjYAQF15ofyhUhr2qdDWrRpJQ69Vf2xdo5nN5S1tWW+Y+89nwZVgMEQIhdaQaC1d6abN8oY2EeS6zt5pWW9WvJYsAzFsrdjVKvt6+738zHqaPlCbYAYwBeAHAWwE0ADwN4Ld8+xQaoBSqPfkmrb7ggm7BXrGWe6BIYxowp8xzDQymudcmtlByIfy+2CEKQx4CSq8r38C+z1yV1DsSovv5tw2Spn3mezWnhGjNacZebI1+Um/J9Hwzj+PruZQ4sFqmYH8/gMe0358E8EnfPpWYQckn2PVF2S8rZmsbKELrzmY7D7FdpxFAtsFKmx+9/jKbvQKfhpe18XeZpHQB7BPYZs/Hlb1UH/NQk6C4jleG4MpSp/q+rmctJHZCbzR98iFt3EWOlCnY2wDu1n7/IYCjvn0qIdjNwRXbsmNH1HpzHvbe0V8A34xILo3Z1xDMz4fbd5eWkt35TPNCWtNSVkFoCgszk2CS542tYTDL4guMqhoh5qAQE6xLCLueOd3V2aWx90HZq5xgB3AjgBUAK3v27Cm8AoJYWvLPcaq8Y0KEEuNH1bUv9NqlPblMDbpGrXsvhU5MDPQ2XlK0ppbFdr20FDZGERIubytL2f7cSddvE65pGl3bM6eei6RrH3KNvfqmmKQHst2W7mYHDrhTo5bYKg8lLoFjs5XrkZAhwl9pnbYMf764hrIC5bISYjPXMaMtbfuYCo0rrkBXavrVY03qOWUVrktLvU87OMQ29jqAFwFcqg2e/rZvn1IFu62r5nJ3cg1WubwPmHBsXV5bJGCSjdjl6mgLFjIH0dL6nbsoW1Pr5XxJwjBEY/cNJldJwclDuOZ1Twvu4ZTt7tgC8P2Od8xi0valCnbz4TSFtGu9LjQ4lUBvmA+5OaZh88U270et5u7+ml4PPr9q815mvadlampJwV9pBYnLx9ukKhp7KHmMb/Tqrqu20eNfCkiENnoBSi5cXXfdFGDT6FmYZ8OmAdmErc18Ygp8U/C47J16Qi6XD7jtnhct3PMeULX1NtMIXnM/X+pbc6KOqr4TRfaiQo8d4kmTERbsOrrfsq69qYe1IsEHQ4VrUEt3s1OmFlPj3Ls3vu/kZHotVB3TFqiSt795km9+UY1IiKnEN5Do8/TKUu5+NABZ76mvnkLt7yXklmLBrrDdGDVzkmvOw6pqJYOEy2QSknGz0XA3xK5z2UL/XSaKPIRtqCApcqA1SWNPuk5f2Xotd9njD3mc17WvLeFa0nOomwML6PWzYBfCboZRN2ZhodsDxkzTy8I9GzaB6/I5t/lWh2hKtgYkJK9PUeaRXrdLOleI7dfcLq0PuH6MXnMj9TOQr9d76urtmSbByclkYV2wYsiCXYjuG6bPcWq2xGNj2RIZMXadbfp4AAAYE0lEQVRsA3Uue7EtZ3uSLdj2UpYlXHzdd7OH0OvLXoQm6iqruX0vGme/NPYs2MpsM3NNT/f92liwC+F+yGw3TYVVD9pDWXV8WqFZv7aGtRdbsE+42M6dp5aVt2AL9WDxlSf02vJqEKtqzgzt+ajfujnQNj9DHxQ/FuwKJRhc04OZL01VH8pBRr8HoXnWfQ2xzaQQIqxtxw8J5klD3r2FpaW4m6iZXjdPbBr7sDgWpDGdmYno1BSBruRxJcKCXeG6oUtL3dOSVWW2+2HDN9Zh27bX0O2kRtmm/ZqDtMpcl8e1hrz8IYqEzY1OBRXlrYjoNnZfENOgKUAus6yO2Yi6zDBJc8EWCAt2RVJwR5bJiJkwXG55abTZJEESak/Wt3H5HWcNSgsVekllVsdpt+3lLHKw33bPdBfKQTNZhigXPvfRiszryoJdkfQQuibPGDSNpMqk0dh7JU2GR9tgobmU0eVOUjp0DdGmsfcymXYopvaqa+z99HzJgs0TTi+7mdLZ3LYCqUVYsOv4vBeSIvoGRSOpOnoXP+8G0yeEkrD5KifNLJRnuV3PmSk8dXuvsrEvLSXPJZu1fDYbuykAB2VGMbNOzXQVNpfbrF5COcOCPQRT49mzJ3qQfS82a/P9wVXvWTxHbLZXm3dUUffc12MMESppPTXyuI5B1djNhlDPxa8asdDp9PoEC/YkfAMl9XoxodajTshLEtK7MgevstwTX6/N5d+d50BliI29F0+iLIPQtjIMyzvgcqk1g9wqek0s2JNImhqv1QrTDgdFW+k3IULMpjEpbKaJvLw00nrT5HXPQz010kbRuvZJex2+c/dDk83jnLZrGqB3mgW7D7NLZi4uG22S8GHc+MwlrkFMc2DLN/BVJEVpqLb5V00TQC9BMS5BlfY6qiTw8rwHoQ1hBWHB7sPm1tRsyqyCKh+EyQB11wojq1bsGuhzuZm129G+RXvVhJQ/bw3Vdd26G2OzGQ2chg6O5qVpV0ngFd3I2IR9hWzrChbsPvQHVrmNJaUTqJL20g+y2rHn57unY5uejmzINo19ctKutWYNIqoKrgbLrCflpZPG2ycvwVQVAVdmI1OlBs0gVLDXMYrMzgLHjwMnTgAzM/L3oUPA2ppcv7Ym183ORvvMzAD33CPXjY/L38PM8nK8fk6c8NePud/EBLC6Kj/vvFPuUzcet9On5dJsAlNTwM9/DrzwQrT+1VeBL3xBrm80gI0NWfdHjtjPPWjoz6FeX5/9bLRNrQZsbcnv6+vuei+yjFWoa9s7WxSuZ918J6pMiPTPe+m7xi6E387mmum+KtpL0dg0lpDBT31CaX2Q0zdIrS+NhjS/TE11a63KDdXM+aPuh8ttsKq4nqU8fKkrrHEOBL08/yUBNsV4cN2kvBNCDSq+wTeXK6LLdVQ3JYQsalDVvEemqazdjhoNs/GounBPsoG7fOhDs4+OutkwhCQlzVxfkToNFey1fvcY+oKtqwXIrrDq9ur/jxozM9L0AchPZXaanQWOHpXfDx2SXVNA1tPGhvt4110nTS06e/YACwvA5KR9H9X1np+Xn6ur8Xt2//3A5qb8rT4VqlxVxfX8Ad3XrUwhR49218GxY/bjz8xIkxUwGmZDk+Xl+PNpW3/woDTzHTxo307VuTK5DFqdhkj/vJe+a+yukOiKdLf6htJS2u3uATuf1mhz2zO1G1fd6t4yvsFBfTuzB2AmaBpkjT1pv9AUAqNiNjQJqdtete8K1CnYFOMhKflSSO7wYcPlS67MI66EVPPzyXbhJBeykBdGN/fokcF6ruxhsLEnUYHJHipNaNBXiDKhtq2QDGDB7iPUxqmWUdDeXT7VY2PS7dC2Ttfmi86tYgtwGpZJINIw6r3KJELe36Wl+PiMr9dTsbpmwZ5EqFfCqGhGtiyHtZrfq8WMHi1Ss/ENpvZx4oO+UDEtsjR8162vW1ryT4Qemu64IgOmOizYe2UUNXbT1bPdtgcGTU66k6OVVU5dU9fL5tO8mMGnF08i1yQ6tufa1VgMqMY+mgFKPmxBI4MQkJAF3UtjcxPYvVt6BCwvA088EQVlff7z0fZl18nyMnDbbcCZM9JzqdmUwTvb21G5yw7eYcrD9CS67bboObR5GU1MRB5uRPFjzc0Bjz8uA74A4OWXpXeM8kJSlBkUlTch0j/vpdIa+yiSpA31u9vv8pPfuzdKDFYRjYopCNMjSO9N62Y51eNMmsrOZq4ZgDmPEaixk9y2XPbt2ydWVlZKPy/joSrh0rZyHDokfY51ajXggguk1lWvAzffLFMNMMPLtdcCDz/c/f/8vNTQb79daul6Tw4AxsaABx7ofq6Xl4Ebbog090YDuOmmSvfSiei0EGJf4oYh0j/vhTV2xorP113X2Gu17pQDFRjYYgrGN/7lcnpQ7rcubGM1Fe4Bgm3szMDhisg8cUJqUmfOyN9zc/Lz4MHupGxV6Xkw+aPbvM+dA557Lp6kTCXpazSkxr65Kcdibr3Vfcy5uWgcaWwsimL2JbobAFiwM9XBzKA5MREX3ocPy24yYB/YUqHia2vyOOZgGDP4qPup7vOLLwL798v/rrxSfs7NAadOyechKTul6SyhMpEOQtoADyzYmepgCmtTg7/jDqlR6UJbf2nTphZmrX6wUPftpZe6c+bo3luXXRYJaCX4k4S7Wr9//3A8GyH2mrwXtrEzQeg2VVsqA9/2vvlBK+abzASg37dmMx5PYZuIZUjHX8A2dmbgSdtNDvE7DtXqmf5h61Hp9219HWi1gEsvjZ4DXWOfnZWa+hCYVHomRPq7FgDXA/gugG0A+0L3Y42d6Yk8fOrTaOxV8OEfNXyeUb6JcMx7NaT3DmWkFADwTgBvB/AEC/aKMaQPdi6EZpNkk035+PKz8EQ4wYI900QbQohnhRDPZzkGUwAhEwkME0kTK7j4678G3vUuYHGx+zi+yTCY4tAntKjVgCefjO5rmolwen0mhoUQ6Z+0gDX2alHBrHSFkVazduWdb7fjObpDp6Fj8kHvRS0sxFMC1Gryv9B7nWdvq2I9X+SlsRPRSSI6a1k+nKYBIaIbiWiFiFbOnz/fUyPEBDJo03hlwaVZLy/LEPRrr41rbfr2Oo8+GoWWr6/LYChzijqmGMwe5mOPxVMCbG8Dn/mM/G67J6Z2nldva5B7viHSP2kBa+zVo2KaRmEkpWxVGrhtEM5MKKb/3rNn+Ouu36hn1HRXdE1+nuTi2mhEM5/lobFXsOeLMvOxs2Bn+orZiNnyhpjTH+qCXHX1zQySozaBR5mYMQoqR4tLqLvmd7Xda5XxsUwPqpIIFeyZBk+J6DoiegXA+wE8RESPZTkew+TCzIzMF6JoNqUfvDLN3H478MIL0frtbWl6ue8+YGoq+n9rS25rmnMUoz5A1wuqzo4di8wlW1vyHrRawC23RGbERgOYnpb/33ef3Rymmx0Va2tyoPXo0WwmNBUXMYjmuBDpn/fCGjuTGz6/Z31Sct8Uf7rpZWEhyvFuyyKYdF7GjRk96sqZHjoFnv5fqxUNfg/x/QBHnjIjgSuSVM//ce21UdY+Fz/8IXDdddLFTpoX45hRqhzBmowZQWpGj05PRzNi6YP8rsRdtiRvgDzu3JxchiHPSw6wYGcGGzMjpC19byjb23FvDEBO4rG52e1dNDMDfOlLwMaGNBkMs+dRL+hC+EtfAq6+WtahXp8qna5PGOv30WxM9eRfStAfPVraJVYZFuzMYBOSvvfwYeDkSSmE0/Ke9wBXXGEXPGouTXNOTR/DlFnSvBaXEN7YiM98NDYm74m6fpdAP3Ysum/qPo6PR404wL0mFyH2mrwXtrEzhWJzU7PNcanssbY5Mk1vjCTPmxBXuKrb5UNdZJVNW8+waLoY6r/TuC7Oz7v3bbXiZax6fRYAynR3TLuwYGcKJSSRlPJ5VutMX2p9mZ7unjJNj47U/eR9VNAv+l/x1ZmZXMsmdG2pcs1BTXMgWq2fnpaL2s410F2vd9ez3hiMQNxGqGBnUwwzfJjmGUC62M3M2NP6KnfFZlMO6jWbUpRsbMjB1NOno2OvrQG33QacPRvZ43/1K2k2UOd24RoPqAKuaE1zsNIWuTs+DrzjHbJOFBMT0SCoMtFMTETTG546FU2cYrK5Kc01W1vSzCVE9P+xY/E6NmdU4pmzJCHSP++FNXbGShHRsr7uuqlR6lq8ma9EX2zukKHmgKpGBNvqyezFqHLrLouqvsxtWy3/OVyBSGaAkTlpue24Ve4J5QzKCFBimNwoKi+HL5fMwYNyUE/liNnYkJM3zM525yvREUJq8ibKUyMNVQlyMoNxAODxx6P1uufPlVfKoKF7741cDM+dSz6Hfi9Ulkadel0e9/jxaB7TD30oCjZrNKKJzHVGKTdSKCHSP++FNXami6K0LpfG7gpFV7Zfl7auFjO3jG4H1m33vrJkyWliDiKG5JdP01Mw60ddk6v8pl3cNhGGPuBaq8XreOfOaB8910+jEWY/r2pPKGfAg6fMQFGkh4MrWtE1kGoKtVpNCnJlglHb+wS/7RrM4/Y6N6dZ9qSIyyRzlE0g2ia1MMs/Oem+/qkp+4BrvR4fdJ6ejnvXhJp1RhQW7MzgUbbW5TqfLoiUUFfCh0i6R4Zo9bpwE6JbWKbV2F3ZEJMaB1uD4nMX1FMw6NPQLS11e7iENG623pGvcWPB7oQFO8MkkZSTxOaqp4SdGmBVgrpelw2AKexVhkhTa9WFZajvuK6l21wCQzR2fVvX4Kg5sKk3FuZgprns2GE/pj7g6vN/Vw2O3rAMuXklDaGCnd0dmdHElnfEdKM7cSIaWNXZ3IxcJj/0oWjg8JlnugdcVYbIN785Gjjc3JTZB9NEoZqRnJOTwMsvR+unpoAjR+zHUQOji4uRS6I6lh7JqSJG9YHNej0+GLl7t7+cb3878Oyz8piNBvDSS/J/5SqpuzzOzcmy7d/fHcE6Nha5PTLpCZH+eS+ssTN9x2We0DG1bJcZJe0yNibNOT4zjC8wqNeBV1NzVoFB+viC6ZJo9izabfd1qYFOFXRk2v5Dx1FGyH0xLWBTDMN4cJknfAJW/+6zG4cspi+8Lrza7Wh9UhRokhlnYUE2WrqAbrXiphwzctbXqNjGFdrt7lQAtpS8oQJ7BFMFhBIq2NkUw4wmLvOEmUjKlkJWmRRUpKqiXpciy+ajbSJE9L3ZjMwdi4vA/fdH63zJrfQUwvpvxeIi8KlPye/qGo8ckdvrEZ/r6/Fz6Ne8vCyPo0w3Nt/+tTWZVfHQIft2ug98SOStLbEbkwoW7MzoYoajm8LGlr1QbdtsRgJc2aHPnYunHwhFF/JmoBKRPLYrF7lvnMA81vKyFOwTEzLASglfNcPUoUPyc3U1qgd1fEWjIVMo6GVW6CkTdK6+OipXqMB25WRnwghR6/Ne2BTDVIokP/cQtz09+Ma3mIE5+v7z80IcOBD//8ABWR6bD/mePd3H0GeNMoOoTO8cQG5jGy+wec2ocQg93YIZkGV6E7EpJVfANnaGyYAr9a8r0jLJDVAXjGYUpj4wa9rep6f9qWx1O7lv+r92235dSjDb9pmeTs6zYwYX6etHIBK0bFiwM0wWXAN4riChVisu9CcnwzxfWq3kQCe1vxoItTUYSZGwKsjH5qeuR3/qS1I4P3uvlA4LdobJSlIAky2wxmbS8Wmuabxr1HH0oCkVwKPnV/EJdiHigVXqGErbd00w7aof9l4pFRbsDFM0SUI7NDmXL4WtrrHrIf9qH91VUfUAWq1uW72ZlMun4euulkkJuNjkUios2Bmmn6TRZnW/dXOxZYsMMYG4tlGCWLfZ64nE1BISRMWUTqhg53zsDFMEZh54V5725WWZE14IGT5fNzyQNzejHPGKkPzjtm30nPd33iknh56fB+67T+ZWn5qK9t/aAp56yp7Lnqk87MfOMEUwMwPcfXcUwHTypBSstmAnffKJVivafmNDCvuJifg+IQE8tm30AKK1NemvfvRofL/rr5fnBYBXX42CsHgCi8EiRK3Pe2FTDDMS2LInmoSkzs3LDBJiHrJ5+7ANvTKAUwowTJ+ZmwOeeMIfQu/SvldXo7B/X1qBNIRo+maZVQZGZqAg2QiUy759+8TKykrp52WY0kmTmtfcT091YEsXUFQulSKPzWSCiE4LIfYlbseCnWEqikvAJgl9ZmgJFexsimGYquJKhGV63ORhpmGGCnZ3ZJhBI8TdkRlpWGNnmEGD85UzCbBgZ5hBhPOVMx7YFMMwDDNkZBLsRPTfieg5IvoOET1IRBflVTCGYRimN7Jq7F8HMCWEeDeA7wP4ZPYiMQzDMFnIJNiFECeEEGpW3CcBvCV7kRiGYZgs5Glj/y8AHnGtJKIbiWiFiFbOnz+f42kZhmEYnUSvGCI6CeANllWLQoilzjaLADYBfNV1HCHEFwF8EZCRpz2VlmEYhkkkc0oBIvoYgDkAVwkh1gL3OQ/g5R5OtxPAv/SwX9FUtVxAdctW1XIB1S1bVcsFVLdsVS0X0FvZJoUQu5I2yiTYiegaAJ8D8B+EEIXbV4hoJSRPQtlUtVxAdctW1XIB1S1bVcsFVLdsVS0XUGzZstrYjwJ4HYCvE9EZIrorhzIxDMMwGcgUeSqEeFteBWEYhmHyYdAiT7/Y7wI4qGq5gOqWrarlAqpbtqqWC6hu2apaLqDAsvUlHzvDMAxTHIOmsTMMwzAJVFawh+ahIaJriOh5IvoBEX2ipLJdT0TfJaJtInKOahPRPxPRM52B5VKmjEpRtlLrjYguJqKvE9E/dT5/3bHdVqe+zhDRcoHl8V4/ETWJ6G86608R0SVFlaWHsn2MiM5r9fRfSyrXXxLRT4norGM9EdHnO+X+DhG9tyLlupKIVrX6urWMcnXO/RtE9A0i+l7nvfxTyzb511vIjNf9WADMAKh3vt8O4HbLNmMAXgDwVgANAE8D+K0SyvZOAG8H8ASAfZ7t/hnAzpLrLbFs/ag3AHcA+ETn+yds97Oz7rUS6ijx+gH8NwB3db7/JwB/U9L9CynbxwAcLfO56pz33wN4L4CzjvUtyOhzAnAFgFMVKdeVAP627PrqnPuNAN7b+f46yJxa5v3Mvd4qq7GLsDw07wPwAyHEi0KIDQD/C8CHSyjbs0KI54s+Ty8Elq0f9fZhAF/pfP8KgI8UfD4fIdevl/d+AFcREVWkbH1BCPFNAD/zbPJhAH8lJE8CuIiI3liBcvUNIcSPhRDf7nz/JYBnAbzZ2Cz3equsYDdw5aF5M4Afab9fQXel9RMB4AQRnSaiG/tdGI1+1NtuIcSPO99/AmC3Y7sLOzmFniSiooR/yPX/6zYdBWMVwL8rqDxpywYAH+102+8not8ooVwhVPl9fD8RPU1EjxDRb/ejAB1z3uUAThmrcq+3vs6glFcemn6VLYAPCiFeJaLXQwZxPdfRLqpQttzxlUv/IYQQRORyx5rs1NlbAfwdET0jhHgh77IOOP8HwHEhxDoRzUH2LH6nz2WqMt+GfK5eI6IWgP8N4DfLLAAR7QDwNQCHhRC/KPp8fRXsQoirfes7eWh+FzIPjU0QvApA11be0vmv8LIFHuPVzudPiehByG52ZsGeQ9kKqTdfuYjoHBG9UQjx404386eOY6g6e5GInoDUcPIW7CHXr7Z5hYjqACYA/L+cy9FT2YQQejnuhhy/qAKFvY9Z0AWpEOJhIvofRLRTCFFKDhkiugBSqH9VCPGAZZPc662yphiSeWhuBjAr3MnFngLwm0R0KRE1IAe5CvOkSAMR/Vsiep36DjkYbB217wP9qLdlAH/U+f5HALp6FkT060TU7HzfCeADAL5XQFlCrl8vbxvA3zmUi9LLZthfZyHttlVgGcB/7nh5XAFgVTO/9Q0ieoMaHyGi90HKvTIaaXTO+2UAzwohPufYLP9668dIceBo8g8g7U5nOovyUHgTgIeNEeXvQ2p1iyWV7TpIO9g6gHMAHjPLBunV8HRn+W6VytaPeoO0Tz8O4J8AnARwcef/fQDu7nw/AOCZTp09A+CPCyxP1/UDuA1SkQCACwHc13kO/y+At5Zx/wLL9unOM/U0gG8AeEdJ5ToO4McAftV5xv4YwJ8A+JPOegLwhU65n4HHY6zkch3S6utJAAdKvJcfhBxr+44my1pF1xtHnjIMwwwZlTXFMAzDML3Bgp1hGGbIYMHOMAwzZLBgZxiGGTJYsDMMwwwZLNgZhmGGDBbsDMMwQwYLdoZhmCHj/wM/G/YHeUkjkgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(y[:, 0], y[:, 1], s=10, color='red')\n",
    "xlim, ylim = [-3, 3], [-3, 3]\n",
    "# plt.xlim(xlim)\n",
    "# plt.ylim(ylim)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Animation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "from matplotlib import animation, rc\n",
    "from IPython.display import HTML, Image"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "# !brew install imagemagick"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "x = sample_n01(1000)\n",
    "values = [x]\n",
    "for p, config in zip(ps, cs):\n",
    "  shift_log_scale_fn, flip = config\n",
    "  x = nvp_forward(p, shift_log_scale_fn, x, flip=flip)\n",
    "  values.append(x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJztnX+Irsd137+zu951ruTeICQcN/bK6k2w4yrBzl5CKkLqYnkr1u3aCtehorSkKUiFq9aiFDt6BQ1pMVQxBEMlqEpj0YC5QcEN7xIr+NqJS5s/rPgqsV35J06MiU2wnRbfVCzo9krTP54d7+x5z8ycmWeen+/5wMvuu/v8mJlnnu+cOXNmxlhroSiKosyHjaEToCiKotRFhV1RFGVmqLAriqLMDBV2RVGUmaHCriiKMjNU2BVFUWZGa2E3xrzaGPPHxpjPG2O+aIz5tRoJUxRFUcowbePYjTEGwC3W2heNMa8C8EcA3met/UyNBCqKoih5bLW9gG1ahhdPvr7q5KOznhRFUQaitbADgDFmE8DzAH4MwJPW2ueYYx4E8CAA3HLLLXtvfvOba9xaURRlbXj++ef/ylp7R+q41q6YMxcz5ocB/C6Af2mtfSF03MWLF+21a9eq3VdRFGUdMMY8b629mDqualSMtfb7AD4N4L6a11WUQTg6Ah5+uPmpKBOiRlTMHSeWOowxPwTgnQC+0va6ijIoR0fAAw8ATz7Z/FRxVyZEDYv9dQA+bYz5AoDPAviktfb3KlxXUYbj6lXg+Lj5/fi4+a4oE6FGVMwXALytQloUZTzs7wNPP92I+rlzzXdFmQhVomIUZXYcHgJXrjSW+v5+811RJoIKu6KEODxUQVcmia4VoyiKMjNU2BVFUWaGCruiKMrMUGFXFEWZGSrsiqIoM0OFXVEUZWaosCtKTXR9GWUEqLArSi10fRllJKiwK0otdH0ZZSSosLdFu96KY3+/WVcG0PVllEHRJQXa4Lrex8fNglFXrugU9HVG15fhOTrSMukZtdjboF1vhXJ4CDzxxKmArXuPTscdBkGFvQ3a9VZi+KL23vcC73pXP8I2psZEjZ9BUGFvg+t6X76sbhhlFV/UbtwAnn22e6t1bBayGj+DoMLeFtr1VvIYyrrs476+qDm6tlrHZiGr8TMM1treP3t7e1ZR7HJp7blz1gLNz+VyfvddLq09OLB2Z6e/+w1RpkovALhmBRqrFrsyHENZl33e9/AQ+PjHgWee6cdqVQtZgYY7KkMy1L6iQ9y3z92YdOentUeFXalDSazyUHHfGm+uzBzTuG365eLFi/batWu931fpCH+i1rlz6gJQVhnbJKWxpUeIMeZ5a+3F1HHqY1faM7ZIDGVcjC0Ec2zp6QAVdqU9Gqvcnhrhl2OamOQztoZ/bOnpABV2pT0aiSEjJLxSCzIm3EdHzexWN8t1TOI+toZ/bOnpAklMZO2PxrErs2K5tPby5XjMeCy+/PLl5u/uc/ly3vnWNrHy/jUODurkrRaSMuqTsaVHCDSOXVF6QGptx7r/Egty6u6Dsc3QHlt6KtNa2I0xbzDGfNoY8yVjzBeNMe+rkTBF6YTafmip4IbE20VnPPJI3JWVEv+HHgJ2dprfd3aa78r6IjHrYx8ArwPw0ye/vwbA1wC8JXaOumJGwkS7oyypvHQ1td93kWxvN/eIpcFPo3/u1pa1i0X6Xqk8zuV5KiwQumKq+88BLAG8M3aMCvsImNOaIqm8+P9P+bJL71/SaFDf+sbGtJ+D0jlSYa/qYzfGvBHA2wA8V/O6SgcM7bOt6RJJ5cX/v6NmNMThIXDXXcBLL4XTwLG/D2x4r+ArrwBPPVUnTXNmrGGdI6KasBtjbgXwMQCPWGv/mvn/g8aYa8aYa9/73vdq3VYpZciQr9oTRFJ58f+/vQ0cHPC+7DaCkUoDd+3DQ+Btb8u/1zqzBpOLqiAx61MfAK8C8AkA/1pyvLpiRsJQPllJeF8ubf3PIXdOThmFjo25ipbLUxfOzk43z8JP19T98CV1Z+p59kBfPnYABsBvAfiw9BwV9h6oIUhd3MsdPzb/PicYtdKZEqMuhdfPw85OM8A7pnLPJfeZjLGutaBPYf85ABbAFwB87uRzEDtHhb1jcipz24ovPZ+LCBmTFcXlo0bPImdgtQsRonnoYvC4b3LqThe9wwHpTdhLPirsHZNTmUsrvnu56IzHklmTXVHSkwiFI7Zt9FKhkNZ256Kak8WeS+26N7BBosK+znRtsVPBSlmjQ1hNXL5KXso2LzLN98FBmZ+/LZyrZ7EYV4+pS2qJ8QjcOirs604bH3vq3CEEK/fl5NLY90uZ2wC6c7oW3JoCNTaXWpeMwK2jwr7ulL5w/ku/ucnPhiy18kvTkzv5x53jux0kLqMuyHFZ9UUtgRqBBdsrI8ivCvs606YC0pd+a4s/v2/LUipE1Ke8t9eI6mIx7EtZ+ky6KOe+on3miPrYVdjF1K4sbV64xSJPSLuEi+hICRHXMLnzhvYrx54z978uLcQadW4EFuy6ocI+Fbp4OdpYh9RC7vKFTYlLblQJPWdzcxwNVIrQ85qCRbxOPvYRoMI+Fbp6eUteOJqWu+9OR9SUvtTSxifXyvX/PrT7hUtTTtSQWsQKQYV9Kozp5e06TNKnbYPWtmHo09JMpTX2f7WIFQ8V9ikxppdXmpa+hLmL+/fdmEq3vpOU+5jqSi5jSPsY0tACFXalW4aITafnunDG7e08Qezbd82VFc27pCz6bJBqC+AYeqZjSENLVNiV7hnS+lkuTyNeQiGZ/rH+Cz2E790vK9ooSdPTV4PUhQCOYSB4DGloiVTYdTPrudHnJgRDbgj81FPAzZvN7zdvxjeooBtxXL/erMce22O0Nn5ZPfUUcONG8/cbN4CPfrTdvqm16WITliHX/x9TGnpChb1r+hTaOW1CULPcuBd6TLvU3367THAOD8sbpJzy5Mqr7fNok/ZajCENfSEx62t/1sYVM8ZBuikgKbc2PvauyBkApZtrdJm+knpIXUcT9033Qg91DOpjHwFjGKSbItJ4+sWi+R+3nk0JfcTl596rhli0rYdzMRi6pKd3T4V9DAwhtLXFbgj8cgvNgK1dtovF6UzV2uGXpeIcymPu9dqW1VwMhi7pqfFTYR8LY5oIMyWWy6aBCr0sbePY/WeyWFhrTPx6qecYE+FaC7K5+/exmFhJOOY6oxb7mgl7n0yxyxwTjFj8d2nIIhf66MIm3Wdzs6x3wOWl9kQq6fW6ciupwIdRH7sKeydwojX0S5gr3LHza+SP89/73wFrb7/9rCurr1muXFlxlrO73tZWk87YMTXdSnPqEU4UFfZ1pa1FW+PeUoHJFUx6/O5u+wgPv5w2Ns5e/9KldD5SDZd7FjVnlfrjAdzOTKVLGKTqTs2lEZQiVNjnRu4LU2pl1hzoS6UhJWZcQ+FCHN0nFeoYc+dQ//Hu7qr1HvMxx3oQudZtzvOix+b64UNlkuoNlVxXqYoK+5woeWH6OsdROtAXakh837d/Lt1iLlcEDw7Cx+ZuMkKv7UfV5G6FJ2mAuGNDe6nGGmjuWYWeX8o1lLquUhUV9jnRl/Wdcx+p26Uk/HK5DG+S4bshOIs9ZeW7CUEhLl06PVbiE/f93bQBaTspiE5iih3bNvwx5qLqwzBQRKiwz4m+XhjpfULWdK1BPNrAuEW+/OttbDT7mYZcBf6OSzHrOWWVhtw29DsnhG38zTTNsZ5GCal89+3KU0SosE+N1AvR1wsjSYd0y7k24uBEcnPz1NpPXY/zPcesz1yfccqKjZVdqoGgdCHsOXUoVjZd1MXa15xpA6PCPiVCQjZGQtY0x2JxGmkiWc/FJ2RRSoWYiiJ3vdyGIjZhKiXquW4O34VUYy2ckt5TyTMoofY1+3YJ9diIqLBPiRyxHBLn1nC+31gj5PuIa+Yp5bN3afR93iG/ei2LnTvOd+XQBiHWQND7SwRDImTcIHKJGHUxQFr7mn0O4vbciKiwT4kc90bX6ZBYnb7/OgTnFmmbJ8lL5PKwt1dHPOn/uYYlFB3jR6z4rqFLl9qtS0ORxpe7stvZOe0N5N6/a4tdUreGSGOIniOBehV2AB8B8F0AL0iOV2FnCA1I9kXqZcitwNRiz3XFcKRC8mqJVyxPKYudRsdQS923+t3M0dw0lLpH3Lm5oZipNNRwRdDeYNtn1ueY1FwtdgA/D+CnVdhbMuSAT9vJRBzuZW1rgYXSQF0j3IDj5cuNhVxjxctYL8A9Oz9NklmhUpeIK8tQY1VrYDSXmteaahz8nH3sAN6owi6glnXTRVRC7mSiPip07J6cSNI80ElHMXH3xZnmi16HjhnQnkMsXDK3V+GfUyp8Jc9OckxNMe570HOCjE7YATwI4BqAa7u7u50XwCjhKm6uOHZZ+Yey+krvISlPbuAydS9fvF1DQK/j1+Hcssh1iYTCOKX1p7S3JXXv1PaPzzBMsRajE3b/s7YWu8TCzL1GaHehrumy25wjfCkhkFrsnHgCp0v4xq5TI15fKp47O6fiKT2/JH0559T2jytBVNjHCH0RSwaxOMtyjIOtNa4bWgclF8myBiF3h/9cQtfJcfdw943NAI39TSq+9FmlxhxKhHqq/vGJocI+Vqj/tUQcl0t5LHSXxKzl0i51rXjrEheXG+jNDUdMRevkkBMdlVN/XHr8tXC4RqjUtaL+8V7oOyrmCoC/BPD/AHwLwD+PHb/Wwk4pFYAxv0htGiyppRgqN3cNJ47SWZt0lqf0ubj7+YOgbXZ3yp3PkFt/UmMObSzvPv3ja+qL1wlK68BYK3eJOORYiqGGg9vmDmjCFGPlVLouS8j3zfU6JPQxAznlNhqrwVCjpzsDpMK+BWW6HB42nxocHQFXrwL7++2vub8PPP00cHwMnDvXfE/d8+rV5ngAuHEDuOuu1XS447/xjdNjj4+bvwHA448DL7+8ep8//VPg+eebNF25Uq/M/DS/9NLp386fB7a3m3wAwB/8QZP21H3Pnz/7/T3vqZdWxwc/2Px06XHfHYeHTRnVqgs1ODoCHnigKeunnwbe/vbV5z+GdI4JifrX/qjFPjJqWkDOskptB0fv6bsvuDVofGtcMvHHj2qRRNbE1j6P5dWdt7V1mr6dHfmSBj6h8QXp1nq5vTd6zlR6gCXRZDMB6opRxNSKaEg1EL5wcC9rbJ13KtD+oCr1zW9uNoOEdCYoN1vz4KAR4b29vIFCd23nWzfmbPouXAi7PGLjA35DQbcBjE0caxunPsQeuVK4/I21EeoYFXYp61ZBQqF0NV7qWAPh38MJL93ZKBT+GfM9S3zzoTyH9k9NRftI1oW5887VhihW1nTQl26sHfPX14hTH0OUVYx1e08DqLBLKAkXq90d7pNYfttE50gGtaiQbGxYe8cdq6IVEj33d7pzUmlvI+S6SXXzaeNDexJO8BcLfhG0UGgkF0PPNRw7O6uuGWk9do2H6x1NxWJXfoAKu4TSCR61usNSckU35OeW5DfnXjldZM6dwrkZYufv7a3Gl8cs4JSPn1rFGxurlrbvOuLEmsaFA006reV7IJw/P7RkgHt+1K3DxdhL8kvTzjUQYzZMFBV2EVIxlghiVzPvchsMzvrzXQy5a6/EkIb1OcG4dIm3QiXLInANgytjbhBQ0hBzws6Jq+9CovnlRNnfa5UO8vouIO65UHfScnm2zGga2/RQ+nK39NVgrEHDpMIuRepiGcpiz52JGXIx3Hln2iLmJq+kyiW15jotl8WibF2RnBhvSYNDj6GDny7v9Dgn7s7apmXg3DAhoQ4ZAKHnQq3+CxfK6llqdcqu6LInO8R9BkaFvTZD+dj9CitZ5nWxCLs8Qla6P3NSco5PaA30UiELQYXpnnvkx25srN7HLyc3mJvy729v8zNafd91yI/O3VciQFz5ltSz0klTLn+l9bqrnuxQ9xkYFfY54V6s1KJhdJCRE/hQpIr/ue02+UtChY9zPXQRceNWXZQcS4XMT5O/NG+sNyMpf+76oUYitXuS79fP2dQ6dr3YPIHQOW1XbVSLvSoq7HMkVXk5qyz2YsYG7aQi4AY1d3fDk3Jq9GSorzklrJw/PJRvf4A01KPhIkmkvTdJw0CPp+6rGoOcpQuMtbWE1cfeUCF9Kux90Xdlit0vZimGRMt38/i+YIkIcP5l6iri7i0tM3pcamVCH3rspUvhcgoJNjeG4MYI/DXRU9EoKVeaf59YPD+9nmTPVD99Oe6KUKM/VtEcO5V6FCrsfTDG7l9uQxM6viQSyFnGsbh2aZlxx7URJuq68ddWD103lL+cMY/Q4LdvhVMh39s76966cKHpEXFpjbmkaIPBNWCSRr/NzkiS+jgmS7urtFQaA1Bh74M5D9hIBJizaGOunpg7gr5Q9NzbbosPcHK9gpDrRmqxLxbppQJSz98vI+fWovenLiw/XJK6lLiQ0VC9Sw28pp5xW5GT1qGxGEddpkUtdq8gxtKKhxhTpewCyTNYLJpQSn82qH8+FVCuIVguVwcIQ0vw+mvAuJ+hZxCKQJH42EM+5ltvPU0rtdj9NPni6efDWb9UcLkBUq7H4MraxbPH6l1qOeKuDRPJ9Ut2EeuKrstj7X3soS78GIV+rOmK4dIsWTUwdR2JRRby8TqhibkiOGHzB4RT1mvImk9FmoTmAzhBpz52v4Fxn3PnVvPm8h2q337jwFnsfrhlyj2SymdXFjtXJvT6rkdCG725WuyVmLawcwIw8gIfDTmDeb4AlZRproUTenGo+NFp/ZywhT5cz4FLR2qZ3pDFHsprqCHwGyHaS+EaHHpP38deshRwqhGX+Ni5yJzY/VLnhco2J7a+K0ZuqE1b2GnlGFN3bcxILI6QAJWGseVEaLhzaFQJFax77jlrqW5sNC4Qzr/NCX2qoZI2SC6tvjC761+6dCq21F9N0xFzV8XSRdOXa1FKLGdpWqSTqkoH3XPCR9eYaQu7tXmDPEoDfWG4JQHorExapiFLMvS33I2ffd+576OWfPzjnSuCc4GkGirORRETDmr10lBKNz2fs45z6q6kN0UbRsm1aOPHLUsRG2OQDtbS8yQDplKXUsn7P8PGYPrCTpnSQxoqrRJhoOJvzNmZl/QF4gYg/eNyNl9eLuMrPEo+IVGiG0qnhCIVBkiP9+9J19SJ5Ts0rsClifrYJS6SUD5DFjcXohm6Zsjil+yOVTPEsWRQc6bG4PyEfSoMXaEWi8ZtEXoROHENxW1zg3eXL68eJ53NyImNbzlT37pbe10i2FJLlktHbPEzzrrPWVDLP98Ja6lYc2mPuZE4MebcmpJrtnHrtKXknZppKLIKu6Nv67mLCiXNAxURX4x8QjM46QtEfd/u2NKBNbrwlrsW57aQruHitrbj1qiJnUPDFEM9nVDI4GJx1sceIzVGFLPqYy6SEr+0/ww4d1athiUFfe41LXx3vFrsMxX2IR5u7Xu2seZCXf/Yy+m/QFyoXmn4aUi0ueNKGzE/z1zUi58GuoTCcslvD5eKBY+lz+WDW96AijVn1adcJLl1i7rQ/J5WadRLSf2mz4FbMbMGU3LfClFht3a47ljNCpWTB0m33x0neTmXSz7qJNZVl/p8c0L1uGuGIki4noYT41QauHLxLVzaSITS5g8Qb26uzl5119vaOhV5LpSRptfF8JfWq1hsfmqdf1pObep3Kh1KEBV2a+fRHcvNg9TXLO0Kc3t8+sdxro02jYn0eNrobG2dnYrvp3lvT+4fDrmG/L1Wl8vVwVf/nNgAMRV590kteUyPrWEp+2MHbZbmpeWWkw66I5S/6ffMrO0aqLA75lBBusxDqrvviyF1oYSEJxX5IcmHxLp2YrSxsZouJ1pbW2fDK6kI5/QyQj0Cyf9cOi9cWBWzmNXMuYja9D65Br3NPBHXkOY2DH79om6z2i6oGaHCrsjgIkToYOI99/AvVKhLnWNRhl7W1MtNxYjG7LvjQm4Z3+8t6WX4n93dsHiHrPULF84uBsY1CjGXWGoLwljZhkIouca71N/eptGhdYBr1OfQ+65Ar8IO4D4AXwXwdQC/kjpehX1EcFvpcX51LkIl5IOXdqkl7pZQRIjvPoiJIzf4yU3ScmLE3ZNapP6uRtJPaJmEzc3TNKViw/0NPyRRJJzo+nngxLxNhEwt0eXqxUzDF3PpTdgBbAL4MwB/C8A2gM8DeEvsnFEJ+5y7dzFrmL7EId+vL0ypwUX/I5n4k/uycmGBXCQLzSuNQ+dENuYC4MqSCr4bCKWNjd9YxhrBnBBGP7IlFF0UG6AEwuWWO1jvN7RtBna5a5eGec6YPoX97wD4hPf9UQCPxs4ZjbDPubJwgub+7vvMpRanG5yMWbnU7SERXUkIJHe8ExJ/SQMXwhjbSo7rZbjdlWpE7nAhou7DuXBi5cpBjw1NjvKfv8RizxFQzk/fx7vTxb0mZtj1KeyXAPwX7/s/AfAEc9yDAK4BuLa7u9tDEQiYQ/cuVDFD8dc0z/6EISoAt9++eo2U64R271MiEVozPZRPt5wttYz9uGxfwFK9BN/IyBW1UDo5d5ErD5pG16BJw0+5HgHXYPqzhkPunlSPrsR9NiUmmJfRCbv/UYs9kaaSsDGa/pCwc+JLoxOcKIW2Uovt+Zn6TpE0rr6gcXuCxj6SWHVp+UtdNbHfuYgY9//Uc5f4tEOusVrx4XMwhhxD5KVlD0FdMVLG1BXL2UXe2njFpGIoFSEaAx/zdUr2/JTk2c8D546RbMTBrRQZEt+cDUb8c0KuJangc8IcWjsmFSlEZ866/6cGs2N5TPn2cwdXx07fhl2F+/Up7FsA/hzAXd7g6d+OnTMqYR8Ly+WqpZWaDZhjfXIvb0lFiw3KSS3P2PU4C5sOdvp7d3Jx6ZzLoTTWOhZZ4q4RykPMNbW5eXadGalw+ptwS59NLDzST2NoMJbm49KlcBqmRp+GXYUeQt/hjgcAvnYSHfNY6ngVdgbO9x2yhlOCTQkJeOlgYchiz7HmXLpDC5LRe0lEKpXe3LyGhHJ3d7VXIy3fVAMbWwZZ0pDTvKYGpWkaucHY0LjM1C32vpmSxV7yGY2wj8kNQ19ubmJNTlc45gKQuBBSaeUaFmlD4buJYkJGr+esxJBVzj3LkDhfuBAPz3NWPhdZwrnLQm6sNqGJ1E0jKd9YjH4on6k19f18hDbdqPEujel97Iqp+NhLPqMQ9r79a9I0+YLJRVdINraIuQBS8dlt0s65O7hGgDZabkCRS1soH+4jyV/IncKJJz1nZye+B2uqwfXHTUKx3qE0Ur+4pM5KLHD6zC9dOp3HkIoGCg2oc2GrOXVrjO/jCFFhTzGF0X33YnCxyzkvt8T3nSvwsUbIj5hxf9/YOE03nQx1xx1h6znWI/At+ZTFyYVIho6n97nzzrNzAlw+trbOrknDibIkksVPY2wwmusVcOUV6yVwjaVvNEj3rfXrCueqyR1kncL7OAJU2FNMyULgXsbYCy7NW8wKy0lPaIJNbKA1JHaxBoazbKW+fSeKe3tnRThlsfvHuVhwJ740dNH/zi1oJRGuUP5pLyIl7qEy5BqtWC9BwnK52rClGlvuGjXfx9Ke6MjdQSrsEmo+xK4rRMhCTnWdJZYdtaBTLyEVB3/ikm/1cS+8/7nttlWRkSz/6/vYU3uE0rzShjG2tygnTpyrw103Ni4iicrJEWQ/L7k9Lb8sQouk5UInmpWERdZ0C5aOHY3c2FNhr4WksvVdIWp0W0PWtER0uPwuFqcWq5vc5KJe/I0lfH865zNvY+WFejKclbq3l95P1c+T81Vz9+TKJbRpdKinlapDoUYyZ95D6Fn6PZA2IYz0XZEKdW2jqPT9mIA7SIW9BlLB7qNC5Frskutx3Wepa4eKVGxANLYGOhfFIc0XLXfOFULTH/vQTaxp+UgGBd3/OPdU7Lml6tByGV7HvW29C02M68Mt0YVRlDvRr8u0VEaFvQY54Xs1KoTEvyrxRUuJrdPiT4QJDcb61mlsdciY4HADb9JVAv00pLbwWy7j0S20QeF6NDlrHHHPLDaonbsSpjFNb6hGvUs935p1TnLvNvjplg4G0/PVxz5zYc8R7LYVoo011wYu3dSK5gQkNjBqDG+xc/7gWi8iXSPdH7x0/w9Fxdx6Ky9stLGK1VuuHDnXRCh8M+Qv54R/c7PZ/OTuu1c3xi6BW9YhNo5S05oN1fvS92kC7pQ2qLDXoq8WnHuR/DT00UV0eaWrOu7u8iLlh/r5flon6hsbZxsFahXT+7bNF3UP+V3y0CCui7vm/My+j92Y5ntIwKWNk3++xP3CCT+dresvTVCCxGJvs31eilQDmJO3CbhT2qDCPnZiQhmyONv2COj96CqNoSgZ7nktl2cXGXPrh1y4cPbckPujRBhiZUH/x61yGFr1kPrT/Uk9fuMQWvSMi5SRWJ7UnZUa9HVlxu2B2tYFmLKaa1vWMdpa3ak0jdzdEkOFfcyEXhJuqd0aIWOc5UcbEXpvGg1C4Xzj/nl+o8ANXEo21aCDrNx66y4Shk7957b7Wyya9NBxBc6tFJtgRAVnueRnBEvEcLFoGj+/Z+H3GPxrug1BQgPOnMUvFTDJse4Y14jX8vFz92ljscdcYhO36FXYx0DoZQlZJKnubyhqxb8fV2np/TiLj0a1pHy3MQvf39HIvUw5E2G4hii0SXTo704k3THOjcK5TDhhdu4Nv9GILVNMLXvXG6J5DvUAaJ64BseJvhNXv0HY2DjrkpFuYCLFF/VUQ8edl3v/kvO4+t+nS6kHVNjbUKOrFvO7xqyGkHVBrdQ2DYbEYpdYfzErmh4b2viDQ9IQST60cYqNYywWZxuora3VzT1811jItUIbAuriSYUs+uXPNTg0jI82em4OQWixrhKom9D/xNaY6cM6To1Z9DUI3BMq7KW0rYyugqX2ppQ2HlRwfCu1tMFw30M+9tjEm1RPgIuAcffyF8PK6XlwQkU/3P/pevb+uIA7x/drx3ZmkgojLRMacUOfoz8A7U/w4gaAOXcP5yqiDaHrAZQSK5dYzy7WkNaAxqunFijrMmyzJ1TYU+S6SSTX86eNb2+vWmclVhMXv5zbYEgrsjuOviCxTaz94/xZpu7juxtydrLnfOyctbu11bgkFoutKxnqAAAU8klEQVRVNwfXkHAixYkAt1+ptAy5MQUq3r57iDaAoRUyQ8IV29JQskl4ClpmFy7INtpYLs82uJK19HPqamhco7T+TwAV9hgpy7ZkwJJ7mf0Fp0q7fdSlw4lnrMEoyU9O99UXnZRVTdOa+8ItFk3o5d5eeP0Vd83QGiipgUcqtDFBjblkqNXsr4PjniMX4ZOyckPCFdrSULpufQy/8eAWTYudl2Pc5NTVWCTSjFFhjyGJH86p/Fx3uGa3j7NepQ1GSQ+kpPsaix7hLLaQWITuw6UplreQTz/2rKRlKG3o/AaZNjS0EXRlE7LAS+pQjYW5aL5yzuEip2rVVWrwtO2VTAQV9hglVqz0ejnuhrb3zHUN5OQ190X277O1dRoSF7I+OeGNpZV76WORH7GGI/dZ0XTRtN9552nEix+XfuHC2Qgh32LnGj5XPvT6oTGVVJrbLqXLXbNNnZC4bnLq6oxcLFJU2FPUrhRjrmSl3e/cGGiu2x96UTlhj1ls9FoSP3LM+s8VqNBAs/9JRbxwuy1xwk4FscTlwFnLtBeRY+X2aSDEXFxjfcd6QoVdWaULKz91bExc/UFVKpjcgKV/rdrutNB5ofxx0Uoxa5zrUfiRQtR37TeWOWMq/vlcw8dFkkjKiSvvkL+/ZPAydUypn39mqLArZ5H4JN3LlTOJozSKyL+f/5K6VSVTsxpjL3ob65Kel7s71NbW6rIKsbVcJKLXZhlaem1u0FFy3TahhalnETrGT3/OPIgZo8KunCUVReC/XLEZlpSa4xW08UlZ5NTiD+W1NAbdt6jpffy0UsudzrItCXOleY1t0FHq95b2BDijgCtj6d8ood6AX69q7fQ0caTCvgFFxtER8PDDzc+pcXQEfOMbwIb3uG/eBK5ePf1+9SpwfNz8/tJLwL33AgcHwNvfHr/24SFw5Qpw+XLz8/CQv/+73tV8XPlx5emn4eWXz17j/Pmz369eBW7caH6/ceNsXvb3gXPnmt/PnWu+S6Dnfec7TTk5fvInT/Pn5/vRR0/P29kBvv3t03N2duT35zg6Ah54AHj22SaPzz23+r8nn2x+Suqmn+777we2tk7zG0onfS7Xr/NlLP0bhTvGv+fxMfDa1wLb28337W3goYfSeV1nJOpf+zOIxd5m4KWmVZqDv9lFKdTX6q+THnNftAmNo/enyxeErh0bNOQs9hLfviS9oe6/X285XzznxuIsy5gFTomFDLZ1g0kHUiWuklC5hP7G3SM1tqGDp+qKOUNbYW7zApXCbX5QAudekEQd1Moz54uOhd35g4YhVwuX3hiSgbnQAC9d64WGNeYs6xC7Zki4lsuwy6RNvabPRbLAXNuJTrkMKeQjbURU2H3ailRfFrtfmaj43X133jX8NVtqDSTG0hu7DrXY6QqBoYHcGlEQqYHHlDVKfbucbzpnHaDQuit0CQM/HX5DEIqeKemZ0JBLaf3o8n0YQlAlvYWRoMLuU+NBdV3haBol4he7hv+ylsxcDLkL/IZDOvmEXovrRVBCjXHOc4hZu7H7+OUYWlqXDphK0xdbUIu7Xqys2tZJzoBIzfbMjZrKTU/fgsrdc4geupBehB3AewF8EcArAC5Kz5ucj70PuBd4b69ZG0XqhgmF4NXqocR84DkrCPrWOHVHxNKQ++KnIoGsXXV57e3xfnJ/OV7X8NKegLSX4/JuzNmFykIWe+pZtB0DkVwn1NilenM5RsUQgiqJyBmRZvQl7D8B4E0A/vvgwj524U7hV6bQS+6Oi/nIS7vXPqEXjP6ds1qlzyG2ymIorznjBe58bi/TWF658g9txp1KXyx80J3n/+7H8Kf82TVFUPLMcsdqSurhWCx29/cR6kmvrpjBhX3ELWwWqa5ujt+7zR6pUivRX2kyN5KG+t6lroBcqzHlp+dEyBeu3DyV1sNcMSy5Vxuxyrlfm57jGHzsI2Z0wg7gQQDXAFzb3d2tm9sR+8SKCL1EfeYzVNnp3/2GJNdf67tjpBZTjp835J/m8sAtAcxdI2WF0w1GQoJBQ1klYhgqe6moh4RZep2c42r0HJUVqgk7gE8BeIH5vNs7Ri322nAv0VjzWfIix8LtJCIUi4X3hZVGlMQs8FCZ+1Y/1xPilselm674YwhcKGuqDNs++9hgdBd1qkbPUVlhdBa7/1EfewucddnH0sBSON97bMs0l4fQVnkSEdrePh3oDIkft4BWaLp/rIfgb19HxZrbYCS1EFgolDUmhl2F7M6ttztz1k/Y14UxWO2cSyC16TKXft+SlvRMOHdF7P/cJti+OEuiWUIukpCIcys+UuHkInFSlm2N5z6lXqDC0ldUzP0AvgXgJQDfAfAJyXkq7C0Y2sIKCUFq0+VQ+kP58K1X10NZLFb98gC/jjkdzKVCGvNvU782F63ENWJugxHaq3Dn0OVzd3dXrX2/seFEuIteKXVfqQtltOgEpbnSxsKqIQwpEUxFksR8yVxPgG4Tx80EDW2t51w+e3tnlydw1/Q3lI6VKXdN3yXDuYT8fHA7PXEuGb+hGsKK1kHP0aPCPmdKBLpWl1t6nVgaOauQuy5n3e/t8VY4Z/FT652LpnGTlnLLlAp8rCy46BwuD+46qYifriz3VG+qi/uuy9hYJVTY50DNSt/VhJYaaeSsV3dN6vpwkSlUVKnrIDYAS11FqQWwQuSEQvruGTqA69xCsUaOu1ZtazrVm6p9X/XvZ6PCPnVqV/ouXDiha+aIPRe5Apz1RacsSOluPv49qW87tIJkqlx8wY5FKlELP+XqCDWeXY+x+O4rPz9d3Hfo8aIJosI+dbqo9CUDciGXRkhkUtYm9aFTgaX5pffg1qThXB2pAUCuwSjZlaeNS4brKXBlFBoQ7sLKjTXWtaNy1GLPRoV9CqT80DmVviu/OxWj1L6XORNhQj7dkKCEVpH0j5Fu68fdu3S7tdAM11g6Q4O0kvLkxihqEctLqdvNNX7cc1EfexYq7GNHIqrSSl9q+Uh6BZzflROZlBWWsu5dqCC9lkTAnHC4j6Sns1yeXauGW/lRCh0LiLl1Ys+UK6PQevJdWbu0XNqshU/TuW5ulw4aLRX2sVPT1VJ6LanbhPNJ50bE5PrjpcJFrXU/bjwleH6D0Pblq7FGOed2Ca1O2aV/uuZ666lJZXOlo4ZXhX3s1Hzwba4lEWH6opdGkeRYMFLh4oQjNZDZhtCgZi3Xmf/32P6ptEGrmd+u6mbtdI6ZjhpeFfahyBGvml01ziVSem1aKf3Y8dLY9dx8SHeQ6rOrT0WKrnVT23WW2hg75rsuyVut+pO69jqgFvuM6Mrv2Xc6OD9rX2uZ+D53qVDXFLgYocFeTnRzrhMbB0itJ9/WMuyr7NaRAX3sG5g7R0fAww83P7vm6lXg+Lj5/fi4+T4EbdNxeAjce+/p95deAq5fB554ovlfjXtyz8W/xs2bwOZm8/u5c8D+fjy9H/848Du/A1y+DFy5Ek5nG86fB7a2mt9d2krY32/yBMTzdngIPPNMk6dnnuHzJL0Wx9ER8MADwLPPAjduNH8bst7OjcPD+DvTJRL1r/3pzWLv24Kei8Veco2cwc5Q3Dc3eJjyQ9dAcj0/bZubjasoZU23vWfX11rXgc0JA3XF2GFmto3Fn1jT3y29Rur4mD/cnduHy6fkerF48qGfdSnrOrA5YVTYrR2PBa00hCzEnNmUtRtryWxQa+dbl6beOK0ZUmGft4/98LDxt3bpd1Xk+P7g7W3g4KB5Ltevy/3zbX3K1K+/v9+kxfHCC43fmY7JzLUuSf3Aftn1OW6lFGGaRqBfLl68aK9du9b7fSfP0VEjevv70xUWLg9uEO/4uBHrlHCWlIN/j50d4B3vAB56qPnfL/5iM0Dsc/lyI3hdMIXn6KcROFt21jaDrZJnpVTFGPO8tfZi8kCJWV/7M9twxy7p0xUwRHx91y6BkBuIxop3Xb59Pcc25UnTyG1u0ue4lfIDoD72mdHXQHBXsw77WJlQmhY6+Yp+7zJdfTzHts+QpnFj4/T3rS2NeR8QqbDP28c+J9r4lnNoEwNPfa/0WkdH/cb5++lxPvKDg1Of+rlzwGtfe/acn/1Z3rXgrvXYY6v+Zcn/3N9Kn2OOX7vtPAY/jZubwCuvnP1f1/MFlPZI1L/2Ry32QvqIYCi19rjzSi32WqGakgXOJPnlrH2Xn9jmGaFrl4SRxpYu5txdNeYxuPDTOUYDTRSoK0YppkRYY+uw56xBUksMc1weqWuHlhIIbRLi7lfL7UKv4282Uqu8YmhI5GhQYVf6pZZvPrVuu/TadKZoaBExCZzFzm3nJ7XYS+4fWjtniEl4ymBIhX1rYE+Q0pYhQ+fova9caZ+W8+dXv3M+49T1Dw+BRx4Bfv3Xm3VnPvzh5u/Xr4fTFypLP2/nzzfXOH++uebxceOzv/de4K1vXb1+jTI5PATe/37g8ceBl18+65vf3weefvo0VLSrsRdlWkjUv/ZHLfZKdDG9PmfJ4S7uzW3yENpFKAW1ZmPXKO0VSJZBqMVQoaLKaIC6YtaAmt3wXGHr6t7b22cXCPMH73JdKtQdE0tvF7tQTR1tMEaHVNg13HHK1AyBzA2R6+reN240s0JdOJ2/3MDLLzffpfjLAHzgA/H0luZnLEs118bN1H3ySX6JBWXcSNS/9kct9orUsqrauCK6vHftCVOpiJzc/MzVYtdB2VECocXeaq0YY8yHAPxDADcA/BmAf2at/X7qPF0rpmNKB1RL12CpcU7sOm0HiLseYE5dn667kpP3ochdv0fphV7WigGwD2Dr5PfHATwuOU8t9g7pe02ZNmGI67D5iX//nZ3T6fhuvKCP9JX2rLro4YyVieQFfQ+eArgfwEclx6qwd0ifXeiSe/XdxR/apRDbJ3Vri48CqklXDcfQDWZNJpQXqbDXHDz9ZQC/X/F6Sgl9rSlTeq8+0zfE/WL339kBNrxX7ubN03R1lb6uBnfnNGg8p7yckPSxG2M+BeBHmH89Zq1dnhzzGICLAH7BBi5ojHkQwIMAsLu7u/fNb36zTbqVGH34bN093ISdrv3ybRjah+3f/7nnTidNOd81cDZ9NdPbla98Tj74CeVF6mNvvdGGMeaXADwE4B3W2mPJOTp4OnEm9CKMktRAce2y7aphG7rBrMlE8tKLsBtj7gPwGwD+rrX2e9LzVNgnzsMPN/HNji53G1o3tGxPmYjY9olU2Nv62J8A8BoAnzTGfM4Y859aXk+ZAkP7reeMlm2DTpBqRatFwKy1P1YrIcqEqLXgl7KKlm1DycJvyg/QzawVRbv840PHcVikrhhdtldZb3wBefppFZCxoD2XVqiwK+uNdvnHy+GhPotCdHVHZb3RwUplhqjF3gXqs50O2uVXZogOntZGB30URemIvuLYFcoM151QFGVaqLDXRn22iqIMjPrYa6M+W0VRBkaFvQs0TEtRlAFRV4yiTJWjo2bRMF1HRSGosCvKFNFFspQIKuzKsKjVWYZGXykRVNiV4VCrUwbX+Gn0lRJBhV0ZDrU604QaPxd9dfmyToJTVlBhV4ZDrc40scbv8LDZXUlFXSGosA+B+pUb1OpMo42fUoCuFdM3upaMkotkUTldeG4t0I02xoqu/63kkprwppuFKAR1xfSNdq2V2uggtEJQYe8b9SsrtVFjQSGoj11R5oD62NcC9bErytjoUnx14TnFQ10xitIHOstW6REVdkXpAx3gVHpEhV1R+kAHOJUeUR+7ovSB7qyl9IgKu6L0hQ5wKj3RyhVjjPn3xpgvGGM+Z4y5aoz5m7USpiiKopTR1sf+IWvtT1lr3wrg9wD82wppUhRFUVrQStittX/tfb0FQP+znRRFUZQztPaxG2M+COCfArgO4O9FjnsQwIMnX18yxrzQ9t4j5nYAfzV0Ijpkzvmbc94Azd/UeZPkoOSSAsaYTwH4EeZfj1lrl95xjwJ4tbX2V5M3NeaaZFrsVNH8TZc55w3Q/E0daf6SFru19l7hPT8K4FkASWFXFEVRuqNtVMyPe1/fDeAr7ZKjKIqitKWtj/0/GGPeBOAVAN8E8C+E5/3nlvcdO5q/6TLnvAGav6kjyt8gy/YqiqIo3aFrxSiKoswMFXZFUZSZMZiwz3k5AmPMh4wxXznJ3+8aY3546DTVxBjzXmPMF40xrxhjZhNaZoy5zxjzVWPM140xvzJ0empijPmIMea7c50/Yox5gzHm08aYL53UzfcNnaZaGGNebYz5Y2PM50/y9mvJc4bysRtj/oabuWqM+VcA3mKtlQ6+jhpjzD6AP7TW3jTGPA4A1toPDJysahhjfgLNgPlTAP6NtXby+xwaYzYBfA3AOwF8C8BnATxgrf3SoAmrhDHm5wG8COC3rLV3D52e2hhjXgfgddbaPzHGvAbA8wDeM4fnZ4wxAG6x1r5ojHkVgD8C8D5r7WdC5wxmsc95OQJr7VVr7c2Tr58B8Poh01Mba+2XrbVfHTodlfkZAF+31v65tfYGgN9GE8I7C6y1/wPA/xk6HV1hrf1La+2fnPz+fwF8GcCPDpuqOtiGF0++vurkE9XLQX3sxpgPGmP+AsA/xnwXEPtlAL8/dCKUJD8K4C+879/CTIRh3TDGvBHA2wA8N2xK6mGM2TTGfA7AdwF80lobzVunwm6M+ZQx5gXm824AsNY+Zq19A5pZqw93mZbapPJ2csxjAG6iyd+kkORPUcaGMeZWAB8D8AjxCkwaa+3LJ6vovh7Azxhjou60TjfamPNyBKm8GWN+CcA/APAOO8HJAhnPbi58G8AbvO+vP/mbMhFO/M8fA/BRa+1/Gzo9XWCt/b4x5tMA7gMQHAgfMipmtssRGGPuA/B+AIfW2uOh06OI+CyAHzfG3GWM2QbwjwAcDZwmRcjJAONvAviytfY3hk5PTYwxd7jIOmPMD6EZ4I/q5ZBRMR9DswTlD5YjsNbOwkIyxnwdwA6A/33yp8/MJeIHAIwx9wP4jwDuAPB9AJ+z1v79YVPVHmPMAYAPA9gE8BFr7QcHTlI1jDFXALwdzbK23wHwq9ba3xw0URUxxvwcgP8J4H+h0RQAWFhrnx0uVXUwxvwUgP+Kpl5uAHjGWvvvoudM0EugKIqiRNCZp4qiKDNDhV1RFGVmqLAriqLMDBV2RVGUmaHCriiKMjNU2BVFUWaGCruiKMrM+P/yLrEfa6PVHwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# First set up the figure, the axis, and the plot element we want to animate\n",
    "fig, ax = plt.subplots()\n",
    "ax.set_xlim(xlim)\n",
    "ax.set_ylim(ylim)\n",
    "\n",
    "y = values[0]\n",
    "paths = ax.scatter(y[:, 0], y[:, 1], s=10, color='red')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Reverse the list to view animation in reverse.\n",
    "# values = list(reversed(values))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "def animate(i):\n",
    "  l = i//48\n",
    "  t = (float(i%48))/48\n",
    "  y = (1-t)*values[l] + t*values[l+1]\n",
    "  paths.set_offsets(y)\n",
    "  return (paths,)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "anim = animation.FuncAnimation(fig, animate, frames=48*len(cs), interval=1, blit=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "anim.save('anim2.gif', writer='imagemagick', fps=60)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<img src=\"anim2.gif\"/>"
      ],
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Image(url='anim2.gif')"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}