Skip to content

Instantly share code, notes, and snippets.

@smsharma

smsharma/NGD-Pyro.ipynb

Last active September 14, 2020 16:08
Show Gist options
  • Save smsharma/5836cb7556df6ac99cc64a776ea8717d to your computer and use it in GitHub Desktop.
Save smsharma/5836cb7556df6ac99cc64a776ea8717d to your computer and use it in GitHub Desktop.
Download ZIP
Display the source blob
Display the rendered blob
Raw
NGD-Pyro.ipynb
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import random\n",
"import string\n",
"\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"import pandas as pd\n",
"import tqdm\n",
"\n",
"%matplotlib inline\n",
"%load_ext autoreload\n",
"%autoreload 2"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"import torch\n",
"import gpytorch\n",
"import pyro\n",
"from pyro.infer.autoguide import AutoMultivariateNormal, init_to_mean\n",
"from pyro.infer import SVI, Trace_ELBO, Predictive, MCMC, NUTS\n",
"import pyro.optim as optim\n",
"import pyro.distributions as dist\n",
"from torch.distributions import constraints\n",
"\n",
"pyro.set_rng_seed(101)"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Text(0.5, 1.0, 'Samples from poiss. distrib.')"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAmgAAADfCAYAAABCvy+tAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOydd5hU1fnHP+8uSwdXYFEXWJaOIAIKWFBjhahRscbekmCJMRIlwRIlRhHFgi2JmJ9JUGOswYKKETA2FFFEVLo0F4KrsEhZYMv7++PekWGYO2Vn5t6Z3ffzPPPszsy5935vO3PuOd/zvqKqGIZhGIZhGNlDXtACDMMwDMMwjF2xBpphGIZhGEaWYQ00wzAMwzCMLMMaaIZhGIZhGFmGNdAMwzAMwzCyjEZBCzCyBxEJNdhVbXqvYRiGYQSG9aAZ4ZwGTAFqRWSOiPxFRNomurCINBOREekSIyLnJFCmkYjc5b4+FZFD0rX9iO3sIyJHRnz2FxF5LRPbM4xcQESOEZGnROQhEXlMRH4nIqN93H5r975P6wOliFzq3t/TReTKdK47k4hIOxFZLiIX+LCtq0Vku4iMdd/3EZFVInJ0nOUS+p0QkQIRmSci17vv/09EVER+kpYdyAGsgWb8gKo+B9zrvn1FVS9X1e+SWMXBQFoaaCLSDBibQNFhwGXAGGAukJ+O7Ufhx8CREZ/t4b4Mo8EhInsATwJXqOpVwM+BQUALvzSo6vfAb9K5ThFpDTwKTAKmpnPdPtAIaA00z/SGVPUB4Iuwj5rg1IfN4iya6O9Enru+Vu77XySrMdexIU4jLYhIO+A6oDwN62oG/A7nho9He2CDqtYCl6S6bQ89HYFfAq9EfHUuIJnYpmHkAD1w7r8OQIWq1orIOOAUn3XUpnl9bXAaB+tV9d54hbMJVf2fiLRX1RqfNvnDsVfVuSLSJta2k/mdUNXtItIttD73+kqH5pzBGmhGUojIr3F6rT4BDgFuV9WZwKXAgcBmEfk78FdVfVdEegF3APOBnjg9c0+KyFPA2TiNqr2BHwErVPUK4DicJ6x27rreVdW/RtHSD+epPVRuIk4P4FHu6yvgTaCHqoqIHAE8CywAHgO6AOcA57iViwDXAIOBRa6GU4Ergd5AIxEpBcYD1cBrQAFQ6uo5C+fHaTHQB/izqr7lHrM7XW0VQFdgf2C4qm5K+iQYRnawGKgE5rhD/e8DTwP3A4jI4cA9wL+BI4CpqvqQiPQE/gMozoNPJ2Ak8ACwA6c+GAmcq6qfiMg/ce7Tf+PUIwOB/wGjVHVLpCgROQGnPvoS6ItTR33iDpVtBpoCFwH9IxsTIlIA3O2+vVtEnsO5p88GLsbpSa9U1UtFpAtO3bYYaIvT6LgN5/5OeP+i6O8JvAG0w+nF2+4ev4dU9Sm3jFddczlwr4jcpapjReRkYACwCadevQuY7epc6mrpp6qnRuqIhogMAv4IfIjzAL132HdjgVtE5BJV/buIjARa4jzEng6cB5xJxO8E8I17vJbjXD8jgd/j1K3/JyIvqerFYTIuFpH9cervKuA6Vd2aiP6cQ1XtZa8fXjjDeAqM9fj+D8BN7v9jgVVh370F/D3sfQFOQ+d6931LnIqiCGcoUoFH3O9Odt+3CVv3igT0XhxeDqdyVOBI9/2lzmX+w/d/BDYAfdz3nwAPuP+fi/ODExoiWAUc4n63IvKYuOte4f4/BOdpsth9fwCwBShx378DfIrT/d/a1Xha0OfbXvZK5QUcD6xxr2d17+8h7ncnAVPc/0P1Sjf3/c/c+2O4+/5e9/1+7vvw+zLPXfY6931jd5vjw9ft/t8dp0ETum+HA3Pce25ZWP3yeyDPY59K3e2Vuu9DddVYnIe219zP5gE3uGUE+Ay4Ntn989BwqbvNdu77nwA1wH4J1DVvh+oq4CPggLDzcQBwIvCcexwFuDnBc12I05i6OOy4rCWsXsSpM0PfrwA6uf+PDNP7FmG/E2HHawNO4/tx4Ofu55PZ9TdFgV+HHfMFwLig74NMvcyDZiREmOl0AlAuIqOAjjhPh170xuk16y4i1+D0ds0A9tSdT65z3b/fu39bpyg1crgj8n0NsF1Vvwzbbmibp+M0uL5X1a2qWqKqsxLc1uk4jbu17vvFOI2848O2u1BVK9XxzUDq+2oYgSEiTXAePDri9IqMB7bhPMSB49+a7E4aGOJ+FqovanDuj7fd998Dm1X187D3rcEZ2nI/+9Z9vwOnByeaWfx4nIbHkW6dMwj4GqcBswX4RkRmA2Vh641JWF21SFWnqOrxOPXa/ji9UKHW4WLgjGT3z4Nad73fuu/fwWmo/pj4dU34fs0GPhaRhcBhwOfAQpxRkO+AaTijDIlwKM7D9SeuthqcxvBuul2+AFaKyDxgzzC90ajB6dn8TFUv0J0jJtHO0UZ3+4rTAXBSgvpzDmugGXERkb2A41xT8Jc4ldP9wLsxlhmGc9OBM0Q5UVUn4jyBLgkrui2B7Q+rq3aiTxrw2qbi9PrF0tJaRA6O8lUBu95PIWNr4wS2axi5yCHAVapaq6qfqOr1OIb9Evf7F3AsA38Gnom2AlWtDHtbGa2MB3k4NoNIQnXOv9w653Ycm4LgNEouxel1+bOIHJTE9sAZVg0Rqici7/nw+z2V/YsktJ1qEqtrQtwCnIDT63cxzvn4Fqcn7kacOu9VEWkMICKxJlmFvmuaoOZLgLOAWe62zowsEFG3f6PJe+e2sPOcx9Ofc1gDzUiEa3CGFEpxnoCfc58+SwBkp3NzM9BURNoAx+I82c0Hhoat67fueuIR8oqAc5MnSsiTEqqsDnU1JnKtvwiUikixu0xfEflxhJ7uOMMEkUxxvy913/fG8ZtYGA6jPnOBa/wOsQ/OEBY4vtLXVHUzO++LVFzebeGHnrshOPdcJK8CW3HrHPcH+yGceutBVZ2sqhfheJ9aRVk+FuGhPBbg1G/7hm2nJ45PLm2EhTk6HMdv9RrJ1TUvAtNVdRRwFU6v3anASar6gKoOx/F+FYhIf2C9iJzoIed9nGHIXq625ji2Fa+69RXgeVW9HMdaEuoxjPydCJFoqJQW7vYFp6E5xX1/GU4P6V4JrifrsUkCxg+ISMgICzBCRLrhGEz3A67G6R5/GLhRRKbjNL5WAtNF5FjgQeBmHGPwfapa7ZpU7xORv+F0TX+uqstF5FZ3O5eKyDtAKHbSrcCFOE/cI0RkEk73fjS9pThPxO1E5EFV/ZWqlovI/e56B+KY8gH+4X5+slv+Cpwu/l5ABxEZDjyBMyPtSRH5yF0uNFxzF45PYh/genfyQ2jbl6rqY+4w8L0ishynIXeyqi4VkfPCtnM0Oxusl4rIW6q6ItZ5MYwspQanMXarG5qiGudHdpT7/W+B89wf0jKcH9KHROQMnHsHEbkNeIkY96WqTnPXN8yd4b0fTniPO0SkFW7dISK3qurNbgPjRnFiIubjeJqqgY7uLNONOENzoeHHSEJ1053u8OzP3fejRWStqi5S1SoROcnVMBGn8fgv4K6weinZ/YvG70RkE47P7lRVXQQsilHXXOquu7U7EavS1fgdTt01Aaexd7SIFLnH5R+qukVEqnCGXiOHLQFQ1e9E5BRXUw+chpIAvxCR/+CEz2iHU6+9gtPTd6eIrAeKcRrFEPE7EVGP36mqv3OP3cnsbGgPVdX3gNVAS/c4dseZXHCHu95KV3+6Z/UGhrhmO8PIWdwfgLw6dI8bhpEDiBOI9hJV/btP2xOcnqFaVWcKePh7H7Z/MfA3VW0wcSWsHt8d60Ezch63wrSb2jCMtBBZp1gdk3nsGO+OedAMwzCMrCXCDlEaoBRfiBgivTVmYaNeY0OchmEYhmEYWYb1oBmGYRiGYWQZ1kAzDMMwDMPIMurVJIF27dppaWlp0DIMw/CRjz/++FtVLQpaR6pY/WUYDY9Y9Ve9aqCVlpYyZ86coGUYhuEjIrIyaA3pwOovw2h4xKq/bIjTMAzDMAwjy7AGmmEYhmEYRpZRr4Y4DcMwjNxgytwyJkxbxJqKSooLmzF6eC9GDOwQtCzDyBqsgWYYhmH4ypS5ZVz/wnwqq5zA8WUVlVz/wnwAa6QZhos10OoRW3dU8+mqCpZ9u4W1FZVUVFaRL0JBfh7FhU3p0q4F/TrsQfvWTYOWahhGA2bCtEU/NM5CVFbVMGHaImugGYaLNdBynG83b+fleWt4ed4aPvt6I9W1TmaIRnlCYfMCahW2VdWwdcfOyrDXXq04sncRZw3qRLeilkFJNwyjgbKmojKpzw2jIRJIA01EjgR+DPwD2AjUAt+r6taIck2BlsAOYItlud/J0m828/DMpbw0bw01tUqffVoz8oiuDO7Shj77tKZdyybk5wkAqsqGrVV8Vb6Zj1Zs4N2l5fzfO8t55L9fMaRLG355VHeO6NEOEQl4rwzDaAgUFzajLEpjrLiwWQBqDCM7yXgDTZxf/YuBbkAJ8AhQCvzOfYHTQDsMmBW2XAfga/dtGXAcsCDTerOdbzZt445XFzLl0zKaNsrnokNK+engTvTau5XnMiJCmxaNadOiDYNK23DFkd0o37Sd5z7+mic+WMlFj81mSJc23HDCvgzoVOjj3hiG0RAZPbzXLh40gGYF+Ywe3itAVYaRXfjRgzYMGAd0Bc4D3gTGAufg9J4dCGxX1VkRyxUAdwAvAZ+q6jYftGYttbXKEx+uZMLri9heXcvII7oy8vCutG3ZpE7rK2rVhCuO7Malh5Xy9EereXDGUk7903tccmgXrhvek+aNbfTbMIzMEPKZ2SxOw/DGj1/hj4DbgG3AfKAp8JSqrhKRJsCJwK88lu0GDAaOEJFHVHVjZAERGQmMBCgpKcmA/OD5dvN2rn1mHv9dXM5h3dtx6yl96Zom71iTRvlceEgppw7swF2vL+Kx95bz5oJ1PHzuAfTruEdatmEYhhHJiIEdrEFmGDHIeKBaVV2vqg+rqgIXAver6ir362uBd93vIqkFPlHVB4HmwOMe65+kqoNUdVBRUc6n49uNj1as5/j732HWV9/xxxH78fjPhqStcRZOq6YF/HHEfjxz2SFU1dRy+p/f54kPVhL91BiGYRiGkUl8yyQgIucAS4CJIlLsfnwZMMNjkZHA3u7/a4GjM6sw+/j33K8579EPadWkES/+cigXHNw540b+IV3aMPXqwzm4W1tumvI5t7z0BTW11kgzDMMwDD/xpYEmImcAJwCbgLuAJiKyD86kgS1h5Y4VkVHu2y3A2+7/ewNf+qE1G1BV7n9zCaOenseBnffk31cOZd99Wvu2/TYtGvP3iwfzi8O7MHnWSi5/4mMqd9gEWsMwDMPwi4w30ERkCPAUcD7wV+A0YDXQBPgKx5sWYj/gcPf/h4HhIvJHYCDOTNB6j6pyx2sLue/NxZx+QEf+cekQ9mhe4LuOvDzhxhP7MPakPry5YB0X/W02W7ZX+67DMHIVETlbREaLyI0iYgEHDcNIioxPElDV2TgzMiNZgTMJILzsRGCi+//3wOWZ1pdNqCq3vPQFk2et5MJDOjP2pL7k5QUbm+zioV1o07IJo57+lIv/Npu/XTKElk1shqdhRENEzgU6A+8AP8OZxf4z4B4cS4dhGEZC+OZBM2IT6jmbPGslvzi8C384OfjGWYiT+xfzwNkD+WRVBZf8bbYNdxpGFESkE3A3zujAecASdwLUcuA8EckPUp9hGLmFNdCyhD+9tYxJb3/FhYd05oYT9s26qP4n7r8P9589gDkrN3DVPz+huqY2aEmGkTWISB7wC3ZOeuoOhELlbwBaAHsFIM0wjBzFxqqygGc+Wu0kCR5QzNiT+mZd4yzET/YvZsPWKn4/5XOuf2E+d52xf9ZqNQyfuQgnFNCN7vvw3rIdUT4DGkYcx3hMmVvGhGmLKKuoJF+EGlU6WOBaw7AGWtC8v+xbbvj3fA7v0Y4JZ/bPmmFNLy44uDPfbtrO/dOX0Lltc646ukfQkgwjUESkF1CrqkvCHljKgVBiySaAAusjl1XVScAkgEGDBjW4eDZT5pbtkvKpxo27WFZRyfUvzAewRprRYLEhzgD5qnwzVzzxCV3ateDh8w6gID83Tsc1x/ZgxIBi7n5jMa9/vjZoOYYRNCcAXUVkLHAwcCTwBRCK91gCfKCqW6Iu3YCZMG3RLvk4w6msqmHCtEU+KzKM7MF60AJiy/ZqfjF5Dvl5wmMXD6Z1U/9DadQVEWH86fuz4rutjHp6Hp3aNKdvsaWFMhomqnofgIj8Eqe3rBaYDvQQkRtw8g2PDE5h9rKmojKl7w2jPpMbXTb1DFVlzAvzWf7tFh46dyCd2jQPWlLSNC3IZ9KFB7JHswKufPITNm2rClqSYQTNozihg4YBH6rqBao6TlVPV9XPA9aWlRQXNkvpe8Ooz1gDLQCe+GAlL89bw7XDenFot3ZBy6kz7Vs15cFzB/L1hkrGPD/f8nYaDRpV3aGqtaparaoW1TkBRg/vRbOC6NFHmhXkM3p4L58VGUb2YA00n/l0dQW3vvIlR/duzxU/6hZ/gSxncGkbrhvWi6nz1/LEByuDlmMYRg4xYmAH7jitHx3cnrJ8d5JFh8Jm3HFaP5sgYDRozIPmIxu3VvHLJz+hfaum3HtW9s/YTJTLjujK7OXf8cdXFjCg057062h+NMMwEmPEwA7WEDOMKFgPmo/8/sXPWff9Nh4+7wAKmzcOWk7ayMsT7j1rAO1aNuaX//zEcnYahmEYRooE1kATkSYi0k5ECkUkak+eiBwnIjeIyM0iktNRuF+et4aX5q3h6mN6MKBTYdBy0s6eLRpz/zkDWb1hK7e/uiBoOYZhGIaR02S8gSYOl4jIbSIyWUSGul9txgnm+BVwaJTlugMTgDuB94HJmdaaKf63cRs3TfmcAZ0KufLI3PedeTG4tA0jD+/KPz9cxX8XlwctxzAMwzByFj960IYB44DbgbeBN93esMeBw4HOqvp2lOVOB9aoag1OsuFhItLeB71pRVX57fOfsaO6lnvP6k+jHAlGW1dGHdeTHu1b8tvn5rFxq4XeMAzDMIy64Edr4SPgNmAbMB9oihPMsRgYDFwlIsVRlotMNhz6LKd44sNVvL24nBtO3JeuRS2DlpNxmhbkc+9ZA/hu8w5ueclCPxlGfWDK3DKGjp9BlzFTGTp+BlPmlgUtaReyXZ9h1IWMN9BUdb2qPqxOkKwLgfuBr4HFwERgLfC67J51O+FkwyIyR0TmlJdn17DamopKxr+6gMN7tOP8gxpOIuR+HffgqqO7M+XTNbz++f+ClmMYRgqE8mWWVVSi7MyTmS2NoGzXZxh1xbfxNhE5B1iC0yi7CejkNtrWAv2AthGLRCYbDn22C6o6SVUHqeqgoqKijGivC6rKzS9+Tq3CuFP7sXv7s37zy6O602ef1tzy0ud8b1kGDCNniZYvM5vyZGa7PsOoK7400ETkDJyEwpuAu4AqHD8awN7AOmC9iAwQkdvdz6cDxW7PWglQhtPrlhO8Ov9/vLngG35zXM+cTOWUKgX5edxxWj++2bSdu62iNIycxSsfZrbkycx2fYZRV/yYxTkEeAo4H/grcBrO7MxSEbkVOBU4Q1Vr2ZnHDlV9A6eRdhvwO+ACt0zWs3FrFbe89AX9OuzBJUNLg5YTGP07FXLRIaU8/sFKPlm1If4ChmFkHV75MIPOkxnynXklmAtan2GkSsYzCajqbKAgyle/jlL2eeD5sPfXZlBaxrjjtQVs2LqDv18yuN7P2ozHtcN68vrn/+OGF+bz8q8Oo6CBHw/DyDVGD+/F9S/M32UYMeg8mSHfWeTQZoig9RlGOrBfyzQzZ8V6/vXRan52WBf262Apj1o1LeAPp/Rl4f828dd3lgctxzCMJAnPlylkR57MaL6zENmgzzDSgeXiTCM1tcrNL37BPns05ZpjewQtJ2sY3ndvjuuzF/dPX8zJA4p/SIxsGIb/TJlbxoRpi1hTUUlxYTNGD+8VtzETVL5ML61e/jIB3htzdMrrN4xswHrQ0sg/Z6/iy7Xfc+OJ+9K8sbV9w7nlpD6owjhLA2UYgZFLISliaU2HLy6XjoXRMLEGWppYv2UHd09bxCFd23Jiv32ClpN1dNyzOVcc2Y2pn63l/WXfBi3HMBokuRSSIpbW0cN70axg17CYyfrOculYGA0Ta6CliQnTFrF5ezV/OKVvg4t5liiX/6gbHfdsxh9e+pLqmpyYkGsY9YpcCkkRS2s6fHG5dCyMhomNw6WBz76u4F8freLSoV3ouVeroOVkLU0L8rnpxD5c/sTHPPnhKi46tDRoSYZRb0jET1Vc2IyyKA2QbAxJEU9rqr64wuYFbIiSLzgbj4XRMLEetBRRVW59+UvatmjMr21iQFyG992Lw7q34543FvHd5u1ByzGMekGifqp0DA36RSa1TplbxuZt1bt9XpAvWXksjIaJNdBSZNoX/2POyg385rhetG4aLdybEY6IMPbkPmzZUcPEN5cELccw6gWJ+qmyMWSGF5nUOmHaIqpqdw9x26Jxo6w8FkbDxIY4U2BHdS3jX1tIz71actagjkHLyRm6t2/FuUNK+OfsVVw8tJRuRS2DlmQYOU0yfqqgQmbUhUxp9TpeGystb7CRPVgPWgo88cFKVny3letP2LfBZwxIll8f24NmBfnc+drCoKUYRs6TLemYQumXuoyZytDxM+ocsiJd6/EiW46XYcTCWhV1ZOPWKh6YsYTDe7TjyJ5FQcvJOdq1bMIVR3bjjS/XMXv5+qDlGEZOkw3esnTFFfMjPlk2HC/DiIc10OrIQzOXsLGyiuuP39fCatSRS4d2Ye/WTbn91QWoeqU8NgwjHtngLUtXXDE/4pNlw/EyjHhktQdNRJoArYBqYLOq7j7tJgBWfbeVf7y/kjMP7Eif4tZBy8lZmjXO57rhvbju2Xm88tlaTupfHLQkw8hZgvaWpSuumF/xyYI+XoYRj4w30MTpXroY6AaUAI8AHwGXAR2A9sA9qvpFlMU3uxo3ACOAtzOtNxHunLaQ/Dzh2mHWHZ4qpw7swP+9u5w7X1/IsL570aRRfvyFDMMIlGgx1+LFFZsyt4yxL31BhWvE37N5Abec1JcRAzvssr48EWqi9Kh7+cMsn6ZRX/FjiHMYMA64HaeB9SYwCviZqo4B1gAzRCRaY/Fx4HCgs6pmReNs7qoNTP1sLSOP6MperZsGLSfnyc8TbjihN19vqOTxWSuDlmMYRhyiecRGPzsv6gzIUFyxKXPLGP3svB8aZwAbtlYx+rl53DRl/i7ri9Y48/KHWT5Noz7jRwPtI+A2YBswH2gK/Au4x/1+PtAcaBxl2WJgMHCViEQd/xKRkSIyR0TmlJeXp1v7btz9xiLatmjMyCO6ZnxbDYXDexRxRM8iHpi+hI1RnsANw8geonnEqmqVKGHFfogr5hV3rKpGeerD1butDyBfJK4/zPJpGvWZjDfQVHW9qj6sjgv8QuB+VV2pqo+7w5/nAzer6tbw5UQkD1gMTATWAq9LFDe+qk5S1UGqOqioKLOzKd9b+i3vLf2OXx7VnRZNstq+l3OM+XFvvt9WzaR3lgUtxTCSRkQGich1InKPiHR3PztbREaLyI0iUm+C/SXjBQv1qsVaJlqPGUCtKsvHn8h7Y472HLK0fJpGfca3WZwicg6wBJgY1hs2CngCZ4gzsgK7FOjkNuzWAv2Atn7pjURVuWvaIor3aMq5B5UEJaPe0qe4NSf1L+axd1dQvslSQBm5g4gU4Ngx/gR8CzwvIocBPwPuBtaxc8QgMNIVWyyZWGF5IkyZWxZzmXyPWfCJbKewefTsLZmIZ5bp2GyGEYkvDTQROQM4AdgE3AU0EZFrgC44szTvAraLyAARud1dbAs7JwXsjVPJBRYw6z9frmPe6gquObYnTQvMyJ4JfnNcT3bU1PLwzKVBSzGMZGgO9MKpzzbhTHw6D1jiPmAuB84TkcAqjnR6taLFECvIEwryd29o1ahy/QvzOap3EQV5u39fkC+cc1CnOsUk8zOfpnndjCDIeANNRIYAT+EMZf4VOA3HV3YfcBXwKNBLVatwZnoOcxd9FigVkVuBU4EzVLU203qjUVOr3PPGYrq2a8FpB9jsoEzRpV0LzhrUkX9+uIqvN2yNv4BhZAGquhHoCnwJHIBTd3UHQuNsG4AWwF6BCCS9Xq1oMcQmnNmfCWf0j9obVllVw8yF5Uw4sz+FzXb2eO3ZvIAJZ/TnthH96hSTzM98muZ1M4Ig40YqVZ0NROuHfiZK2eeB593/q4FfZ1ZdYrw8bw2L1m3ioXMHWkqnDPOro3vw/CdlPDB9CXed0T9oOYaRKCuBMcAhOA+hD4d9t8P9u1sPmoiMBEYClJRkzjqRbq9WKIZYKHTGNU9/GrN8WUUlo57+lOLCZow9ue9uDahkY5JNmVtGWYL5NGOF90gU87oZQWCtjTjsqK7l3v8sps8+rTlhv32CllPvKS5sxgUHd+a5j79mWfnmoOUYRlxEpDHQQVXvAO4AZgMbgZARqgmgRLFo+DXJKRO5J6OFzohFuoYGQ8ONXoTvU6zwHslosNydRhBYAy0Oz8xZzar1Wxn9417kRfFQGOnnyiO70awgn3v/szhoKYaRCD8H5rj/fw+0BN7CCRMEToDuD1R1i//SHDKRe9JriDEeqQ4NRhtuDBG5T7HCeySjwXJ3GkFgDbQYbKuq4YHpSxhcuqclRPeRti2b8LPDujD1s7V8XrYxaDmGEY8XgWdE5Cac4crfAg8Bm0TkBhz/7cgA9WUk92Qqw3uZWjZyn2KVTUaD5e40gsCCecVg8qwVfLNpOw+de4AlRPeZnx/RlckfrOSeNxbxt0uGBC3HMDxR1TLg6ihfXeC3llikO/dkcWEzTx8YOI0YIGqZVIYGvbbbobDZbvvnlX6qLhoSOX7Jpp2yNFVGLKwHzYNN26r481vL+FHPIoZ0aRO0nAZH66YFjDyiKzMXlTN31Yag5RiGEcHo4b2ihs6AnaEuMjE0mOg6vcJwhOtLJ8mG4rDQHUY8rIHmweRZK9mwtYprh/UMWkqD5cJDStmzeQH3T18StBTDMCIYMbDDbqEzYGf4jFCPU7qHBhNdp5f/TIQf9KWTZENxWOgOIx42xCmFyk0AACAASURBVBmFTduqePSdrzimd3v271gYtJwGS8smjRh5RDfufH0hn6zawAElewYtyTCMMBIZ9kv30Gqi6/T0mCkZGUZMNhSHhe4w4mENtChMnrWSiq1V/PrYHkFLafBceEhnHn3nKya+uYTJl5oXzTAyQTpihWWKuvq0vPxnsbxnoW2VVVSSL0KNKh0S3KaXN664sFnUfYhV3shu/PIO2hBnBNZ7ll20aNKIy47oytuLy/l4pXnRDCPdpCtWWCaoq0+rLmmgwrcFO5O4J7pNL2/cUb2Lou7DUb2LLHRHDuKnd9AaaBFY71n2ccEhnWnbojET37S4aIaRbtIVKywT1NWnVZc0ULHiqyWyTS9v3MyF5VH3YebCcgvdkYP46R20Ic4wrPcsO2neuBGX/agr415dyJwV6xlUarNqDSNdpCtWWCaoq0/L6/vINFDpWGc40bxxozzSYK2pqMyIP8/ILH56B62BFob1nmUv5x/cmUlvf8V9by7myZ8fHLQcw6g3xIoVlifClLllKTciIj07R/UuYubC8rgeHi+flgKlY6Z6+sPq4u+KF9ctT4QuY6buojcRL5J5zXKXoL2DgQ1xishxInKDiNwsInvVtUy6CPWeHbuv9Z5lI80bN+LyH3XjvaXfMXv5bikNDcOoA7FihYHjw0pX7sxwz84TH6xKyMMTzdcVjteydYm/Fm9bNaq76L1pyvyEvEiWJio38fKa+ekdzHgDTRwuEZHbRGSyiAwVke7ABOBO4H1gcpTl4pZJJ/94f4XTe3aMxT3LVs47qDPtWjbhPsvRaRhpIZF8mpnMnRlvGyFfV36MTC7Rlq1L/LXwZYAfthlt25VVNTz14eqEvEiWJio38fKa+ekd9GOIcxgwDugKnAe8CYwF1qhqjYgsB4aJSHtV/SZsudMTKJMWnN6z5Ry7b3v6ddwj3as30kSzxvlccWQ3/vjKl3zw1Xcc3LVt0JIMI6dJ1DfjR95Nr3IjBnbw9HHFWrYu/q5oy3QZMzVq2dAsz0xpMYIlltfMr/PpxxDnR8BtwDZgPtAU6A2E9j4UO6F7xHLdEyiDiIwUkTkiMqe8vLxOAv/x/go2VlrvWS5w3kElFLWyXjTDSIYpc8sYOn4GXcZMZej4GT8MwyXqm9mjWYHnOuKR6Dbi+cMSWbauGuuiy6tXL+TbM3Ibr/Pup3cw4w00VV2vqg+rqgIXAvfjeDxD7HD/Rg7+5ydQBlWdpKqDVHVQUVFR0vqs9yy3aFqQzxU/6saHy9cza9l3QcsxjKwnVtymeL6rEN9vc+Ki1SX2UyLbSMUfFlo2U/GpvDxk5xzUKaqmdPj2jODJBu+gb5MEROQcYAkwEagBQs3QJu7fyO6v8gTKpIz1nuUe5x5UQvtWTbjP4qIZRlxixW2K5o9q0Xj3RketOnHRoq0jHtG2cf7BJWnxh4Uvm6n4VF4esttG9PP0x1lOzdwnG7yDvoTZEJEzgBOAt4C7gOnAYBERoAQoAxaLyADgTFW90S0zPLJMOnVZ71lu0rTA8aL94eUvmbXsOw7pZl40wwgRmbbJi5DHJtJP4+W5irUOLx3JpsPxWib8c6/QGl5ayioq6xwqJFLPfT8dsMt6YvnjoulJZP8ymTrISI6gvYMZb6CJyBDgKXdb5+P0np0L9MTxpnUDLlDVWhHphjOp4EZVfUNEhkeWSac26z3LXc4ZUsKf31rGxDcXc0i3Q4KWYxhZQShtU7yZmRDbYxMrHlgi6wgNN4Z6tELDjeCdqNxrmTkr1/P8x2Vx1xVLd7xtp7IPicbFSnX/jIaHHx602apaoKrivhqparWqXquqN6rq2ao60y37vKoODlt2tzLpwnrPcptQL5p50QxjJ4mEzYDYXppEfWmx1lGX4UavZRINZxFLd12GHBPdh0S9Sqnun9HwaLC5OK33LPc5Z4jjRbMcnYbhEC+kRSJemki/lxex1lGXdDhe3yUaziKkO9n1J1vea7vxvEqp7p/R8EhqiFNEOgF9gf2AfkBfVR2UCWGZxHrP6gfmRat/vDp/Ldurazilfwfy8ryDk2YKEckHfqeq43zfeJIkk4YGnIbDe2OOTmjdIe/N0PEzPNcX6uFJJkVTKqmWElnXiIEd+MPLX0RNXRUeiiMRv1esfYi2jnjH1mt9+SJRG2npSrNlZJZM+gfj9qCJyGUi8r6IVOCY9H8OtARewvGS5Ryzl69ny/Zq6z2rB1gvWv1hR3Utt73yJY/PWkmMwPEZRVVrgIPczCenikhWjjLESkNTEKVhW5AvdQoPEGvYMNkUTammWkpkXV6pq0L7n0woDq99OKp3UZ3CeVi4jvpHpkK7hEik8rkeGAUcCLyCE2j2MdcvlpO/isfsuxezrj/Ges/qAeZFqz88+/Fq1mzcxjXH9kSCaqEBqnqKqt4ErAc+EJFrRSSrEvTGSkMz4cz+FDYr+OHzPZsXMOGM/nV6qo833BkvRVMqoTS8iLUuLw9ei8aNkg7F4bUPMxeW18kzZuE66h+ZCu0SIpEhzp+o6ufu/2eKyPHAyyLyd+D+dM+s9IuiVk3iFzJyApvRmfvsqK7lTzOXMbCkkMN7tAtUi4j8BuiAM1rwObAWuFxE5qnqa4GKc/EzDU1ofV3GTCWaWypWiqa6plry2pZAzKFELy0b3ZAjyXrjou1DMmE1Ellf6PNU1msEQ128lskQtwctrHEWev8aMARoA7yXFhWGkQLWi5b7PPfx15RVVAbee+ZyB9AWOEZVh6rqP4EXgR8HK2sniaahSWfqIz9T39R1W/GWK2xeEPN7iH3MpswtI8/j+oynLd65yIbUQkZyZPqc1clfoarbVfX3wEVpUWEYKWJetNxlR3UtD89cysCSQo4IuPfM5VJVvVhV54R9Nhj4OihBkSTi8Uq3P8bP1Dd13Vas5eL50yD2MQt9F83QH09bIuciG1ILGcmR6XOWUqDaXPWgGfUPm9GZu4R6z8ad1i8bes9Q1SejfDY5CC1ehIbJYs0ei5fiKRPbTBd13Vas5YaOnxHTnxZaLpanKPI7cGZhxvPXJXIu/Dy+RnrI9DnzJdWTYfiBedFyjyzsPcsZ4nm8MuGP8TP1TV235bVcPH9arDKxjlmtalydycRUswZZbpHJc2YNNKPeYL1oucfznzi9Z7eful9W9J7lGrFiMHnF3coTocuYqVnRQ+NnDsp4sdlC/rJoQ5ihMsnGdkt02+FYXs7gyZZzkJUxfgyjrpgXLXfYUV3LQzOWMqBTIT/qWRS0nJwjnq/JK65YjWpGYjYlS6ZjSEUSz58Wz1+Wit8o0WX9PibG7mTTOQikgSYiR4rIeBHZV0SKRWRvEWkepVxTEWknIq3dCN+GEROb0Zk7hHrPrjm2h/We1YF4fqnIuFvZFmcr0zGkIokVmy2aFtjVX1aX2G6JbDscv4+JsTvZdA4yPsTpBni8EjgGWKaqI4FS4HfuC6AWOAyYFbZcB3bOmioDjgMWZFqvkfuYFy37aQi9ZyJyHM7sz0bAI6q6Lp3rT8TXFO6P6TJmalLryTSZjiEVjWT9aZH+slT8RoksG8QxMXYlm85BxhtoqloBjBORnkBj9+M84BxgI06Ggu2qOiti0QKceEQvAZ+q6rZMazXqB+ZFy37qk/dMnB24GOgGlACPAOuACTj121HAZGB4OrebbL7LuuTHzCTZoicR71lk+Xj+pLp6mDJxTPz0U2WLdysVvM6BAkPHz/B1nwIZ4lTVx1T1X8AMYG/gbo+i3XCeQK8WEcvLZCRMyIt2/3TzomUboZmb/etP79kwYBxwO/A28CZwOrDGze25HBgmIu3TudFkPVHZFmcrG/QkG9ssEX9SKh6mdB8TP/1U2eTdSoW65p/NBEFPErgWeFc1yt3hDHt+oqoPAs2Bx6OtQERGisgcEZlTXl6eQalGLhHqRfvgK/OiZRsvfPI1X2+oV96zj4DbgG3AfJx8xb2B0GP4Bvdv93RuNFlPVCoeqkyQDXoS8Z7FKx/pT0rFw5TuY+KnnyqbvFupUNf8s5kg6DAbl+H0kEVjJNDC/X8tEDUBm6pOAiYBDBo0KFpDz2ighLxo9083L1q2sKO6lofc3rMj60fvGaq6HngYQEQuBO4HWocV2eH+3e2xXERG4tR1lJSUJLS9VIaRsi3OVtB6EvWexSsf/rlXmWjDZuFEntf7fjog5aFTP/1Ufm4r00Op8XLCllVUMmVuWcav3cB60ERkHxy/xpawz44VkVHu2y04wwXgDIN+6a9CI9exXrTsox72nv2AiJwDLAEmAjVA6BG8ift3ty5+VZ2kqoNUdVBRUfwGa30ZRsoWks2lmMjnXmUEPM9TpoZOcyF/arL4eQ/E0u7HfZfxBpqI5InIWGAAMEBExopII5xK6yucYYEQ+wGHu/8/DAwXkT8CA3FMuIaRFOZFyx52VNfy4Iyl9O+4R73pPQshImcAJwCbgLuAOUCxO4GgBGcmesoXYX0ZRsoWMuHjGz28F9EePRQ8z1Omhk5zIX9qsvh5D8Tyo/lx3/kxi7MWGOu+wlmBMwkgvOxEnKdPVPV74PJM6zPqNzajM3t4+qNVWZVzM12IyBDgKZz69Hyc3rNzgZ443rRuwAVuXZgS2RQCoD6QbC7FRMqPGNiBa57+NOryyZ6/RIZOY537XMifmix+3gMh7cmez3QRtAfNMDKOedGCp3JHDQ/OWMqQ0jb1Luemqs7GCQsUybXp3lY6wzD4ERIhW8MupBIqIxHfXIcYabaieZe8zusezQoYOn4GayoqY4YCibU/uZA/NRmSvQdSvQZDgYy9zk86tuFF0LM4DSPjhHvR3l/2bdByGiSPf7CCbzZt59phPetV75nfpGsYyQ8fT7b65TIdKgNip9mKtp5o5QvyhC07qn/Q4BUK5KjeRVl5nDNFMvdAuq7B0cN7UZC3e721ZUc1N02Zn7Hjbw00o0FwzpAS9m7dlLunLSJ6VBcjU2zeXs2f31rG4T3acVBXG2JOhXSFYfDDx5OtfrlMh8qAnecp0fRa0c5ry6aNqKrZva7KF9nl3M9cWJ6VxzlTJHMPpOsaHDGwAy2b7j7gWFWjPPXh6owdfxviNBoETQvyufqYHtzw7/lMX/ANx/bZK2hJDYbH3l3Ohq1VXDcsmICo9Y10DCP54ePJVr9cpvxekYwY2IFRSXiXIs+rV2quWlWWjz/xh/fJbKO+kOg9kM5rsGJrVdTPo/Vs1nUbkVgDzWgwnDmoI5PeXsbdbyzi6N7tyYvSZW2kl4qtO3j07a8Y1mcv+ncqDFpOvSFVz4sfaZayJZVTtO3H05Uu7V7ryROhy5ipFBc246jeRcxcWL7buUxUg5+erCA8hclsM7JsYfMCNkRpWCV7HmOlA8v3+DzkT0sFG+I0GgwF+XmMOq4nC/+3iZfmrQlaToPgkbe/YvOOan4zrGfQUuoN6fDV+BESIRtSOUUj0VAZ6dAey4sWOndPfLAq6rlMVINfnqwgPIXJbDNa2c3bqinI3/VBPNnzGC8d2DkHdfL0p6V6bKyBZjQoTtq/mH33ac29/1nMjuqUox4YMSjftJ2/v7eCk/YvpvfereMvYCREOnw1fqRZyoZUTnXVlS7tkeuJ5kmLJHQuE9XglycrCE9hMtuMVraqVmnRuFFK5zFeOrDbRvTz9KelemxsiNNoUOTlCaOH9+TSv8/hmTmrOf/gzkFLqrf86a2l7KipZdRx1nuWTtLlq/EjJELQqZy8SERXurSHr8fLVxZJ6FwmqsEPT1YQnsJktulVdmNlFZ/eMiztGsLTgXn501I9NtZAMxocR/Vqz6DOe/LA9CWcfkBHmjWOHinaqDtrKip58oNVnHFAR7q0axF/ASNh0uGPytb4ZLlMIsfU69xFEutcphLDLRFPXOT6QuvymvueSU9hMtd6rLKpXO9+ehYjsSFOo8EhIvz2x735ZtN2/jFrRdBy6iUPTF8CwNXH9ghYSf0jVX9UtsYny2USPaZe8bTCiXUuU43hlognLnx94etKVms6SOZa9yqbapw4Pz2LkVgDzWiQDOnShiN7FfHnt5axsTJ697RRNxav28Qzc1Zz3sEldAh4xl59JFV/VLbGJ8tlEj2mXvG0QsQ7l6nGcEvEExe+Pi//VSJa00Ey17pX2VTjxPnpWYzEhjiNBst1w3rxkwff5dG3v+K6gGeW1SfufG0hLZo04uqjrfcsU6Tij8rW+GS5TDLH1MuvJMB7Y45OeTvxyiTiiQuV9VpXIlrTRTLXerSy6YgT56dnMZzAGmgi0gRoBVQDm1W1OkqZ44DBODofUdV1/qo06jP7ddiDE/ffh/97dzkXHtqZ9q2aBi0p55m17DumL/yGMcf3Zs8WjYOW0yBI1l+Tql/G/Gu7E8vbFZl7M9HjHzrOZRWVP8Ta8oq5lYgfKprPLJ6Wul4r0a4R8Cdpe+T2veKXKTB0/IyUdGT6Xsj4EKeIFIrIDSIyXUQmhX21GSgHvgIOjbJcd2ACcCfwPjA501qNhsd1w3pRVVPLff9ZHLSUnKe2VrnjtQUU79GUiw8tDVpOg6AufrJU/DLmX4tOMrk3Ezn+kd6vUAPDKxZXPD9UaNlkY63V5VqJdo2MfnYeo5+b58t1E7l9r0j/pKjDj3sh4w00Va1Q1XHAaiD8kfpx4HCgs6q+HWXR04E1qloDLAeGiUj7TOs1GhZd2rXggkM68/RHq1mw9vug5eQ0r8xfy2dfb+TaYb1oGuUHwkg/dfGTpeKXMf9adJLJvZnI8Y/l/YLd83HG8kPF0hRPS12uFa94ZJF5RTN13cSKWxaNuurw414I0oNWjDN8ebiI/ENVI0O7dwdCfasbwj77JryQiIwERgKUlJRkTq1Rb/n1MT144ZMyxr26gMmXDkESCCZp7Mr26hruen0h++7TmlMb+HCXn9TVT5aIXyba8I3513YS7fjUJpiXMd7xj3c8I/NxRpKIzyzUOxcqF9qXUAMjvJEWS2u0odhEycR14zXjtFYVgajhQuqiw497IZBZnCKSBywGJgJrgddl91/F8EfwHVE+A0BVJ6nqIFUdVFRUlBG9Rv2msHljrj6mB+8s+Za3FpcHLScneXzWSr7eUMkNJ/S2HKc+4uUFSjX+ktfwTWHz6PkFg86v6TeZPj7xyiezPq+ygrMf6Ur/BLGHE5PRVlemzC3Dq/YpLmyW1vslU/deOEGF2bgU6KSqitNA6we0jShTDoT2tEnYZ4aRdi44uDOlbZtz+9QFVNdYCqhk2Li1igdnLOXwHu04vIc9JPlJpuIveQ3fqJKV+TX9JtPHx8tHVpf1jR7eK2qjRXH2I93pnxIlE9eNV0BdwTkO6bxf/Mg1G1QDbQsQ8p3tDawD1ovIABG53f18OlDs9qyVAGU4vW6GkXYaN8pjzPH7svSbzTz10eqg5eQUE6cvZtO2Kq4/ft+gpTQ4MhV/KVbanGzMr+k3mT4+4ecVdvqn6rK+EQM7eGYBWFNRmZH0T4mQievGS4+yc6g2XdevH7lmM+5Bc4czbwYGuO/HArcB94jIrcD+wBmqWisi3YBhwI2q+oaIDHfLdgMuUFXr2jAyxvC+ezGkSxsm/mcxpwwopnXT6MMVxk4Wr9vE5FkrOXtICX2KLSF6EIR7hEJ+oFFPf5rStP9Y4RWyNb+mn/hxfOrqE4y2TIc44TK8fFulrn9tz+YF3HJS3zqnropGpNctHXjp6eCmexr70hdUuIHJmxfksXVHNaOe/pSxL32BCGzYWvWDj65DAvdPpu8FP2Zx1qrqWFUd4L7Gqmq1qv5aVW9W1RGq+q5b9nlVHRy27LWqeqOqnq2qMzOt1WjYiAg3nbgv323ZwcMzlwYtJ+tRVW59+UtaNM7numENa4grG0nntH8/hm9ymWw4Psmc71h6vYZTw3vdNmytYvRz86KmrvIaio1HJsJSxEr3NPrZeT80zgC2VtWyYWsVClRUVrHBDSAc8tFlQwgZS/VkGGHs37GQ0w7owGPvLuer8s1By8lq3vhyHe8u/ZbfHNeTNhaUNnDSOe3fj+GbXCYbjk8y5zuW3kTCcgBU1WjU1FVeQ7Hnh6V6S3eICy9ipXuqqk1uAkMm9CWLaJKzLrKZQYMG6Zw5c4KWYeQ432zaxjF3/5cBJYUWdsODbVU1HHfff2lWkM+rVx9Oo/zgnvVE5GNVHRSYgDSRav3VZcxUT4N0rJAMRm6SqfPttd5U1x3k9Rlrn+KRaX2x6i/LxWkYEbRv1ZRRx/Xk1le+ZNoX/+PH++0TtKSs4//eXc7q9ZU8+fODAm2cGTtJNYUTWBqnXMLrfCuOd8zLNxYi8lwf1buImQvLYzZkvNJRJXK9xEpBVTpm6m7prBLxgCWqIRWvXDpSQtUVq1kNIwoXHtKZ3nu34o+vLGDrjt3SxDZoVq/fyoMzlvDjvnsztHu7oOUYLqn6oiyNU24xengvCmLEHPTyjUH0c/3EB6tiNmIK8sUzHVUi10usFFTR/iZy/SWqId6xikdQ94I10AwjCo3y87j1lP0oq6i0CQNhqCo3v/g5+SLccnKfoOUYYaTqi7I0TrnFiIEdaNk09iBYNN8YJB+/bM/mBUw4o3/cdFSxrpdEvW6Jri8ZDSMGdmDCmf0pbLZzZn7zgjySabMFcS/YEKdheDCkSxtOG9iBR99ezmkHdKRbUcugJQXOa5//j5mLyrnpxH3ZZ4+GFT0+F0hl2r+lcco9KrZWxS0T7fwlek5j+a/qcr0kkoIqmfUloyHavZGohkS0ZAJroBlGDMac0Js3F6zj+ufn86+RBzfoNEabtlXxh5e/oM8+rbn40NKg5RhxSNZPlg4Pm+EviXirop2/RD1ZezQrYOj4GVGvoVSulylzy/BMjBlByAMW8siFa0n1mvVa3iunqN/3gg1xGkYM2rdqyk0n9mH2ivX8c/aqoOUEyj1vLOabTdsZd1o/mxiQ5dTFT5YNsb2M5IgXhyzSNxa+XDxPVkGesGVHtec1VNfrZcrcMkY/O49kAkiEe+TCtRzVuyila9ZrH845qFNW3AtWyxpGHM4c1JGh3dsy/rWFrN3YMId7Pl1dweRZK7jg4M4M6FQYtBwjDnXxk2VDbC8jOSLjkIU3uaL5xsKXi+Vf61DYjJZNG1FVs2srKvwaquv1MmHaojrFJIuksqqGmQvLU7pmvfbhthH9suJesDhohpEAq77byrCJ/+Ww7u149MJBDSo22raqGk584B0qd9Tw+qgjsi4FVrbGQROR44DBOFaSR1R1Xazy6ay/LCaaEY9410gQcdaSpT5czxYHzTBSpKRtc649rhe3v7qAlz9by8n9i4OW5Bv3/Wcxy8q38PjPhmRd4ywbEJFC4ErgGGCZqo4Uke7ABOBA4ChgMjDcL03mJzPA24c4ZW4ZeXF8VrF8agNvfcMzb2X4NvdoVoCIM5khnm8MvL1fngg/DLmG59mMFwMuHfgRMzCrhzhFpImItBORQhGxxqQRKJcMLaV/xz245cXP+eb7bUHL8YWPV65n0jtfce5BJRzeoyhoOVmJqlao6jhgNRDKeXU6sEZVa4DlwDARae+XJvOTGV4+xJumzOf6F+ZHbQiFXyOxcnR65a0MrTu0zVCOy0jfWDT/W0G+RPV+xUIVfvPMp/zmmU93ybMZKwZcOvArZmDGG2hu4+oGEZkuIpPczxqLyK9EZLyIPCYifT0W3wyUA18Bh2Zaq2HEolF+HvecNYDKqhque+4z6pM9IBqVO2q47tnPKN6jGTecsG/QcnKN7kCom2BD2Ge+YH4yw8uH+NSHq6PGQMsX2eUaCV1DicQri7fu8DIzF5bvFpMs5JcL936FNIFz/bZoHL3hVqvOKxKvGHDpwK+YgRnvlVLVCmCciPRk59PlJcDPVHWAiNwGzBCRDqoaGbL9ceAxYJ6qbsq0VsOIR/f2LbnxhH35/YtfMHnWSi6qx+Em7nx9Icu/3cI/f34QLZtYB3aShP+a7IjyGQAiMhIYCVBSUpJWAanERDOyiylzy5IewvMaRvQaQqxV3W19IwZ2YNTTnyasM5HhybKKSiZMW8TGyqrdhkbDQ3pM/OmAXfQkG7MM4sctq+swpV8xA4Ma4nwduMf9fz7QnJ2Nt3CKcUy2V4lIwzH9GFnN+Qd35qheRYx7dQFL1tXP54b/fLmOv7+/gosPLeVQS+dUF8qBkOGrSdhnu6Cqk1R1kKoOKiqyIWRjd0JhKZIZwpsytwyvfi+vHjEvf2IyvsVEetsEPIddYw0Z1sU/GWuZVIYp03GsEiGQBpqqrlTVx8WZCnc+cLOqbg0vIyJ5wGJgIrAWeF2iTJ0TkZEiMkdE5pSX71b/GUbaERHuOqM/LZs04pqnP2VHdW3QktJKWUUl1z07j/06tOb6E3oHLSdXmQ4Uu3VWCVCGU58ZRlJ4haWINYQ3YdoizxmYycb4ihdrLXwd8Txk0WLTeg2NRg4ZJptP0ysGXIhUhin98ngGPUlgFPAEzhBnZB6dS4FO6hh91gL9gLaRK7AnUCMIilo14c7T9+eLNd9z29Qvg5aTNqpravn1U3OprqnlwXMOoEmjxA27DRURyRORscAAYID7/wycRtptwO+AC1S1frXkDV9IV6ojcBpHycb4ivQzFjYrYM/mjn8s3CMWLX5YqGxoO14DoF5Do+H7ES2fphexYsBFW3cin4fjl8czMGOJiFwDdAG+B+4CfiIiA4AzVfVGYAvwtlt8b2AdsD4IrYYRjWP77MVlR3Tlkbe/YkCnQk47oGPQklLmvjcXM2flBu4/ewBd2rUIWk5O4Da8xrqvcK71XYxR74gVliJPhC5jpiachilkvo/lT/TyZSUT/NWr7NDxM+qUWilS09iT+zJh2iLPfXxvzNFxdXodo1BqqXh+ND88nn7M4tzt6VJEzgLuA64CHgV6qWoV0A0Y5i76LFAqIrcCpwJn2BOokW2MHt6Lg7u24YZ/z+fLNd8HLSclpi9Yx5/eWsZPB3XilAFmLjeMtS1ZpAAACzNJREFUbCDW0F6NatrTMGUyfERdUit5acpEmqcQmQqbkSwZb6Cpaq2qjlXVAe5rrKo+o6oS9ip1yz6vqoPd/6tV9deqerOqjlDVdzOt1TCSpVF+Hg+ecwCFzRpz+RMfs3FrVfyFspCl32zi1//6lL7FrRl7slfUG8Mw/Cba0F609lq60jBlMnxEXVIreWlKZ5qnaGQibEayWKonw0gDH6/cwNmTZnFIt3Y8dtGgnEomvrGyihEPv8embVW8eNVhnhVWtpKtqZ6SxeovI1H8TsMUZEolPzQFud+x6q/c+RUxjCzmwM578sdT9uPtxeX8/sUvciaI7fbqGi5//GO+3rCVP59/YM41zgyjIZKpMA9+hY9IBj80ZeN+gzXQDCNtnD2khCuO7MZTs1fxl/9+FbScuNTWKqOf/YxZX33HXWfsz+DSNkFLMgwjATIV5iEbU4T5oSkb9xssWbphpJXRw3qxev1W7nx9IW1aFPDTwemNDp8uVJU7XlvAS/PW8Nsf9+LUgbk/A9UwGgohn1W6k3Vnar3Zrikb9xvMg2YYaWd7dQ2/mPwx7ywpZ+JPB2TdjEhV5e43FvHwzGVcdEhnxp7clygxoHMG86AZhpGrxKq/rAfNMNJMk0b5PHL+gVz0t9n85pl55OcJP9k/OzKVqSoT31zCwzOXcfbgTtxyUm43zgzD8J+65rD0m1zR6YV50AwjAzRrnM9jFw/mgJJCfvXUXJ6avSpoSdTWKrdPXcD905dwxoEdGXdqP/KSSJ1iGIaR6Vhp6SJXdMbCGmiGkSFaNmnE5EsP4kc9i7j+hfk8PHNpYLM7q2pqufbZefz13eVcfGgpd52+vzXODMNImkzHSksXuaIzFtZAM4wM0qxxPpMuGMQpA4qZMG0R1zz9KdsiKo1M882mbZz36If8e24Zo4f34paT+ljjzDCMOpFKDks/yRWdsTAPmmFkmMaN8pj40wH03KsVd7+xiGXlm3ng7IF0LWqZ8W1/vHI9Vz75CRsrq7j/7OybsGAYRm7hlcMy6JhhkeSKzlhYD5ph+ICI8MujuvPoBYP4ekMlJz7wLk9+uDJjQ57bqmoY/9pCzvzLLJo0yuffVw61xplhGCmTrTHDIskVnbGwBpph+MixffZi2jVHMKh0T2789+ec+ZdZfF62MW3rV1VmLvqGnzz4Ln/57zJ+OrgTU68+jH33aZ22bRiG0XCpa55Pv8kVnbEILA6aiBwHDMYZZn1EVdfVpUw4FkfIyBVqa5VnP17NXa8vYv3WHZzcv5jLjuhGn+K6NaRUlVnLvuPht5by3tLvKG3bnLEn9+XIXu3TrDz7sDhohlG/yPXwGMkQaBw0ESkErgSOAZap6kgR6Q5MAA4EjgImA8MjlotbxjBylbw84aeDSzi+3z48PHMpT8xayYufruGQrm05qX8xw/vuRduWTeKuZ8W3W3hzwTqembOaxes206ZFY8ae1IdzD+pM40bWQW4YRm4RCo8RmoEZCo8B1NtGmhcZb6CpagUwTkR6Ao3dj08H1qhqjYgsB4aJSHtV/SZs0UTKGEZO07ppAdcfvy9XHtmdJz9cybNzvuaGf8/nxinz6dquBf07FtKpTXPatmxM00b5bNlRzYatVSxZt4kv137Pyu+2ArBfh9ZMOGN/TupfTNMI34VhGEauECs8hjXQ/KE7EJpesSHss2+SLIOIjARGApSUZGfeQ8OIxx7NCrjyyO5c8aNuLFi7iTcXrGPe6greWfot5Zu271JWBErbtmDfvVtz8aGlHNN7L0raNg9IuWEYRvqoD+Ex0kVQDbTwR/wdUT5LtAyqOgmYBI6HI10CDSMIRIQ+xa138aJV1dRSsbWK7dU1tGzSiOaNG9nwpWEY9ZL6EB4jXQRVy5cDoaPdJOyzZMsYRr2nID+PolZN6LhncwqbN7bGmWEY9Zb6EB4jXQRV008HisXJ0lwClAGLRWSAiNweq0wgag3DMAzDyDj1ITxGuvBjFmcecDMwwH0/FrgNpwF2G9ANuEBVa0WkGzAMuFFV3xCR4ZFlMq3XMAzDMIzgGDGwQ4NskEXixyzOWmCs+wrn2ihlnweeD3u/WxnDMAzDMIz6jplZDMMwDMMwsgxroBmGYRiGYWQZgaV6ygQiUg6sTLB4O+DbDMrJJLmq3XT7S0PR3VlVizIlxi+SrL+g4ZzfbMF0+0tD0e1Zf9WrBloyiMicXM3fl6vaTbe/mO76Ta4eJ9PtL6bbX9Kp24Y4DcMwDMMwsgxroBmGYRiGYWQZDbmBNiloASmQq9pNt7+Y7vpNrh4n0+0vpttf0qa7wXrQDMMwDMMwspWG3INmGIZhGIaRlWQ8k0A2ICKDgBMAAZ5U1aXJfB8U2aorHvVVt4g0BX7lfv+lqr7iv8rdSUD3gcDhQDPgNVX91H+VyZOr11G6sfrLX+qrbqu//Cfla0lV6/ULaAMsAZoD3YHPgYJEv89i3QJcgpOrdDIwNGjNyR5P4Fzg+qA1J3G8pwIXAL2AhUFrTlB3AfAP9/8mwJygNbtaCoEbcHLyTkrlOqrPL6u/skt3RFmrvzKvOyvrL1dPxuuwhjDEeRxQpapbgRVAX2BgEt8HRTxdw4BxwO3A28CbIrKX3yKjkNDxFJFOwN04N102EE/3QW6Zp4ANwHN+C/Qgnu5WwFki8nOgFqjwXWEUVLVCVccBq4HGUYpk633pN1Z/+YvVX/6Sk/UX+FOHNYQGWnegEkBVq4FN7meJfh8U8XR9hPP0uQ2YDzQlOyqLuMdTRPKAXwAzfFfnTTzd+wPrgTOBy4Fn/BboQUzdqroeeBF4FKeSmOC/xDqRrfel31j95S9Wf/lLfa2/IA33ZkNooOVHvN8R8Vm874Mipi5VXa+qD6vTl3ohcL+qrvJToAeJHM+LgMeBal8UJUY83e2AvYAPgbdwnvhb+iMtJjF1i0hzYAvwZxz9D4lINvwQxiNb70u/sfrLX6z+8pf6Wn9BGu7NhtBAK8cxF4Zo4n6W6PdBkZAuETkHZ5x7oogU+6QtFjF1i0gvoFZVl/gtLA7xjvc6nCfQ5e7/RUBv39R5E0/3qcA3qnolcBiO7n7+yasz2Xpf+o3VX/5i9Ze/1Nf6C9JwbzaEBtpbQFsRKRCRIqAG+FJE/uR2VUf7fnZganfyFrF1IyJn4MwQ2QTcRXYMEbxFbN0nAF1FZCxwMHCkiBwWmNqdvEVs3W/gmD3zgD1wnp6z4Yn/LWLrbo57XajqB275NQFpjYmIlOTAfek3b2H1l5+8hdVffvIW9aT+gvTXYQ0iUK2IjAK6AC2BacBK4Emgn6pujfxeVZ8OTGwYsXQD+wHvsTNUSg3Q1B3rDpR4x9st80vgOpwnut+r6nsByf2BBK6Ti4DBQFvgP6r6WGBiw4hznSjO8MBWHIPtElX9W0BSf8CtwG4GRrgfTfn/9u7YtoooCgLojIUs2nEDFICLcIhICV2FA8dIjmiCFhy5BEsOoIhH8BEZXkf778rnxBu8YHc0T3elm+RnkocM/y73Jr/2Jb/2dcT8SvbJsHdR0Jin7WVOt7iL5N9PlADjyS/2oKABAAzzHv5BAwA4FAUNAGAYBQ0AYBgFDQBgGAUNAGAYBQ0AYBgFDQBgGAWNw2p70/Z329u2V22f/u72AxhNfrHlw/YjMNNa66HtryT3Oa0v+bzWej7zsQA2yS+22CTA4bX9kdMC4OvlhQYORH7xP0acHFrbqySPOd1Av5z5OABvJr94jREnh9X2a5JvST4leUnyve3HtdbdeU8G8Dr5xRYjTgCAYYw4AQCGUdAAAIZR0AAAhlHQAACGUdAAAIZR0AAAhlHQAACGUdAAAIb5A2Pi807QkhLiAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 720x216 with 2 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"# Specify a 'true' latent function lambda\n",
"scale = lambda x: 9 * np.sin(2 * np.pi * x) + 1 + 18.\n",
"\n",
"# Generate synthetic data\n",
"# here we generate some synthetic samples\n",
"NSamp = 200\n",
"\n",
"X = np.linspace(0, 1, NSamp)\n",
"\n",
"fig, (lambdaf, samples) = plt.subplots(1, 2, figsize=(10, 3))\n",
"\n",
"lambdaf.plot(X,scale(X))\n",
"lambdaf.set_xlabel('x')\n",
"lambdaf.set_ylabel('$\\lambda$')\n",
"lambdaf.set_title('Latent function')\n",
"\n",
"Y = np.zeros_like(X)\n",
"for i,x in enumerate(X):\n",
" Y[i] = np.random.poisson(scale(x))\n",
"samples.scatter(X,Y)\n",
"samples.set_xlabel('x')\n",
"samples.set_ylabel('y')\n",
"samples.set_title('Samples from poiss. distrib.')"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"train_x = torch.tensor(X).float()\n",
"train_y = torch.tensor(Y).float()"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"class PVGPRegressionModel(gpytorch.models.ApproximateGP):\n",
" def __init__(self, num_inducing=8, name_prefix=\"test\"):\n",
" \n",
" self.name_prefix = name_prefix\n",
" \n",
" # Define all the variational stuff\n",
" inducing_points = torch.linspace(0., 1., num_inducing)\n",
" variational_strategy = gpytorch.variational.VariationalStrategy(\n",
" self, inducing_points,\n",
" gpytorch.variational.NaturalVariationalDistribution(num_inducing_points=num_inducing),\n",
" )\n",
"\n",
" # Standard initializtation\n",
" super().__init__(variational_strategy)\n",
" \n",
" # Mean, covar, likelihood\n",
" self.mean_module = gpytorch.means.ConstantMean()\n",
" self.covar_module = gpytorch.kernels.ScaleKernel(gpytorch.kernels.RBFKernel())\n",
"\n",
" def forward(self, x):\n",
" \n",
" mean = self.mean_module(x)\n",
" covar = self.covar_module(x)\n",
" return gpytorch.distributions.MultivariateNormal(mean, covar)\n",
"\n",
" def guide(self, x, y):\n",
"\n",
" # Get q(f) - variational (guide) distribution of latent function\n",
" function_dist = self.pyro_guide(x)\n",
" \n",
" # Use a plate here to mark conditional independencies\n",
" with pyro.plate(self.name_prefix + \".data_plate\", dim=-1):\n",
" \n",
" # Sample from latent function distribution\n",
" pyro.sample(self.name_prefix + \".f(x)\", function_dist)\n",
"\n",
" def model(self, x, y):\n",
" pyro.module(self.name_prefix + \".gp\", self)\n",
" \n",
" # Get p(f) - prior distribution of latent function\n",
" function_dist = self.pyro_model(x)\n",
" \n",
" # Use a plate here to mark conditional independencies\n",
" with pyro.plate(self.name_prefix + \".data_plate\", dim=-1):\n",
" \n",
" # Sample from latent function distribution\n",
" function_samples = pyro.sample(self.name_prefix + \".f(x)\", function_dist)\n",
" \n",
" # Use the exponential link function to convert GP samples into scale samples\n",
" scale_samples = function_samples.exp()\n",
"\n",
" # Sample from observed distribution\n",
" return pyro.sample(\n",
" self.name_prefix + \".y\",\n",
" pyro.distributions.Poisson(scale_samples), \n",
" obs=y\n",
" )\n"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "49236c69d8e84835b4a855b17ce451df",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"HBox(children=(IntProgress(value=0, max=3000), HTML(value='')))"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"ename": "RuntimeError",
"evalue": "Non-negative-definite natural covariance. You probably updated it using an optimizer other than gpytorch.optim.NGD (such as Adam). This is not supported.\n Trace Shapes: \n Param Sites: \n Sample Sites: \nnum_particles_vectorized dist | \n value 32 | \n .u dist 32 1 | 20\n value 32 1 | 20",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mRuntimeError\u001b[0m Traceback (most recent call last)",
"\u001b[0;32m<ipython-input-6-42de373c70c1>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m 26\u001b[0m \u001b[0mhyperparameter_optimizer\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mzero_grad\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 27\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 28\u001b[0;31m \u001b[0mloss\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0melbo\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdifferentiable_loss\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmodel\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmodel\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmodel\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mguide\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtrain_x\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtrain_y\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 29\u001b[0m \u001b[0mloss\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mbackward\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 30\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/pyro/infer/trace_elbo.py\u001b[0m in \u001b[0;36mdifferentiable_loss\u001b[0;34m(self, model, guide, *args, **kwargs)\u001b[0m\n\u001b[1;32m 107\u001b[0m \u001b[0mloss\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;36m0.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 108\u001b[0m \u001b[0msurrogate_loss\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;36m0.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 109\u001b[0;31m \u001b[0;32mfor\u001b[0m \u001b[0mmodel_trace\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mguide_trace\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_get_traces\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmodel\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mguide\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 110\u001b[0m \u001b[0mloss_particle\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0msurrogate_loss_particle\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_differentiable_loss_particle\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmodel_trace\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mguide_trace\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 111\u001b[0m \u001b[0msurrogate_loss\u001b[0m \u001b[0;34m+=\u001b[0m \u001b[0msurrogate_loss_particle\u001b[0m \u001b[0;34m/\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mnum_particles\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/pyro/infer/elbo.py\u001b[0m in \u001b[0;36m_get_traces\u001b[0;34m(self, model, guide, args, kwargs)\u001b[0m\n\u001b[1;32m 165\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmax_plate_nesting\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0mfloat\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'inf'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 166\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_guess_max_plate_nesting\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmodel\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mguide\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 167\u001b[0;31m \u001b[0;32myield\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_get_vectorized_trace\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmodel\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mguide\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 168\u001b[0m \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 169\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0mi\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mrange\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mnum_particles\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/pyro/infer/elbo.py\u001b[0m in \u001b[0;36m_get_vectorized_trace\u001b[0;34m(self, model, guide, args, kwargs)\u001b[0m\n\u001b[1;32m 147\u001b[0m return self._get_trace(self._vectorized_num_particles(model),\n\u001b[1;32m 148\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_vectorized_num_particles\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mguide\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 149\u001b[0;31m args, kwargs)\n\u001b[0m\u001b[1;32m 150\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 151\u001b[0m \u001b[0;34m@\u001b[0m\u001b[0mabstractmethod\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/pyro/infer/trace_elbo.py\u001b[0m in \u001b[0;36m_get_trace\u001b[0;34m(self, model, guide, args, kwargs)\u001b[0m\n\u001b[1;32m 51\u001b[0m \"\"\"\n\u001b[1;32m 52\u001b[0m model_trace, guide_trace = get_importance_trace(\n\u001b[0;32m---> 53\u001b[0;31m \"flat\", self.max_plate_nesting, model, guide, args, kwargs)\n\u001b[0m\u001b[1;32m 54\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mis_validation_enabled\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 55\u001b[0m \u001b[0mcheck_if_enumerated\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mguide_trace\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/pyro/infer/enum.py\u001b[0m in \u001b[0;36mget_importance_trace\u001b[0;34m(graph_type, max_plate_nesting, model, guide, args, kwargs, detach)\u001b[0m\n\u001b[1;32m 42\u001b[0m \u001b[0;32mand\u001b[0m \u001b[0mthe\u001b[0m \u001b[0mmodel\u001b[0m \u001b[0mthat\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0mrun\u001b[0m \u001b[0magainst\u001b[0m \u001b[0mit\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 43\u001b[0m \"\"\"\n\u001b[0;32m---> 44\u001b[0;31m \u001b[0mguide_trace\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mpoutine\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtrace\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mguide\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mgraph_type\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mgraph_type\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_trace\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 45\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mdetach\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 46\u001b[0m \u001b[0mguide_trace\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdetach_\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/pyro/poutine/trace_messenger.py\u001b[0m in \u001b[0;36mget_trace\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m 185\u001b[0m \u001b[0mCalls\u001b[0m \u001b[0mthis\u001b[0m \u001b[0mpoutine\u001b[0m \u001b[0;32mand\u001b[0m \u001b[0mreturns\u001b[0m \u001b[0mits\u001b[0m \u001b[0mtrace\u001b[0m \u001b[0minstead\u001b[0m \u001b[0mof\u001b[0m \u001b[0mthe\u001b[0m \u001b[0mfunction\u001b[0m\u001b[0;31m'\u001b[0m\u001b[0ms\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mvalue\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 186\u001b[0m \"\"\"\n\u001b[0;32m--> 187\u001b[0;31m \u001b[0mself\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 188\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmsngr\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_trace\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/pyro/poutine/trace_messenger.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m 169\u001b[0m \u001b[0mexc\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mexc_type\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34mu\"{}\\n{}\"\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mformat\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mexc_value\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mshapes\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 170\u001b[0m \u001b[0mexc\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mexc\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mwith_traceback\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mtraceback\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 171\u001b[0;31m \u001b[0;32mraise\u001b[0m \u001b[0mexc\u001b[0m \u001b[0;32mfrom\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 172\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmsngr\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtrace\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0madd_node\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"_RETURN\"\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mname\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m\"_RETURN\"\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtype\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m\"return\"\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mvalue\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mret\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 173\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mret\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/pyro/poutine/trace_messenger.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m 163\u001b[0m args=args, kwargs=kwargs)\n\u001b[1;32m 164\u001b[0m \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 165\u001b[0;31m \u001b[0mret\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mfn\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 166\u001b[0m \u001b[0;32mexcept\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mValueError\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mRuntimeError\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 167\u001b[0m \u001b[0mexc_type\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mexc_value\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtraceback\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0msys\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mexc_info\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/pyro/infer/elbo.py\u001b[0m in \u001b[0;36mwrapped_fn\u001b[0;34m(*args, **kwargs)\u001b[0m\n\u001b[1;32m 135\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mfn\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 136\u001b[0m \u001b[0;32mwith\u001b[0m \u001b[0mpyro\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mplate\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"num_particles_vectorized\"\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mnum_particles\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mdim\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m-\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmax_plate_nesting\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 137\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mfn\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 138\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 139\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mwrapped_fn\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m<ipython-input-5-145b1ab1dfab>\u001b[0m in \u001b[0;36mguide\u001b[0;34m(self, x, y)\u001b[0m\n\u001b[1;32m 27\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 28\u001b[0m \u001b[0;31m# Get q(f) - variational (guide) distribution of latent function\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 29\u001b[0;31m \u001b[0mfunction_dist\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mpyro_guide\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 30\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 31\u001b[0m \u001b[0;31m# Use a plate here to mark conditional independencies\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/gpytorch/models/approximate_gp.py\u001b[0m in \u001b[0;36mpyro_guide\u001b[0;34m(self, input, beta, name_prefix)\u001b[0m\n\u001b[1;32m 60\u001b[0m \u001b[0;34m:\u001b[0m\u001b[0mparam\u001b[0m \u001b[0mstr\u001b[0m \u001b[0mname_prefix\u001b[0m\u001b[0;34m:\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mdefault\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m\"\"\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0mA\u001b[0m \u001b[0mname\u001b[0m \u001b[0mprefix\u001b[0m \u001b[0mto\u001b[0m \u001b[0mprepend\u001b[0m \u001b[0mto\u001b[0m \u001b[0mpyro\u001b[0m \u001b[0msample\u001b[0m \u001b[0msites\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 61\u001b[0m \"\"\"\n\u001b[0;32m---> 62\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0msuper\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mpyro_guide\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0minput\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mbeta\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mbeta\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mname_prefix\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mname_prefix\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 63\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 64\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mpyro_model\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0minput\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mbeta\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;36m1.0\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mname_prefix\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m\"\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/gpytorch/models/pyro/_pyro_mixin.py\u001b[0m in \u001b[0;36mpyro_guide\u001b[0;34m(self, input, beta, name_prefix)\u001b[0m\n\u001b[1;32m 15\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 16\u001b[0m \u001b[0;31m# Draw samples from q(f)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 17\u001b[0;31m \u001b[0mfunction_dist\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0minput\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mprior\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mFalse\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 18\u001b[0m function_dist = pyro.distributions.Normal(loc=function_dist.mean, scale=function_dist.stddev).to_event(\n\u001b[1;32m 19\u001b[0m \u001b[0mlen\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfunction_dist\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mevent_shape\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;34m-\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/gpytorch/models/approximate_gp.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, inputs, prior, **kwargs)\u001b[0m\n\u001b[1;32m 79\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0minputs\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdim\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 80\u001b[0m \u001b[0minputs\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0minputs\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0munsqueeze\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m-\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 81\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mvariational_strategy\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0minputs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mprior\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mprior\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/gpytorch/variational/variational_strategy.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, x, prior)\u001b[0m\n\u001b[1;32m 164\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mupdated_strategy\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mfill_\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;32mTrue\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 165\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 166\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0msuper\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m__call__\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mprior\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mprior\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/gpytorch/variational/_variational_strategy.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, x, prior)\u001b[0m\n\u001b[1;32m 121\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 122\u001b[0m \u001b[0;31m# Get p(u)/q(u)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 123\u001b[0;31m \u001b[0mvariational_dist_u\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mvariational_distribution\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 124\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 125\u001b[0m \u001b[0;31m# Get q(f)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/gpytorch/utils/memoize.py\u001b[0m in \u001b[0;36mg\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m 49\u001b[0m \u001b[0mkwargs_pkl\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mpickle\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdumps\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 50\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0;32mnot\u001b[0m \u001b[0m_is_in_cache\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcache_name\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mkwargs_pkl\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mkwargs_pkl\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 51\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0m_add_to_cache\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcache_name\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmethod\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m**\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mkwargs_pkl\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mkwargs_pkl\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 52\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0m_get_from_cache\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcache_name\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mkwargs_pkl\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mkwargs_pkl\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 53\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/gpytorch/variational/_variational_strategy.py\u001b[0m in \u001b[0;36mvariational_distribution\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 61\u001b[0m \u001b[0;34m@\u001b[0m\u001b[0mcached\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mname\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m\"variational_distribution_memo\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 62\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mvariational_distribution\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 63\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_variational_distribution\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 64\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 65\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mforward\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mx\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0minducing_points\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0minducing_values\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mvariational_inducing_covar\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/gpytorch/variational/_variational_distribution.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 40\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0m__call__\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 41\u001b[0m \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 42\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mforward\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 43\u001b[0m \u001b[0;31m# Remove after 1.0\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 44\u001b[0m \u001b[0;32mexcept\u001b[0m \u001b[0mNotImplementedError\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/gpytorch/variational/natural_variational_distribution.py\u001b[0m in \u001b[0;36mforward\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 58\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 59\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mforward\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 60\u001b[0;31m \u001b[0mmean\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mchol_covar\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0m_NaturalToMuVarSqrt\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mapply\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mnatural_vec\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mnatural_mat\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 61\u001b[0m \u001b[0mres\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mMultivariateNormal\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmean\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mCholLazyTensor\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mTriangularLazyTensor\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mchol_covar\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 62\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mres\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/gpytorch/variational/natural_variational_distribution.py\u001b[0m in \u001b[0;36mforward\u001b[0;34m(ctx, nat_mean, nat_covar)\u001b[0m\n\u001b[1;32m 115\u001b[0m \u001b[0;34m@\u001b[0m\u001b[0mstaticmethod\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 116\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mforward\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mctx\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mnat_mean\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mnat_covar\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 117\u001b[0;31m \u001b[0mmu\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mL\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0m_NaturalToMuVarSqrt\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_forward\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mnat_mean\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mnat_covar\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 118\u001b[0m \u001b[0mctx\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0msave_for_backward\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmu\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mL\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 119\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mmu\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mL\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/lib/python3.7/site-packages/gpytorch/variational/natural_variational_distribution.py\u001b[0m in \u001b[0;36m_forward\u001b[0;34m(nat_mean, nat_covar)\u001b[0m\n\u001b[1;32m 100\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mstr\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0me\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstartswith\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"cholesky\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 101\u001b[0m raise RuntimeError(\n\u001b[0;32m--> 102\u001b[0;31m \u001b[0;34m\"Non-negative-definite natural covariance. You probably \"\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 103\u001b[0m \u001b[0;34m\"updated it using an optimizer other than gpytorch.optim.NGD (such as Adam). \"\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 104\u001b[0m \u001b[0;34m\"This is not supported.\"\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;31mRuntimeError\u001b[0m: Non-negative-definite natural covariance. You probably updated it using an optimizer other than gpytorch.optim.NGD (such as Adam). This is not supported.\n Trace Shapes: \n Param Sites: \n Sample Sites: \nnum_particles_vectorized dist | \n value 32 | \n .u dist 32 1 | 20\n value 32 1 | 20"
]
}
],
"source": [
"pyro.clear_param_store()\n",
"\n",
"num_iter = 3000\n",
"num_particles = 32\n",
"\n",
"# Instantiate model\n",
"model = PVGPRegressionModel(name_prefix=\"test\", num_inducing=20)\n",
"\n",
"# Define optimizers separately for the variational and hyper-parameters\n",
"variational_ngd_optimizer = gpytorch.optim.NGD(model.variational_parameters(), num_data=train_y.size(0), lr=0.1)\n",
"hyperparameter_optimizer = torch.optim.Adam(model.hyperparameters(), lr=0.01)\n",
"\n",
"elbo = pyro.infer.Trace_ELBO(num_particles=num_particles, vectorize_particles=True, retain_graph=True)\n",
"\n",
"model.train()\n",
"\n",
"# Training loop \n",
"\n",
"iterator = tqdm.notebook.tqdm(range(num_iter))\n",
"\n",
"for i in iterator:\n",
" \n",
" model.zero_grad()\n",
" \n",
" variational_ngd_optimizer.zero_grad()\n",
" hyperparameter_optimizer.zero_grad()\n",
" \n",
" loss = elbo.differentiable_loss(model.model, model.guide, train_x, train_y)\n",
" loss.backward()\n",
" \n",
" variational_ngd_optimizer.step()\n",
" hyperparameter_optimizer.step()\n",
" \n",
" iterator.set_postfix(loss=loss, lengthscale=model.covar_module.base_kernel.lengthscale.item())"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.7.3"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment