add variational inference post notebook and rename the notebook of the chase the posterior blog post
This commit is contained in:
ritchie46 committed
1 parent
057084bc79
commit
81f46a8b7b
2 files changed
+387
No files matched your search
File renamed without changes.
@@ -0,0 +1,387 @@
|
||||
{
|
||||
"cells": [
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 1,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"import pandas as pd\n",
|
||||
"import matplotlib.pyplot as plt\n",
|
||||
"\n",
|
||||
"import numpy as np\n",
|
||||
"import torch\n",
|
||||
"from torch import nn\n",
|
||||
"from sklearn import datasets\n",
|
||||
"from torch.utils.data import DataLoader, TensorDataset\n",
|
||||
"from sklearn.preprocessing import StandardScaler"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 394,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"<matplotlib.collections.PathCollection at 0x7f6d652fbfd0>"
|
||||
]
|
||||
},
|
||||
"execution_count": 394,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
},
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA6AAAAFlCAYAAAD1f4WfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3df4xd6Vkf8OfN7JDchKIJzVLiSZY1VWogMYnbUZRqqgqWtE5JAOM2LZRWqK20qkQriqhbb1MpoVLkkaxSEKJSI0JBagoE2AwRW2oCE0pZNYFZZsMSNm4jUsLOpmQpTAvNkHidt3/MjD0e33t978z58Z5zPh9ptfaZ8b3n3nPuuef7/njelHMOAAAAqNsL2t4BAAAAhkEABQAAoBECKAAAAI0QQAEAAGiEAAoAAEAjBFAAAAAacV8bT/qyl70sP/jgg208NQAAADV74oknfj/nfP/R7a0E0AcffDA2NzfbeGoAAABqllL6nXHbDcEFAACgEQIoAAAAjRBAAQAAaIQACgAAQCMEUAAAABohgAIAANAIARQAAIBGCKAAAAA0QgAFAACgEQIoAAAAjbiv7R0AAADog/Wt7bh67Xo8u7Mbp5ZGcen8mbhwbrnt3SqKAAoAAHBC61vb8cijT8XujZsREbG9sxuPPPpURIQQeoghuAAAACd09dr1W+HzwO6Nm3H12vWW9qhMAigAAMAJPbuzO9f2oRJAAQAATujU0miu7UMlgAIAAJzQpfNnYrS4cMe20eJCXDp/pqU9KpMiRAAAACd0UGhIFdzpBFAAAIAKXDi3LHDegyG4AAAANEIABQAAoBECKAAAAI0QQAEAAGiEAAoAAEAjBFAAAAAaIYACAADQCAEUAACARgigAAAANEIABQAAoBECKAAAAI0QQAEAAGiEAAoAAEAjBFAAAAAaIYACAADQiBMH0JTSi1JKv5pS+khK6aMppe+pYscAAADol/sqeIzPRsRDOec/TiktRsSvpJR+Luf8oQoeGwAAgJ44cQDNOeeI+OP9vy7u/5dP+rgAAAD0SyVzQFNKCymlJyPi0xHxgZzzh8f8zsMppc2U0uZzzz1XxdMCAADQIZUE0JzzzZzz6yLiFRHx+pTSa8b8zrtyzis555X777+/iqcFAACgQyqtgptz3omIX4qIN1X5uAAAAHRfFVVw708pLe3/eRQRb4yIj530cQEAAOiXKqrgvjwifjSltBB7gfa9OeefreBxAQAA6JEqquD+RkScq2BfAAAA6LFK54ACAADAJAIoAAAAjRBAAQAAaIQACgAAQCMEUAAAABohgAIAANAIARQAAIBGCKAAAAA0QgAFAACgEQIoAAAAjRBAAQAAaIQACgAAQCMEUAAAABohgAIAANAIARQAAIBG3Nf2DgAAAP21vrUdV69dj2d3duPU0igunT8TF84tt71btEQABQAAarG+tR2PPPpU7N64GRER2zu78cijT0VECKEDJYACAAC1uHrt+q3weWD3xs24eu364AKonuA9AigAAFCLZ3d259reV3qCb1OECAAAqMWppdFc2/tqWk/w0AigAABALS6dPxOjxYU7to0WF+LS+TMt7VE79ATfJoACAAC1uHBuOa5cPBvLS6NIEbG8NIorF88ObtipnuDbzAEFAABqc+Hc8uAC51GXzp+5Yw5oxDB7giMEUAAAgFodBHBVcAVQAACA2ukJ3mMOKAAAAI0QQAEAAGiEAAoAAEAjBFAAAAAaoQgRAACVWt/aVu0TGEsABQCgMutb23esd7i9sxuPPPpURIQQChiCCwBAda5eu34rfB7YvXEzrl673tIeASURQAEAqMyzO7tzbQeGRQAFAKAyp5ZGc20HhkUABQBo0frWdqyubcTpy4/F6tpGrG9tt71LJ3Lp/JkYLS7csW20uBCXzp9paY+AkihCBADQkj4W7DnYb1VwgXEEUACAlkwr2NPlwHbh3HKn9x+ojyG4AAAtUbAHGBoBFACgJQr2AEMjgAIAtETBHmBozAEFAGiJgj3A0AigAAAtUrAHGBIBFAAAaMX61rYRAAMjgAIAAI3r4zq43JsACgAA1O5ob+dnPvd8L9fBZToBFAAAqNW43s5JrIPbb5ZhAQAAanX12vW7ejsnsQ5uvwmgAABAbda3tqf2eB5mHdz+MwQXAAAqprrrnoOht5MsjRbjJS+8b/Dv05AIoAAAUKF5q7v2OaxOG3o7WlyId3zjq3vzWpnNiYfgppRemVL6YErp6ZTSR1NK31nFjgEAQBeNC10H1V2POgir2zu7keN2WF3f2m5ob+s1raDQlYtnhc8BqmIO6PMR8d0556+MiDdExHeklL6qgscFAIDOmRS6xm2fJ6x20aSCQstLI+FzoE4cQHPOn8o5//r+n/8oIp6OCGcTAABTrW9tx+raRpy+/Fisrm30ptdvUugat32esNpFl86fidHiwh3bFBoatkqr4KaUHoyIcxHx4TE/eziltJlS2nzuueeqfFoAADqmz0NP5wld84TVLrpwbjmuXDwby0ujSLHX82no7bClnHM1D5TSF0bEf4mId+acH532uysrK3lzc7OS5wUAoHtW1zbGLs2xvDSKxy8/1MIeVWvWwkJHCxZF7IVVIY2uSyk9kXNeObq9kiq4KaXFiPjpiHjPvcInAAD0fejphXPLMwXIg9/paxVcOOrEATSllCLi3RHxdM75e0++SwAA9N2ppdHYHtC+DD2dx6xhFfqgijmgqxHxdyPioZTSk/v/fX0FjwsAQE8pTgPDdOIe0Jzzr0REqmBfAAAYCENPYZgqmQMKAADzMvQUhqfSZVgAAABgEj2gAABA62ZduoZuE0ABAIBWHV0PdXtnNx559KmICCG0ZwRQAACgVVevXb8VPg/s3rgZV69d72QA1Zs7mQAKAAC06tkxa8JO214yvbnTKUIEAAC06tTSaK7tJZvWm4sACgAAHNP61nasrm3E6cuPxeraRqxvbR/rcS6dPxOjxYU7to0WF+LS+TNV7Gaj+tSbWwdDcAEAqIy5b8NR5VDTg9/vw7lzamkU22PCZhd7c+sggAIAUAlz34al6sJBF84t9+I8uXT+zB2fg4ju9ubWwRBcAAAqYe7bsBhqOt6Fc8tx5eLZWF4aRYqI5aVRXLl4thfhugp6QAEAqIRAMiyGmk7Wl97cOugBBQCgEn2qZMq99alwEM0RQAEAqIRAMix1DDWtqqou5TIEFwCASvSpkimzqXKoqSJWwyCAAgBQGXPfOK6qq+pSJkNwAQCA1iliNQwCKAAA0DpFrIZBAAUAAFqniNUwmAMKAAC0ThGrYRBAAQCAIihi1X+G4AIAANAIARQAAIBGCKAAAAA0whxQAKAI61vbio8A9JwACgC0bn1rOx559KnYvXEzIiK2d3bjkUefiogQQgF6RAAFAFp39dr1W+HzwO6Nm3H12nUBtIf0dnPAuTA8AigA0Lpnd3bn2k536e3mgHNhmBQhAgBad2ppNNd2umtabzfD4lwYJgEUAGjdpfNnYrS4cMe20eJCXDp/pqU9oi56uzngXBgmARQAaN2Fc8tx5eLZWF4aRYqI5aVRXLl41jC8HtLbzQHnwjCZAwoAFOHCueVeBE5FVaa7dP7MHfP+IvR2D5VzYZgEUACAiiiqcm8H74OQjnNhmFLOufEnXVlZyZubm40/LwBAnVbXNmJ7zPy15aVRPH75oRb2CKAdKaUncs4rR7ebAwoAUBFFVQCmE0ABACqiqArAdAIoAEBFLCcDMJ0iRAAAFVFUBWA6ARQAoEJ9WU4GoA6G4AIAANAIARQAAIBGCKAAAAA0whxQAABqtb61rTATveX8no8ACgDQoKHdrK5vbccjjz4VuzduRkTE9s5uPPLoUxERvX7dDIPze36G4AIANOTgZnV7Zzdy3L5ZXd/abnvXanP12vVbN+cHdm/cjKvXrre0R1Ad5/f8BFAAgIYM8Wb12Z3dubZDlzi/5yeAAgA0ZIg3q6eWRnNthy5xfs9PAAUAGrG+tR2raxtx+vJjsbq20ethp5MM8Wb10vkzMVpcuGPbaHEhLp0/09IeQXWc3/MTQAGA2g1x7uM4Q7xZvXBuOa5cPBvLS6NIEbG8NIorF88q0EIvOL/nl3LOjT/pyspK3tzcbPx5AYB2rK5txPaYYabLS6N4/PJDLexRe4ZWBRcYppTSEznnlaPbK1mGJaX0wxHxloj4dM75NVU8JgDQH0Oc+zjJhXPLAicwWFUNwf2RiHhTRY8FAPTMEOc+DoW5vcA8KgmgOedfjog/qOKxAID+GeLcxyEwtxeYV2NFiFJKD6eUNlNKm88991xTTwsAFEChjv5Z39qO737vRwa3rilwMpXMAZ1FzvldEfGuiL0iRE09LwBQBnMf++Og5/PmhGKWQ5zbC8zGMiwAAMzl6rXrd/V8HmZuLzCJAAoAwFym9XCa2wtMU0kATSn9WET8t4g4k1J6JqX0D6p4XAAAyjOph3MhJXN7gamqqoL7rTnnl+ecF3POr8g5v7uKxwUAoDyTqhr/67/5WuETmKqxIkQAAAfWt7bj6rXr8ezObpxaGsWl82cElw45OFaOITAvARQAaNRBBdWDIjYHa0dGhABTqEkNBo4XMC8BFABo1LgKqgdrRwo05dFg0G9GI9A0VXABgEZNqqBq7cgyTWswoNsOGhe2d3Yjx+3GhfWt7bZ3jR4TQAGARk2qoGrtyDJpMOgvjQu0QQAFABo1qYKqtSPLpMGgvzQu0AYBFABo1IVzy3Hl4tlYXhpFiojlpZG1IwumwaC/NC7QBkWIAIDGqaDaHZZc6a9L58/cUWAqQuMC9RNAAQCYSoPBbX2qGqtxgTYIoAAAMIM+LkmjcYGmCaAAADPqU+8X87OGLZycAAoAMIM+9n4xH1Vj4eRUwQUAmIE1E1E1Fk5OAAUAmIHeLyxJQ9vWt7ZjdW0jTl9+LFbXNmJ9a7vtXZqbIbgAADM4tTSK7TFLine truncated
|
||||
"text/plain": [
|
||||
"<Figure size 1152x432 with 1 Axes>"
|
||||
]
|
||||
},
|
||||
"metadata": {
|
||||
"needs_background": "light"
|
||||
},
|
||||
"output_type": "display_data"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"w0 = 0.125\n",
|
||||
"b0 = 5.\n",
|
||||
"x_range = [-20, 60]\n",
|
||||
"\n",
|
||||
"def load_dataset(n=150, n_tst=150):\n",
|
||||
" np.random.seed(43)\n",
|
||||
" \n",
|
||||
" def s(x):\n",
|
||||
" g = (x - x_range[0]) / (x_range[1] - x_range[0])\n",
|
||||
" return 3 * (0.25 + g**2.)\n",
|
||||
" \n",
|
||||
" x = (x_range[1] - x_range[0]) * np.random.rand(n) + x_range[0]\n",
|
||||
" eps = np.random.randn(n) * s(x)\n",
|
||||
" y = (w0 * x * (1. + np.sin(x)) + b0) + eps\n",
|
||||
" y = (y - y.mean()) / y.std()\n",
|
||||
" idx = np.argsort(x)\n",
|
||||
" x = x[idx]\n",
|
||||
" y = y[idx]\n",
|
||||
" return y[:, None], x[:, None]\n",
|
||||
"\n",
|
||||
"Y, X = load_dataset()\n",
|
||||
"\n",
|
||||
"plt.figure(figsize=(16, 6))\n",
|
||||
"plt.scatter(X, Y)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 395,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"X = torch.tensor(X, dtype=torch.float)\n",
|
||||
"Y = torch.tensor(Y, dtype=torch.float)\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 396,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"class MaximumLikelihood(nn.Module):\n",
|
||||
" def __init__(self):\n",
|
||||
" super().__init__()\n",
|
||||
" self.out = nn.Sequential(\n",
|
||||
" nn.Linear(1, 20),\n",
|
||||
" nn.ReLU(),\n",
|
||||
" nn.Linear(20, 1)\n",
|
||||
" )\n",
|
||||
" def forward(self, x):\n",
|
||||
" return self.out(x)\n",
|
||||
" \n",
|
||||
"epochs = 200\n",
|
||||
"m = MaximumLikelihood()\n",
|
||||
"optim = torch.optim.Adam(m.parameters(), lr=0.01)\n",
|
||||
"\n",
|
||||
"for epoch in range(epochs):\n",
|
||||
" optim.zero_grad()\n",
|
||||
" y_pred = m(X)\n",
|
||||
" loss = (0.5 * (y_pred - Y)**2).mean()\n",
|
||||
" loss.backward()\n",
|
||||
" optim.step()\n",
|
||||
" "
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 397,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA6AAAAFlCAYAAAD1f4WfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3de3Cc2Xnf+d/pC9ANNNAgLg2C4HBIApwZjqDRaIS5SEOPZWkkUZZiWY7lUVLrcu2uw/VWokq01mbt9R9OtiqVrdokE682uxWW5dhJHLvsOIpvMe2xFIeekTUSRxpJ1EAmAWI4HAyIK9FAA91AX87+0Y0GGsQd3e+l+/upYgFovOw+3edF9/uc85znGGutAAAAAACotYDbDQAAAAAANAYCUAAAAACAIwhAAQAAAACOIAAFAAAAADiCABQAAAAA4AgCUAAAAACAI0JuPGh3d7c9ffq0Gw8NAAAAAKix1157bdZa27P1dlcC0NOnT+vatWtuPDQAAAAAoMaMMbe3u50UXAAAAACAIwhAAQAAAACOIAAFAAAAADiCABQAAAAA4AgCUAAAAACAIwhAAQAAAACOIAAFAAAAADiCABQAAAAA4AgCUAAAAACAIwhAAQAAAACOIAAFAAAAADgi5HYDAAAAAMCvRiaTunJ9ShMLafV3RHVxqFfn++JuN8uzmAEFAAAAgEMYmUzq8tVxJdNZ9cUjSqazunx1XCOTSbeb5lkEoAAAAABwCFeuTykeDSseDStgTPn7K9en3G6aZxGAAgAAAMAhTCyk1RapXNXYFglpYiHtUou8jwAUAAAAAA6hvyOqpUyu4ralTE79HVGXWuR9BKAAAAAAcAgXh3qVTGeVTGdVsLb8/cWhXreb5lkEoAAAAABwCOf74rr03BnFo2FNJjOKR8O69NwZquDugm1YAAAAAOCQzvfFCTgPgBlQAAAAAIAjCEABAAAAAI4gAAUAAAAAOIIAFAAAAADgCAJQAAAAAIAjCEABAAAAAI4gAAUAAAAAOIIAFAAAAADgCAJQAAAAAIAjCEABAAAAAI4gAAUAAAAAOIIAFAAAAADgCAJQAAAAAIAjCEABAAAAAI4gAAUAAAAAOOLIAagxJmKM+YYx5jvGmO8bY/5xNRoGAAAAAKgvoSrcx6qkD1lrU8aYsKSXjTF/Yq39ehXuGwAAAABQJ44cgFprraRU6cdw6Z896v0CAAAAAOpLVdaAGmOCxpjXJU1Lesla++o2x1wyxlwzxlybmZmpxsMCAAAAAHykKgGotTZvrX1c0klJTxljhrY55rK1dthaO9zT01ONhwUAAAAA+EhVq+Baaxck/YWki9W8XwAAAACA/1WjCm6PMaaj9H1U0vOSfnDU+wUAAAAA1JdqVMHtk/QbxpigigHt71hr/6gK9wsAAAAAqCPVqIL7XUnvrUJbAAAAAAB1rKprQAEAAAAA2AkBKAAAAADAEQSgAAAAAABHEIACAAAAABxBAAoAAAAAcAQBKAAAAADAEQSgAAAAAABHEIACAAAAABxBAAoAAAAAcAQBKAAAAADAEQSgAAAAAABHEIACAAAAABxBAAoAAAAAcAQBKAAAAADAESG3GwAAAACg/oxMJnXl+pQmFtLq74jq4lCvzvfF3W4WXMYMKAAAAICqGplM6vLVcSXTWfXFI0qms7p8dVwjk0m3mwaXEYACAAAAqKor16cUj4YVj4YVMKb8/ZXrU243DS4jBRcAAABAVU0spNUXj1Tc1hYJaWIh7VKL3EMqciVmQAEAAABUVX9HVEuZXMVtS5mc+juiLrXIHaQi348AFAAAAEBVXRzqVTKdVTKdVcHa8vcXh3rdbpqjSEW+HwEoAAAAgKo63xfXpefOKB4NazKZUTwa1qXnzjRc6unEQlptkcpVj42airyONaAAAAAAqu58X7zhAs6t+juiSqazikfD5dsaMRV5M2ZAAQAAAKAGSEW+HwEoAAAAANQAqcj3IwUXAAAAAGqEVORKzIACAAAAABxBAAoAAAAAcAQBKAAAAADAEawBBQAAQNWMTCZ15fqUJhbS6u+I6uJQL+vfAJQxAwoAAICqGJlM6vLVcSXTWfXFI0qms7p8dVwjk0m3mwbAIwhAAQAAUBVXrk8pHg0rHg0rYEz5+yvXp9xuGgCPIAAFAABAVUwspNUWqVzh1RYJaWIh7VKLAHgNASgAAACqor8jqqVMruK2pUxO/R1Rl1oEwGsIQAEAAFAVF4d6lUxnlUxnVbC2/P3FoV63mwbAI6iCCwAA4JJ6qxh7vi+uS8+dqXhOLzx50tfPCUB1EYACAAC4YL1ibDwarqgYe+m5M74O2M73xX3dfgC1RQouAACAC6gYC6AREYACAAC4gIqxABoRASgAAIALqBgLoBERgAIAALiAirEAGhEBKAAAgAvWK8bGo2FNJjOKR8O+L0AEAHuhCi4AAIBLqBgLoNEQgAIAAABwVL3tgYv9IwUXAAAAgGPW98BNprMVe+COTCbdbhocwAwoAAAAgJrZOts5u5Qp73srqfz1yvUpZkEbADOgAAAAAGpiu9nOvxyd02qucgsi9sBtHASgAAAAAGriyvWp8mxnwBjFo2EdawnrjXeWKo5jD9zGQQAKAAAAoOpGJpN66Y0pff3WrL5+a06zqYwk6Xxfm+6tsAduo2INKAAAAFAlB63uWq/VYNdTb8NBI9mAMtm8Xru9oPc92KFIOKQLg12KR8Pl5/3Ckyfr4nljb8Zae7Q7MOYBSf9W0nFJBUmXrbW/stv/GR4etteuXTvS4wIAAABesh50xaNhtUVCWsrklExndem5M9sGVwc93k9efOmGkums1nJ5feutBTWHArLWyhijsz2xuniO2J0x5jVr7fDW26uRgpuT9PPW2vOSnpH0d40xj1bhfgEAAADf2G69Yzwa1pXrU1U53k8mFtJqi4TU0xbRE6c61BwOKpu3WssXCD4b3JFTcK21k5ImS98vGWNGJPVLeuOo9w0AAID6VI+ppxMLafXFIxW37Vbd9aDH+0l/R1TJdFbxaFg9bRH1tEXKP/u9n3E0VS1CZIw5Lem9kl6t5v0CAACgfmy3Ncflq+MamUy63bQj6e+IailTub3IbtVdD3q8n1wc6i0XF6LQEDarWgBqjIlJ+j1J/8Bau7jN7y8ZY64ZY67NzMxU62EBAADgM/WaenrQoKueg7TzfXFdeu6M4tGwJpMZxaNhUm8hqQpFiCTJGBOW9EeS/tRa+y/2Op4iRAAAAI3rC7/7HfXFIwoYU76tYK0mkxn9s8+8x8WWHR1VcIGinYoQHXkNqDHGSPqSpJH9BJ8AAABobJvXB66rl9TT833xAwWQBz0e8LtqpOA+K+mnJX3IGPN66d+PVuF+AQAAUIfqOfUUwO6qUQX3ZUlmzwMBAAAAbawP3Jx6+sKTJ5kJBBrAkQNQAAAA4KBIPQUaU1W3YQEAAAAAYCcEoAAAAAAAR5CCCwAAAMAVbEPTeAhAAQAAADhuZDKpy1fHFY+G1RePKJnO6vLVcV167owvg1CC6f0hBRcAAACA465cn1I8GlY8GlbAmPL3V65Pud20A1sPppPpbEUwPTKZdLtpnkMACgAAAMBxEwtptUUqEzLbIiFNLKRdatHh1VMwXWuk4AIAAAA4kGqkm/Z3RJVMZxWPhsu3LWVy6u+IVru5NTexkFZfPFJxm1+D6VpjBhQAAABHNjKZ1Isv3dAXfvc7evGlG6Qe1rFqpZteHOpVMp1VMp1Vwdry9xeHemvU8trp74hqKZOruM2vwXStEYACAADgSFj/1liqlW56vi+uS8+dUTwa1mQyo3g07NsCRPUUTNcaKbgAAAA4ks0BiaTy1yvXp3wZTGB31Uw3Pd8Xr4tzZD2Y3pyW/MKTJ+viuVUbASgAAACOhPVvjaWe1m5WU70E07VGCi4AAACOhPVvjYV0UxwFM6AAAAA4kotDvbp8dVxSceZzKZNTMp3VC0+edLllqIVqp5tWo6Iu/MNYax1/0OHhYXvt2jXHHxcAAAC1QRCBw1gvYBWPhisGL/xajAgbjDGvWWuHt97ODCgAAACOjPVvOAwKWDUe1oACAAAAcMXEQlptkco5MQpY1TcCUAAAAACuoIBV4yEABQAAAOAKKuo2HgJQAAAAAK5Yr6gbj4Y1mcwoHg1TgKjOUYQIAAAAgGsoYNVYmAEFAAAAADiCABQAAAAA4AgCUAAAAACAI1gDCgAAXDcymdSV61OaWEirvyOqi0O9rAkDgDrEDCgAAHDVyGRSl6+OK5nOqi8eUTKd1eWr4xqZTLrdNABAlTEDCgAAXHXl+pTi0bDi0bAklb9euT7FLGidYaYbEudBo2MGFAAAuGpiIa22SOWYeFskpImFtEstQi0w0w2J8wAEoAAAwGX9HVEtZXIVty1lcurviLrUItTC5pnugDHl769cn3K7aXAQ5wEIQAEAgKsuDvUqmc4qmc6qYG35+4tDvW43DVXETDckzgMQgAIAAJed74vr0nNnFI+GNZnMKB4N69JzZ1gTVmeY6YbEeQCKEAEAAA843xevi4CT4io7uzjUq8tXxyUVZ7yWMjkl01m98ORJl1sGJ3EegBlQAACAKqC4yu6Y6YbEeQBmQAEAAKqC7WT2Vi8z3TgazoPGxgwoAABAFVBcBQD2RgAKAABQBRRXAYC9EYACAABUAdvJAMDeCEABAACqgOIqALA3ihABAABUCcVVAGB3zIACAAAAABxBAAoAAAAAcAQBKAAAAADAEawBBQAAQE2MTCZ15fqUJhbS6u+I6uJQL2tkURecPrfnl9c0mUzrXSf8//dDAAoAAOCQRgrIRiaTunx1XPFoWH3xiJLprC5fHacyMHyvVue2tVaTyYxGp1PFfzOp8vfzy2tqi4T03V/+qIwxVXw2ziMABQAAcECjBWRXrk8pHg0rHg1LUvnrletTdfl80TiOem7n8gXduZfWzamlcpA5Np3S2MyyUquLine truncated
|
||||
"text/plain": [
|
||||
"<Figure size 1152x432 with 1 Axes>"
|
||||
]
|
||||
},
|
||||
"metadata": {
|
||||
"needs_background": "light"
|
||||
},
|
||||
"output_type": "display_data"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"with torch.no_grad():\n",
|
||||
" plt.figure(figsize=(16, 6))\n",
|
||||
" plt.plot(X, m(X))\n",
|
||||
" plt.scatter(X, Y, alpha=0.5)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 398,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"tensor(0.4997)\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"with torch.no_grad():\n",
|
||||
" print(torch.abs(m(X) - Y).mean())"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"$$ \\text{ELBO} = E_{\\theta \\sim Q}[\\log P(D|\\theta)] + D_{KL}(Q(\\theta)||P(\\theta)) $$\n",
|
||||
"\n",
|
||||
"$$ \\text{ELBO} = E_{\\theta \\sim Q}[\\log P(D|\\theta)] + \\int Q(\\theta) \\log \\frac{P(\\theta)}{Q(\\theta)}d\\theta $$"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 399,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"def ll_gaussian(y, mu, log_var):\n",
|
||||
" sigma = torch.exp(0.5 * log_var)\n",
|
||||
" return -0.5 * torch.log(2 * np.pi * sigma**2) - (1 / (2 * sigma**2))* (y-mu)**2"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 451,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"class VI(nn.Module):\n",
|
||||
" def __init__(self):\n",
|
||||
" super().__init__()\n",
|
||||
" \n",
|
||||
" self.q_mu = nn.Sequential(\n",
|
||||
" nn.Linear(1, 20),\n",
|
||||
" nn.ReLU(),\n",
|
||||
" nn.Linear(20, 10),\n",
|
||||
" nn.ReLU(),\n",
|
||||
" nn.Linear(10, 1)\n",
|
||||
" )\n",
|
||||
" self.q_log_var = nn.Sequential(\n",
|
||||
" nn.Linear(1, 20),\n",
|
||||
" nn.ReLU(),\n",
|
||||
" nn.Linear(20, 10),\n",
|
||||
" nn.ReLU(),\n",
|
||||
" nn.Linear(10, 1)\n",
|
||||
" )\n",
|
||||
" \n",
|
||||
" def reparameterize(self, mu, log_var):\n",
|
||||
" # std can not be negative, thats why we use log variance\n",
|
||||
" sigma = torch.exp(0.5 * log_var) + 1e-5\n",
|
||||
" eps = torch.randn_like(sigma)\n",
|
||||
" return mu + sigma * eps\n",
|
||||
" \n",
|
||||
" def forward(self, x):\n",
|
||||
" mu = self.q_mu(x)\n",
|
||||
" log_var = self.q_log_var(x)\n",
|
||||
" return self.reparameterize(mu, log_var), mu, log_var"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 452,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"def elbo(y_pred, y, mu, log_var):\n",
|
||||
" # likelihood of observing y given Variational mu and sigma\n",
|
||||
" likelihood = ll_gaussian(y, mu, log_var)\n",
|
||||
" \n",
|
||||
" # prior probability of y_pred\n",
|
||||
" log_prior = ll_gaussian(y_pred, 0, torch.log(torch.tensor(1.)))\n",
|
||||
" \n",
|
||||
" # variational probability of y_pred\n",
|
||||
" log_p_q = ll_gaussian(y_pred, mu, log_var)\n",
|
||||
" \n",
|
||||
" # by taking the mean we approximate the expectation\n",
|
||||
" return (likelihood + log_prior - log_p_q).mean()\n",
|
||||
"\n",
|
||||
"def det_loss(y_pred, y, mu, log_var):\n",
|
||||
" return -elbo(y_pred, y, mu, log_var)\n",
|
||||
"\n",
|
||||
" \n",
|
||||
" "
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 453,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"epochs = 1500\n",
|
||||
"\n",
|
||||
"m = VI()\n",
|
||||
"optim = torch.optim.Adam(m.parameters(), lr=0.005)\n",
|
||||
"\n",
|
||||
"for epoch in range(epochs):\n",
|
||||
" optim.zero_grad()\n",
|
||||
" y_pred, mu, log_var = m(X)\n",
|
||||
" loss = det_loss(y_pred, Y, mu, log_var)\n",
|
||||
" loss.backward()\n",
|
||||
" optim.step()"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 454,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"tensor(0.5258)\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"with torch.no_grad():\n",
|
||||
" print(torch.abs(m(X)[1] - Y).mean())"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 455,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"# draw samples from Q(theta)\n",
|
||||
"with torch.no_grad():\n",
|
||||
" y_pred = torch.cat([m(X)[0] for _ in range(1000)], dim=1)\n",
|
||||
" \n",
|
||||
"# Get some quantiles\n",
|
||||
"q1, mu, q2 = np.quantile(y_pred, [0.05, 0.5, 0.95], axis=1)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 464,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"<matplotlib.collections.PolyCollection at 0x7f6d648794e0>"
|
||||
]
|
||||
},
|
||||
"execution_count": 464,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
},
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA6AAAAFlCAYAAAD1f4WfAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdd3Tkd33v/+d3ujTqWm3vXntxA5fF2AZibIoNOLkOBEIcfAkQnIIJJDkOcbgnN4SbwDn+JaSY5OIQIHECBC7GITjggL3gGBbjdYF1k7dq+6prevmW3x8fdY369Hk9ztGRNBp95zO7I+n7+r4/n/fH8jwPERERERERkVLzVXoAIiIiIiIi0hgUQEVERERERKQsFEBFRERERESkLBRARUREREREpCwUQEVERERERKQsFEBFRERERESkLAKVeNA1a9Z427dvr8RDi4iIiIiISIk9+eSTg57n9cy+vSIBdPv27ezfv78SDy0iIiIiIiIlZllWX6HbNQVXREREREREykIBVERERERERMpCAVRERERERETKQgFUREREREREykIBVERERERERMpCAVRERERERETKQgFUREREREREykIBVERERERERMpCAVRERERERETKQgFUREREREREyiJQ6QGIiIiIiIjUgweePsXdD/VyejTNxo4m7rxxN7dcvqnSw6oqCqAiIiIiIiKr9MDTp7jr/gOk8w4Ap0bT3HX/AQCF0Gk0BVdERERERGSV7n6odzJ8TkjnHe5+qLdCI6pOCqAiIiIiIiKrdHo0vazbG5UCqIiIiIiIyCpt7Gha1u2NSgFURERERERkle68cTdNQf+M25qCfu68cXeFRlSd1IRIRERERERklSYaDakL7sIUQEVERERERIrglss3KXAuQlNwRUREREREpCwUQEVERERERKQsFEBFRERERESkLBRARUREREREpCwUQEVERERERKQsFEBFRERERESkLBRARUREREREpCwUQEVERERERKQsFEBFRERERESkLBRARUREREREpCwUQEVERERERKQsFEBFRERERESkLBRARUREREREpCwUQEVERERERKQsVh1ALcuKWJb1E8uyfmpZ1nOWZX28GAMTERERERGR+hIowjGywA2e5yUsywoCj1mW9W3P835chGOLiIiIiIhInVh1APU8zwMS458Gx9+81R5XRERERERE6ktR1oBaluW3LOsZoB/4rud5jxe4z+2WZe23LGv/wMBAMR5WREREREREakhRAqjneY7neZcBm4GrLMu6pMB97vU8b4/neXt6enqK8bAiIiIiIiJSQ4raBdfzvFHg+8BNxTyuiIiIiIiI1L5idMHtsSyrY/zjJuANwIurPa6IiIiIiIjUl2J0wd0A/JNlWX5MoP2q53nfKsJxRUREREREpI4Uowvuz4DLizAWERERERERqWNFXQMqIiIiIiIiMh8FUBERERERESkLBVAREREREREpCwVQERERERERKQsFUBERERERESkLBVAREREREREpCwVQERERERERKQsFUBERERERESkLBVAREREREREpCwVQERERERERKQsFUBERERERESkLBVAREREREREpCwVQERERERERKQsFUBERERERESmLQKUHICIiIiIi9euBp09x90O9nB5Ns7GjiTtv3M0tl2+q9LCkQhRARURERESkJB54+hR33X+AdN4B4NRomrvuPwCgENqgFEBFRERERKQk7n6odzJ8TkjnHe5+qLfhAqgqwYYCqIiIiIiIlMTp0fSybq9XqgRPURMiEREREREpiY0dTcu6vV4tVAluNAqgIiIiIiJSEnfeuJumoH/GbU1BP3feuLtCI6oMVYKnKICKiIiIiEhJ3HL5Jj75tkvZ1NGEBWzqaOKTb7u04aadqhI8RWtARURERESkZG65fFPDBc7Z7rxx94w1oNCYlWBQABURERERESmpiQCuLrgKoCIiIiIiIiWnSrChNaAiIiIiIiLlkB6Bc89BNlHpkVSMAqiIiIiISCkkhyB2ptKjkGoRPwsjx8DJwdAhSPRXekQVoSm4IiIiIiLFlkvC2AnAM4GjYytYVqVHJZXgeTB6HNLD02+E2CnzOunYBr5l1AXtLGBBIFTskZaFKqAiIiIiIsXk2KbShWc+Tw/D0GFwnYW+a/XSI5AaXvx+Uj6Obaqd6Xn+XzKjMNgL+czSjxk7zeRrqwapAioiIiIiRfXA06cat9un501Ns5wuF4fBg9C1sziVq4kw6/Ob9/mMqbJ5LiQHoHUDRNpW/ziyck7e/J872YXvZ2dg8CVTJW/qWPi+uaQJrW0bizfOMlMFVERERESK5oGnT3HX/Qc4NZrGA06Nprnr/gM88PSpSg+tPGKnTdgsxE6boJFLrfz4ngeJAeh/fmoNoefBaJ8JnwD5FAwfNlXXfHrljyWGYy/8dc+DzBjYsy46xM8sHj4nj+HAyFHz+vEWqG7WwZpiVUBFREREpGjufqiXdH7mVNN03uHuh3rrvwqaHoHkIo1l3DwMHYTO7RBpX/7xY9NCTWoQWtZB/LQJnbNlYzAQg+ZuUxH1B5f3eHPG7prxO/nx97ap2K32uPPxPHDtqcdznWkf21PB0Oc3b5YffIHxj33TPp5++wrW4eYShSuT+bSZ8pweMWPyBcx6zkibqUivZDp04py5QNG5HfyzolpmbP6LGzVEAVREREREiub0aOGK23y3VyXPg2zcTHXEgo4ti39PPm2mwC7p+C4MH4X2LRDtXvz+ueR4ZXXW1h2ubSqfmdGFvz81ZEJSdK0JrAs1vHFss15xesicCJ1egTWsbn7h6aCeZx47M2qmHy9FLgnDR8zzK7aJYGr5TXAONkGoBcKt84fT6QHUsc3zSQ/PDf2ubSrPrRvGK88rXKeZi5t1oZ3bIRQ1t3leXVQ/QQFURERERIpoY0cTpwqEzY0dTRUYzTJ4ngka6RFTaZoefkJRaO6a/3tdxwTKiSmwS3tAGDtuqpnzBTg7a4LnQgFzsfA5+XAuJM6aYN21c251zfMgOWimjRYKmvNJDZnANTu8ua75WrLfrIe1/Es/pp0tTfgE8+8wsT7XTpsqMefM+CLt5i3cNjOkZxOQHjWhMxNj0WAZL0JQdHJm/ejERYrUsBlvHVAAFREREamgemvYc+eNu7nr/gMzpuE2Bf3ceePuCo5qAdnpoTNf+D5jJ00IDYQLf33k2NLX+s2WODe+Tcu2qRDn2CYsJgcperfTfNJMAe46b2YzpKFDcyusS+GOVwQnArqTN02QkoMzg+xywrkzz/9DKXmOCZjpYVMlDbeZqmeoxTQJGjla/jFNXKTIJ8eDb31QABURERGpkImGPRNhbaJhD1CzIXRi3FUdqnPJ8WmUo/OHzuk8B0b6YM35cyt98bPjVbRVSI+Y0NW53VS6EueWV4VcLjszFUKDERNuVhI+J6SGINhsqp2pYQqHZs9UWZeyBnN2B+Fy81xTWZ6Ygl3pLU9SQ5V9/CJTABURERGpkHpt2HPL5Zuqb/y51NRaxJUEnHzShM22DVO3ZcaKM90STAA892xxjrUUTm48hO40FcvVyCVg4IXF7+c6c6f+zje2qlG7+21WKwVQERERkQqpi4Y91WJ6h1QnP7N7ai618imy0yXOmWY14RazTnGkb/XHrCTXNlNvl7V2dRXme5zU8Mw1tpWYgitlowAqIiIiUiE127CnnKZ3YZ0eKmd/XpYQNb7f5poLxpsOlXCabLmUK3yCaaIzfd0pmCryaJ9ZXzvR8bWqKqBSbAqgIiIiIhVScw17yiE5aNa8TQTMapsC6eSg/4X6CJ/llhqeufep68DoCfPx2Eno2W1u079tXVMAFREREamQmmjYU075tAki1RY6Z1NAWpnMmAmYvvEtWWKnpppA5VOQGDDTm6WuKYCKiIiIVFBVNuypBM+D0eNUffiUVfBMI6joGtN5NzXE3t4B7tvXx0A8Q3drM++84VW8ZWulxyml5Fv8LiIiIiIiJRY/a6pgUt9Sw+C6MHaCvb0D3PPIIfrjGTxgMJ7iHx58jL29q+zKK1VNFVARERERqaxcynSYlbpmqp37ORL/b3a2umTyNll75nTmvJ3nvn19XL+7p0KjlFJTABURERGRytHU24YwUe3M2g4tQH98/vsOxDNlG5eUn6bgioiIiMjqeKsIj6khsz2H1LX79vXNqXbOp6c1UuLRSCUpgIqIiIjI6iQHV/Z9rgPxM8Udi1Sdvb0D9C+xqhkO+Lntmm0lHpFUkqbgioiIiMjqpAYh0gaB8PK+L35mfK/P+jO9u2tPa4Tbrtk277rG5dy31kxMvZ1PazhAUyhQl89dClt1ALUsawvwz8B6wAXu9Tzvr1d7XBERERGpAa4DdgbSo9C6bv77ZeOQS45P1/XM+5VWTqvc9PWOAP3xzGQImx2ulnPfWrTQ1NtwwM/t151XF89Tlq4YU3Bt4Pc9z7sQuBr4oGVZFxXhuCIiIiJS7XJJ8z4zOv99MjEYPmIqnomzpuNtsp96bTxUKHRlbYf79vWt6r61aKGGQnfcsEvhswGtugLqed4Z4Mz4x3HLsl4ANgHPr/bYIiIiIlLF7CzkxxsI5VPm89nTcLNxGDkKnjvn2+t16ul8oavQ7cu5by3qaY0UXP+5tjVSF//XsnxFbUJkWdZ24HLg8WIeV0RERETKLJsAZ4H1mfGzMPgS5BJTt6VnVUGzCVP5nCd83vPIIfrjGTympp7u7R0ozvgraL4uroVuX859a9Ft12wjHPDPuE2Nhhpb0QKoZVktwNeBj3ieFyvw9dsty9pvWdb+gYHa/8UiIiIiUrcyYzB8GIYOgp2b+bWJfTsnGghlp532zZ6GGztdMHxCfU89XU7oqveAdv3uHu64YRdrWyNYmMqnpt42tqJ0wbUsK4gJn//qed79he7jed69wL0Ae/bsqc8J/yIiIiLVyvPAsha/X3oURo4BnmkuNHQLine truncated
|
||||
"text/plain": [
|
||||
"<Figure size 1152x432 with 1 Axes>"
|
||||
]
|
||||
},
|
||||
"metadata": {
|
||||
"needs_background": "light"
|
||||
},
|
||||
"output_type": "display_data"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"plt.figure(figsize=(16, 6))\n",
|
||||
"plt.scatter(X, Y)\n",
|
||||
"plt.plot(X, mu)\n",
|
||||
"plt.fill_between(X.flatten(), q1, q2, alpha=0.2)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 300,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"def det_loss(y, y_pred, mu, log_var): \n",
|
||||
" reconstruction_error = (0.5 * (y - y_pred)**2).sum()\n",
|
||||
" kl_divergence = (-0.5 * torch.sum(1 + log_var - mu**2 - log_var.exp()))\n",
|
||||
"\n",
|
||||
" return (reconstruction_error + kl_divergence).sum()"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 422,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"tensor(0.5044, grad_fn=<MeanBackward0>)"
|
||||
]
|
||||
},
|
||||
"execution_count": 422,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"source": []
|
||||
},
|
||||
{
|
||||
"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": 4
|
||||
}
|
||||
Reference in new issue
Block a user