core team vi training init
This commit is contained in:
ritchie46 committed
1 parent
a755cba5fb
commit
f9e66aa340
20 files changed
+489
No files matched your search
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
File renamed without changes.
@@ -0,0 +1,489 @@
|
|||||||
|
{
|
||||||
|
"cells": [
|
||||||
|
{
|
||||||
|
"cell_type": "markdown",
|
||||||
|
"metadata": {},
|
||||||
|
"source": [
|
||||||
|
"# Expectation maximization\n",
|
||||||
|
"\n",
|
||||||
|
"$$ \\begin{eqnarray} \n",
|
||||||
|
" z_i &\\sim & Multinomial(\\phi) \\\\\n",
|
||||||
|
" x_i|z_i & \\sim & \\mathcal{N}(\\mu_{z_i}, \\sigma_{z_i})) \\\\\n",
|
||||||
|
"\\end{eqnarray}$$ \n"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "code",
|
||||||
|
"execution_count": 379,
|
||||||
|
"metadata": {},
|
||||||
|
"outputs": [
|
||||||
|
{
|
||||||
|
"data": {
|
||||||
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA6IAAAFlCAYAAADxilWiAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdd3RUdf7/8eedNNILKbQgLSF0pEuXSBOpFlC/KK7r6lrBtf/c9csXXdsiWLCCiKIioiIgJXRFQAjSQ5eSACkE0kmd+/sjyIoECJCZO0lej3M4x9yZzPuFejJ5zf3czzVM00RERERERETEWWxWBxAREREREZHqRUVUREREREREnEpFVERERERERJxKRVREREREREScSkVUREREREREnEpFVERERERERJzK3arBoaGhZoMGDawaLyIiIiIiIg60adOmE6ZphpX1mGVFtEGDBsTHx1s1XkRERERERBzIMIzDF3pMS3NFRERERETEqVRERURERERExKlUREVERERERMSpVERFRERERETEqVRERURERERExKlUREVERERERMSpVERFRERERETEqVRERURERERExKlUREVERERERMSpVERFRERERETEqVRERURERERExKlUREVExGnsdrvVEURERMQFqIiKiEiFW7VqFddffz0tWrTgmmuuoWbNmnh5eTFixAgANm/erFIqIiJSjblbHUBERKqOU6dOERwczOHDh7n//vtp2bIlfn5+Z/94eXlht9t57rnn+OSTT/D09CQgIAA3Nzero4uIiIgT6YyoiIhUCNM0GTJkCLt37+buu+9m1KhRtGzZkgYNGhAaGkqNGjUwDAObzcaiRYuIiIhg0qRJtGjRgq+++srq+CIiIuJEKqIiInJVEhISGDVqFEVFRSxfvpyYmJhyf+/48eOZMmUKzz//PC+//LIDU4qIiIgrUREVEZErdvz4cQYMGEDv3r3x8PDA09Pzsr7fMAxiY2NZvXo1M2bMYMKECQ5KKiIiIq5E14iKiMgVycvLY8iQIdx///088MADV/VaderUYdWqVfTp0wcvLy+eeuqpCkopIiIirkhFVERELpvdbmf06NE0b96c5557rkJes1atWqxatYoTJ05QUlKCzWbDMIwKeW0RERFxLVqaKyIil+2FF14gLS2NDz/8sELLYnh4OM2bN2f8+PFMnz69wl5XREREXIthmqYlgzt06GDGx8dbMltERK7Oli1bqFevHqGhoQ55/ZMnT+Lt7Y2Hhwfu7lq8IyIiUhkZhrHJNM0OZT2mM6IiIlJu69at4+mnn6Zt27YOK6EAISEheHp60rlzZ7Zs2eKwOSIiImINfcwsIiLl1qpVq8veGfdKubm5MW7cOO644w7i4+Px8fFxylwRERFxPJ0RFRGRcnnnnXc4duwY7du3d9rMO++8kzZt2mgXXRERkSpGRVRERC4pISGB8ePHExwc7NS5hmHw3nvvMX/+fBYuXOjU2SIiIuI4KqIiInJRpmkyduxYnn/+ecLCwpw+PygoiE8//ZS//vWvpKamOn2+iIiIVDwVURERuaj58+eTmJjIgw8+aFmGXr16MWHCBAoLCy3LICIiIhVHRVRERC6ooKCAxx9/nMmTJ+Ph4WFplnvvvZcaNWqwePFiS3OIiIjI1VMRFRGRC5o2bRrNmzenf//+VkcBID09nQ0bNlgdQ0RERK6SYZqmJYM7dOhgxsfHWzJbRETKp7CwkMzMTEuuDb2YxMREIiMjrY4hIiIiF2EYxibTNDuU9ZjOiIqISJnGjx/P/v37Xa6E5ufn07lzZ7Zs2WJ1FBEREblCKqIiIlKmDh06UK9ePatjnKdGjRr885//ZNy4cVi1qkdERESujoqoiIicJy4ujtjYWAICAqyOUqb77ruPtLQ0vv/+e6ujiIiIyBVQERURkXMcPHiQO+64g6KiIqujXJC7uztvvPEGTz75JAUFBVbHERERkcukIioiIueYPHky9957L/7+/lZHuah+/foRHR3NO++8Y3UUERERuUwqoiIiclZGRgafffYZjzzyiNVRymXixIl4enpaHUNEREQuk7vVAURExHV8+OGHDBo0yCU3KSpLTEwMMTEx7N27l+joaKvjiIiISDnpjKiIiACl9wx96623+Mc//mF1lMuSnZ3N3XffTX5+vtVRREREpJx0RlRERABIT09n9OjRtG3b1uool8Xf35+1a9daHUNEREQug86IiogIpmkSEBDAyy+/bHWUKxYbG8vu3butjiEiIiLloCIqIiLs3LmTAQMGWB3jihmGQWxsLC+++KLVUURERKQcDNM0LRncoUMHMz4+3pLZIiJyvvz8fGrUqGF1jCuWlZVF48aN+fnnn7VxkYiIiAswDGOTaZodynpMZ0RFRKq53bt388ILL1TqEgoQEBDAI488wksvvWR1FBEREbkEFVERkWpu0qRJGIZhdYwK8eijj/LDDz+wf/9+q6OIiIjIRaiIiohUYzk5OcyePZsHHnjA6igVIigoiIceeogZM2ZYHUVEREQuQrdvERGpxr755ht69OhBrVq1rI5SYZ577jk8PDysjiEiIiIXUa4zooZhDDAMY49hGPsNw3jmIs/raBhGiWEYt1RcRBERcZRPPvmEu+++2+oYFcrLy4tff/2VF154weooIiIicgGXLKKGYbgBU4CBQHPgdsMwml/gea8CSyo6pIiIVLyDBw+yfft2brrpJqujVLhGjRrRr18/q2OIiIjIBZTnjGgnYL9pmr+ZplkIzAKGlvG8R4BvgNQKzCciIg7i6+vLjBkz8PLysjpKhQsJCaF169YsWrTI6igiIiJShvIU0bpA4h++Tjpz7CzDMOoCw4H3Ky6aiIg4it1uJz09nUGDBl31a2XnF7E/NYeMvEKsujd1WfLz87nzzjtJSUmxOoqIiIj8SXk2KyprT/8//6YxGXjaNM2Si90CwDCMvwF/A6hfv355M4qISAVLSkpi7NixLF68+Ipv3ZJbUMwnaw/x/uoDZOcXA+DhZhDi60mYnxd1grx5pE8UreoFVmT0cgsLC2PkyJG88847TJgwwZIMIiIiUjbjUp9eG4ZxHfC/pmn2P/P1swCmab78h+cc5L+FNRTIA/5mmubcC71uhw4dzPj4+KtLLyIiTldQXMIXvxzhnRX7Sc8tpEdUKH1iwsk8XcTJ3EJO5haSnlvInuRssk4XMa5vNA/0aoybzfn3Kt23bx9du3bl0KFD+Pr6On2+iIhIdWYYxibTNDuU9Vh5zohuBKIMw2gIHAVGAXf88QmmaTb8w7BPgAUXK6EiImKd3Nxcbr31VubNm4e7e/nv4mWaJnM2JTFp2V6OZeTT/ppgXr2lNa3qln3GM/N0Ea8u2s3rS/awck8qk25rS2SIT0X9NcolKiqKHj16MH36dB5++GGnzhYREZELu+Q1oqZpFgMPU7ob7i5gtmmaOw3DeMAwjKpxB3QRkWrk22+/xWazXVYJBZiycj9PztlGsI8n79x+Le/e2e6CJRQg0NuDl4a35IXBzdl1LIsBb/7IN5uSnH4d6ZNPPskbb7xBcXGxU+eKiIjIhZXrtxDTNBcCC/90rMyNiUzTHHP1sURExFFmzJjBAw9c3ueIs+MT+U/cXga0rMULg5tjK+d1pYZhcGOr2rSNDGL8/AT+8fVWNieeYsLQlld8berluu666xg9ejRZWVmEhIQ4ZaaIiIhcXHl2zRURkSri8OHDbNmyhcGDB5f7e1buTuXZb7fTqWEIzw9qVu4S+kd1grx598523N4pkpnrj/DlhsRLf1MFGj9+PIBL7eorIiJSnamIiohUI5999hkjR44s971DtyRm8ODnv9IkzI9XRrTCw+3K3zbcbAaP9Imic8MQXpi3g62JGVf8Wldi0KBB7N+/36kzRUREpGwqoiIi1YRpmsyYMYMxY8aU6/mHTuRyz/QNBPt6MGlkG3y9Lu+a0rK42QwmDG1JTV8v/v75Jk7mFl71a5bXmjVriIqKcto8ERERuTAVURGRasI0Td588006dChzF/VzpGUXcNfHGzCBN0deS02/8p1BLY9AHw9eHtGKtOwCHpu1mRK7c5bLurm5MXbsWHbu3OmUeSIiInJhKqIiItXEjh076N279yU3CTJNkyfnbCU1K5+Jt7ahfs2Kv+VK8zoBPNGvKT/tO8Gby/ZW+OtfSM2aNZk4caLT5omIiEjZVERFRKqJSZMmcfDgwUs+b2lCCqv2pPFA78a0vMjtWa7W0LZ1uKl1bd5asZ8Vu1McNuePHnzwQebOncvx48edMk9ERETKpiIqIlJNTJ8+nRYtWlz0OacLSxg/P4HGYb7c2r6eQ/MYhsGT/ZvSNMKfsbO2kJKV79B5UHpG9Pbbb+fdd991+CwRERG5MBVREZFq4N///jfz5s275PPeW7WfoxmnebJ/U9yvYofc8qrh4caLw1uSV1jCpKXOWaL72GOP8cEHH3D69GmnzBMREZHzqYiKiFRxpmny4Ycf0qBBg4s+79CJXN7/8Tf6t4jg2vrBzgkH1A/x4eb29Zgdn8jelGyHz4uOjqZTp07MnDnT4bNERESkbCqiIiJV3KZNm/D09KRVq1YXfI5pmoyfvxOPM/f6dLZ7uzXEx9OdVxbtdsq8t956ixEjRjhlloiIiJxPRVREpIqbM2cOt95660V3y122K5WVe9L4a49GhPlX3K1ayivQx4O7u17Dit2prD1wwuHzGjVqxJEjR9i8ebPDZ4mIiMj5VERFRKow0zSZM2cOt9xyywWfk19Uwvj5O2kU6sttHRy7QdHF3NYhkogAL15euAu7E+4teujQIRITEx0+R0RERM6nIioiUoVt2bIF0zRp27btBZ/z3qoDJJ06zRNO2qDoQmp4uPFAr8ZsP5rF/G3HHD5v+PDhDB48mLy8PIfPEhERkXOpiIqIVGERERF8+OGHF1yWm5qVz3urD9CveQTtr3HeBkUXMqBLine truncated
|
||||||
|
"text/plain": [
|
||||||
|
"<Figure size 1152x432 with 1 Axes>"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
"metadata": {
|
||||||
|
"needs_background": "light"
|
||||||
|
},
|
||||||
|
"output_type": "display_data"
|
||||||
|
}
|
||||||
|
],
|
||||||
|
"source": [
|
||||||
|
"import matplotlib.pyplot as plt\n",
|
||||||
|
"import numpy as np\n",
|
||||||
|
"from scipy import stats\n",
|
||||||
|
"\n",
|
||||||
|
"np.random.seed(654)\n",
|
||||||
|
"# Draw samples from two Gaussian w.p. z_i ~ Bernoulli(phi)\n",
|
||||||
|
"generative_m = np.array([stats.norm(2, 1), stats.norm(5, 1.8)])\n",
|
||||||
|
"z = stats.bernoulli(0.75).rvs(100)\n",
|
||||||
|
"x_i = np.array([g.rvs() for g in generative_m[z_i]])\n",
|
||||||
|
"\n",
|
||||||
|
"# plot generated data and the latent distributions\n",
|
||||||
|
"x = np.linspace(-5, 12, 150)\n",
|
||||||
|
"plt.figure(figsize=(16, 6))\n",
|
||||||
|
"plt.plot(x, generative_m[0].pdf(x))\n",
|
||||||
|
"plt.plot(x, generative_m[1].pdf(x))\n",
|
||||||
|
"plt.plot(x, generative_m[0].pdf(x) + generative_m[1].pdf(x), lw=1, ls='-.', color='black')\n",
|
||||||
|
"plt.fill_betweenx(generative_m[0].pdf(x), x, alpha=0.1)\n",
|
||||||
|
"plt.fill_betweenx(generative_m[1].pdf(x), x, alpha=0.1)\n",
|
||||||
|
"plt.vlines(x_i, 0, 0.01, color=np.array(['C0', 'C1'])[z])\n",
|
||||||
|
"plt.show()"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "markdown",
|
||||||
|
"metadata": {},
|
||||||
|
"source": [
|
||||||
|
"# Log likelihood\n",
|
||||||
|
"Let $\\theta = \\{ \\mu, \\sigma \\}$.\n",
|
||||||
|
"$$\\ell(\\theta) = \\sum_{i-1}^n \\log P(x_i; \\theta) $$\n",
|
||||||
|
"\n",
|
||||||
|
"Or in vector form:\n",
|
||||||
|
"\n",
|
||||||
|
"$$\\ell(\\theta) = \\log P(x; \\theta) $$\n",
|
||||||
|
"\n",
|
||||||
|
"This log likelihood function cannot be solved directly as it not includes the latent variable $z$.\n",
|
||||||
|
"\n",
|
||||||
|
"$$ \\ell(\\theta) = \\sum_{k=1}^K \\log P(x | z_k; \\theta) P(z_k) $$\n",
|
||||||
|
"\n",
|
||||||
|
"## Continuous form\n",
|
||||||
|
"\n",
|
||||||
|
"$$ \\begin{eqnarray}\n",
|
||||||
|
"\\ell(\\theta) & = & \\int \\log P(x | \\theta, z) P(z) dz \\\\\n",
|
||||||
|
"\\log P(x; \\theta)& = & \\int \\log P(x | \\theta, z) P(z) dz\n",
|
||||||
|
"\\end{eqnarray}$$\n",
|
||||||
|
"\n",
|
||||||
|
"\n",
|
||||||
|
"\n"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "markdown",
|
||||||
|
"metadata": {},
|
||||||
|
"source": [
|
||||||
|
"**Jensens inequality:**\n",
|
||||||
|
"*Only for concave functions*\n",
|
||||||
|
"$$ E[f(X)] \\le f(E[X]) $$\n",
|
||||||
|
"\n",
|
||||||
|
"**Jensens equality:**\n",
|
||||||
|
"$$ E[f(X)] = f(E[X]) \\iff X = E[X] $$"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "markdown",
|
||||||
|
"metadata": {},
|
||||||
|
"source": [
|
||||||
|
"## Evidence Lower BOund\n",
|
||||||
|
"\n",
|
||||||
|
"\n",
|
||||||
|
"Multiply the **evidence** w/ an arbitrary distribution over $Z, \\quad \\frac{Q(z)}{Q(z)} = 1$ \n",
|
||||||
|
"\n",
|
||||||
|
"$$ \\begin{eqnarray}\n",
|
||||||
|
"\\log P(x; \\theta)& = & \\int \\log Q(z) \\frac{ P(x | \\theta, z) P(z)}{Q(z)} dz \\\\\n",
|
||||||
|
"& = & \\log E_{z \\sim Q}[ \\frac{ P(x | \\theta, z) P(z)}{Q(z)}] \\\\\n",
|
||||||
|
"\\end{eqnarray}$$\n",
|
||||||
|
"\n",
|
||||||
|
"As the log function is a concave function, we can apply Jensen's inequality.\n",
|
||||||
|
"\n",
|
||||||
|
"$$ \\begin{eqnarray}\n",
|
||||||
|
"\\log P(x; \\theta)& = & \\log E_{z \\sim Q}[ \\frac{ P(x | \\theta, z) P(z)}{Q(z)}] \\\\\n",
|
||||||
|
"& \\ge & E_{z \\sim Q}[ \\log \\frac{ P(x | \\theta, z) P(z)}{Q(z)}]\n",
|
||||||
|
"\\end{eqnarray}$$\n",
|
||||||
|
"\n",
|
||||||
|
"\n",
|
||||||
|
"\n",
|
||||||
|
"Resulting in a lower bound on the evidence called the **ELBO**:\n",
|
||||||
|
"\n",
|
||||||
|
"$$ \\begin{eqnarray}\n",
|
||||||
|
"\\log P(x; \\theta)& \\ge & E_{z \\sim Q}[ \\log \\frac{ P(x | \\theta, z) P(z)}{Q(z)}]\n",
|
||||||
|
"\\end{eqnarray}$$\n"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "markdown",
|
||||||
|
"metadata": {},
|
||||||
|
"source": [
|
||||||
|
"# Now what?\n",
|
||||||
|
"\n",
|
||||||
|
"Jensen's inequality holds to equaltiy if and only if $X = E[X]$. This will be the case if we choose a distribution that leads to:\n",
|
||||||
|
"\n",
|
||||||
|
"$$ \\begin{eqnarray}\n",
|
||||||
|
"\\frac{P(x|z) P(z)}{Q(z)} &=& \\frac{P(x, z)}{Q(z)} = 1\n",
|
||||||
|
"\\end{eqnarray}$$\n",
|
||||||
|
"\n",
|
||||||
|
"Wich means that \n",
|
||||||
|
"\n",
|
||||||
|
"$$ P(x, z; \\theta) \\propto Q(z) $$\n",
|
||||||
|
"\n",
|
||||||
|
"A next requirement is that $\\int Q(z) dz = 1$ (requirement of probability distribution).\n",
|
||||||
|
"\n",
|
||||||
|
"$$ \\begin{eqnarray}\n",
|
||||||
|
"\\int \\frac{P(x, z; \\theta)}{\\int P(x, z; \\theta) dz} dz &=& 1 \\\\\n",
|
||||||
|
"\\end{eqnarray}$$\n",
|
||||||
|
"\n",
|
||||||
|
"This will lead to:\n",
|
||||||
|
"\n",
|
||||||
|
"$$ \\begin{eqnarray}\n",
|
||||||
|
"Q(z) &=& \\frac{P(x, z; \\theta)}{\\int P(x, z; \\theta) dz} \\\\\n",
|
||||||
|
"&=& \\frac{P(x, z; \\theta)}{\\int P(x; \\theta)} \\\\\n",
|
||||||
|
"&=& P(z|x; \\theta)\n",
|
||||||
|
"\\end{eqnarray}$$"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "markdown",
|
||||||
|
"metadata": {},
|
||||||
|
"source": [
|
||||||
|
"## Gaussian mixture model\n",
|
||||||
|
"\n",
|
||||||
|
"### Expectation\n",
|
||||||
|
"Given a Gaussian mixture model, observed data points $x$ and current parameters $\\theta_t$ find $Q(z)$"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "code",
|
||||||
|
"execution_count": 377,
|
||||||
|
"metadata": {},
|
||||||
|
"outputs": [
|
||||||
|
{
|
||||||
|
"data": {
|
||||||
|
"text/plain": [
|
||||||
|
"array([ 4.28244273, 4.66972656, 3.35275075, 2.36158676, 3.74751768,\n",
|
||||||
|
" 9.75240814, 5.15756521, 4.80213993, 5.90145058, 1.91271217,\n",
|
||||||
|
" 4.71273852, 2.36241133, 1.61084671, 6.71697406, 3.11442326,\n",
|
||||||
|
" 5.15422491, 5.09580485, 2.11274369, 6.76574344, 1.48354695,\n",
|
||||||
|
" 2.65304672, 6.81892055, 2.77103106, 5.06418788, 5.72517362,\n",
|
||||||
|
" 3.54313443, 5.67858554, 1.03222311, 1.63983108, 5.51107055,\n",
|
||||||
|
" 0.65124788, 4.97843073, 4.32989217, 2.28812021, 4.51797044,\n",
|
||||||
|
" 4.52825271, 6.93340533, 4.94273792, 2.54969132, 0.78634279,\n",
|
||||||
|
" 4.81172184, 5.19576022, 6.11967179, 3.5235398 , 0.75307613,\n",
|
||||||
|
" 4.23428727, 7.15230773, 2.83345813, 8.60608829, 4.27073679,\n",
|
||||||
|
" 7.75388195, 6.79259227, 0.62404892, 4.26260727, 4.09584051,\n",
|
||||||
|
" 3.26993383, 4.86431429, 8.38977374, 6.48961063, 4.09783105,\n",
|
||||||
|
" 4.40498267, 3.27433241, 7.44164281, 7.62286463, 3.21797707,\n",
|
||||||
|
" 4.42538064, 7.06135769, 3.94733606, 1.91385988, 5.23226903,\n",
|
||||||
|
" 3.61125312, 1.86866663, 5.71568321, 2.50687193, 2.26842535,\n",
|
||||||
|
" 3.75523264, 0.71738199, 8.02274291, 2.46995812, 5.35222504,\n",
|
||||||
|
" 2.72716235, 4.15793385, 5.33299762, 1.00091742, 7.04351757,\n",
|
||||||
|
" 6.14452194, 2.62636128, 4.55251749, -0.38101307, 4.27776483,\n",
|
||||||
|
" 1.57442618, 3.75797097, 0.52259715, 1.59952851, 6.13345271,\n",
|
||||||
|
" 1.42972786, 2.82725842, 4.57300419, 2.59096887, 3.45567689])"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
"execution_count": 377,
|
||||||
|
"metadata": {},
|
||||||
|
"output_type": "execute_result"
|
||||||
|
}
|
||||||
|
],
|
||||||
|
"source": [
|
||||||
|
"# Data points x\n",
|
||||||
|
"x_i"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "code",
|
||||||
|
"execution_count": 266,
|
||||||
|
"metadata": {},
|
||||||
|
"outputs": [],
|
||||||
|
"source": [
|
||||||
|
"import torch\n",
|
||||||
|
"from torch import nn\n",
|
||||||
|
"import torch.distributions as dist\n",
|
||||||
|
"import torch.nn.functional as F"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "code",
|
||||||
|
"execution_count": 363,
|
||||||
|
"metadata": {},
|
||||||
|
"outputs": [
|
||||||
|
{
|
||||||
|
"name": "stderr",
|
||||||
|
"output_type": "stream",
|
||||||
|
"text": [
|
||||||
|
"/opt/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:1: UserWarning: To copy construct from a tensor, it is recommended to use sourceTensor.clone().detach() or sourceTensor.clone().detach().requires_grad_(True), rather than torch.tensor(sourceTensor).\n",
|
||||||
|
" \"\"\"Entry point for launching an IPython kernel.\n"
|
||||||
|
]
|
||||||
|
}
|
||||||
|
],
|
||||||
|
"source": [
|
||||||
|
"x_i = torch.tensor(x_i, dtype=torch.float)"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "code",
|
||||||
|
"execution_count": 364,
|
||||||
|
"metadata": {},
|
||||||
|
"outputs": [
|
||||||
|
{
|
||||||
|
"data": {
|
||||||
|
"text/plain": [
|
||||||
|
"<BarContainer object of 100 artists>"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
"execution_count": 364,
|
||||||
|
"metadata": {},
|
||||||
|
"output_type": "execute_result"
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"data": {
|
||||||
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA6IAAAHjCAYAAAA5Y6JnAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAbMElEQVR4nO3df7Dld13f8dfbXULk14aa1YHdyKZOQKJAwCXSagsa1A1aUzrqJOGHUmiMJRZrZ8r2JzvFtiJK0TEQI00DrSXjIC0rWUgdEBlwgGwAA0sauoQtWUPN8sMUQSdG3v3jnqWHm5u9Zzfnfs7eex+PmTt7v9/v55zzPpczWZ77/Z5zq7sDAAAAo3zDogcAAABgcxGiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQA2har691X1czOs+1BVfceImSaPN+tc+6pq3xrOMfR5A7C5CVEA1r2qenhV/UJVfaqqvlRVn6iqn546vj3JC5P8xgx398tJ/s2c5np0VXVV/VlVfaWq7pqOzpOcay6q6qlV9f7JPB+qqm+dHJrb8waA1QhRANa1qnp0kvclOTfJRUkeleQfJHllVb14suynkhzo7j+f4S73J/m+qnrMHMa7IMmx7n5Edz8syc8k+Q9VtfMU5nrQJo97IMmrknxTkjuS/MvJ4Xk+bwA4ISEKwHr3q0n+JMnzu/tIL3l/ktckeelkzcVJ/mD6RlX125Mzlce/uqqu6u6/SHJLkh+cw2wXJLl5avuDkz/PWGmuB5pppTuuql+qqv82tf3qqnpXVT2kql5XVa9b4Wa/kuQ3u3v/JH5vSPL0JJnz8waAE9q66AEA4FRV1eOSPC/J07u7lx3+VJJdk++flOT26YPd/RNT9/PzSV6Q5L9Odt2W5CkrPN7bk3zvA4zzvu7+kWX7nprkQ5PbnpXk32Yp9j690lyrzLTcq5J8qqouSPKMJHuSfG93/2WSf7jC7I9KckmSx0/t/oYkfzG1veLzBoB5E6IArGfPTnJnd394hWM7khydfH9Wki+tdAdV9bIsvU/zou7+wmT3l5Lc7xLVFUJzNRck+XuTx/hikt9P8nemonnFuR5gpuWzfL6qXpvkTUm2ZSlC7znBLBcleUiSW6vq+L6HJnnb1JoVnzcAzJsQBWA9257/H5vLPTfJuyfffzHJI5cvmFz2+uIk39/dn5869Mgkf/pgBquqhyZ5YpJzu/uBZrzfXCeYaSUfSfKKJM/r7jtXWbsryf7ufu7UYx1I8s6pNQ/6eQPALLxHFID17NNJHldVX/f3WVX9QJLdWXpPZJLcmq+/JDVV9TNJfjpLZx0/t+x+n5jkj5Y/WFW9Y9l7OKe/3rFs+Xcm+fIJIvR+c60y0/JZnpTk9UnemOTvn2jtxEOTfGXq9udm6We0f2rNis8bAOZNiAKwnt04+fMXquphVfXQqnp+kjcn+fGps4QHkjzz+I2q6ookVyV5dncfm77DyZnM70rye8sfrLsvnnwC7kpfFy9b/tQkh1aZ/2tznWim5apqR5LfTXJllt4P+qSqetbU8eur6vplN7s5yTOr6rFVdU6W3nv6L45f+nui5w0A8yZEAVi3uvvPsvTexyclOZLkz5P8fJJndveBqaVvSvKcqvrGyfYvJfm2LH3Yz/Ezmi+YHPvRJO/p7rse5HgXJPn4Kmum5zrRTF8z+dChA0leM/n0268keXWWPgjpuHOSvH/ZTd+dpXj9ZJZ+3c1/7u7fnDo+r+cNAKuq+3/IIACsT1X1E1n6dS7fsfxDfqrq3yW5u7tfu8p9fDDJi7t7tYici5OYa1+SdPe+VdadkaXLa588+QTdWecY+rwB2NyEKAAbyuR9loe6+72LnmWeZg1RAFgPhCgArAPH3wPa3e9Z7CQA8OAJUQAAAIbyYUUAAAAMtXVRD3z22Wf3rl27FvXwAAAArKFbbrnlc929faVjCwvRXbt25eDBg4t6eAAAANZQVf3vBzrm0lwAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKFWDdGquq6q7q6qjz/A8aqqX6uqw1V1a1U9bf5jAgAAsFHMckb0+iR7TnD84iTnTb6uSPL6Bz8WAAAAG9WqIdrd703yhRMsuSTJm3rJB5KcVVWPmdeAAAAAbCzzeI/ojiR3Tm0fnewDAACA+5lHiNYK+3rFhVVXVNXBqjp47NixOTz0Gtu3bdETnN5W+/kcP34yP8d5/Myn7+N0+9/wRPNM/7xO5Wd3Mo/9YO/3VG4//bxO5f4eaP5Tvc/lx1faXqvXz8k+31lvczKPtxavMdhkdu29cdEjnJIHmvv4/lme16xrFvkzOtXHPpnbrcXzW6+vq2R9z8548wjRo0nOmdremeSulRZ297Xdvbu7d2/fvn0ODw3A3J3sPx6JWTaZRQcWwEYwjxDdn+SFk0/PfUaSe7r7s3O4XwCAuRGPAKePrastqKo3J3lWkrOr6miSVyR5SJJ09zVJDiR5TpLDSb6S5EVrNSwAAADr36oh2t2XrXK8k7x0bhMBsHns25bsu2fRUwAAg83j0lwAAACYmRAFANYV7/UEWP+EKAAAAEMJUQAAAIYSogCsH35nKQBsCEIUAACAoYQoAAAAQwlRANY3l+sCwLojRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAG5NP0z0t7dp746JHAOA0IEQBAAAYSogCsLE5MwoApx0hCsDGJ0YB4LQiRAHYHMQoAJw2hCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCsDn58CIAWBghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgATNu3bdETAMCGJ0QBAAAYSogCACdt194bFz0CAOuYEAUAAGAoIQoAy3mfKACsKSEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAcD+79t646BEA2MCEKAAAAEMJUQA4Wfu2LXoCAFjXhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAEiS7Np746JHAGCTEKIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoADxY+7YtegIAWFdmCtGq2lNVt1fV4arau8LxbVX1u1X1R1V1qKpeNP9RAQAA2AhWDdGq2pLk6iQXJzk/yWVVdf6yZS9N8onufkqSZyX5lao6Y86zAgAAsAHMckb0wiSHu/uO7r43yQ1JLlm2ppM8sqoqySOSfCHJfXOdFAAAgA1hlhDdkeTOqe2jk33Tfj3JE5PcleRjSV7W3V9dfkdVdUVVHayqg8eOHTvFkQEAAFjPZgnRWmFfL9v+oSQfTfLYJBck+fWqetT9btR9bXfv7u7d27dvP+lhAYD52rX3xkWPAMAmNEuIHk1yztT2ziyd+Zz2oiRv7SWHk3w6ybfPZ0QAAAA2kllC9OYk51XVuZMPILo0yf5laz6T5KIkqapvSfKEJHfMc1AAAAA2hq2rLeju+6rqqiQ3JdmS5LruPlRVV06OX5PklUmur6qPZelS3pd39+fWcG4AAADWqVVDNEm6+0CSA8v2XTP1/V1JfnC+owEAALARzXJpLgAAAMyNEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBYBNYtfeGxc9AgAkEaIAAAAMJkQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAKy1fdsWPQEAnFaEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiALABrZr742LHgEA7keIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAjLRv26InAICFE6IAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQwlRAAAAhhKiAAAADCVEAQAAGEqIAgAAMJQQBQAAYCghCgAAwFBCFAAAgKGEKAAAAEMJUQAAAIYSogAAAAwlRAEAABhKiAIAADCUEAUAAGAoIQoAAMBQQhQAAIChhCgAAABDCVEAAACGEqIAAAAMJUQBAAAYSogCAAAwlBAFAABgKCEKAADAUEIUAACAoYQoAAAAQ80UolW1p6pur6rDVbX3AdY8q6o+WlWHquoP5jsmAAAAG8XW1RZU1ZYkVyf5gSRHk9xcVfu7+xNTa85K8roke7r7M1X1zWs1MAAAAOvbLGdEL0xyuLvv6O57k9yQ5JJlay5P8tbu/kySdPfd8x0TAACAjWKWEN2R5M6p7aOTfdMen+TRVfWeqrqlql640h1V1RVVdbCqDh47duzUJgYAAGBdmyVEa4V9vWx7a5LvSvLDSX4oyb+qqsff70bd13b37u7evX379pMeFgAAgPVv1feIZukM6DlT2zuT3LXCms9195eTfLmq3pvkKUk+OZcpAQAA2DBmOSN6c5LzqurcqjojyaVJ9i9b87Ykf6uqtlbVw5J8d5Lb5jsqAAAAG8GqZ0S7+76quirJTUm2JLmuuw9V1ZWT49d0921V9c4ktyb5apI3dPfH13JwAAAA1qdZLs1Ndx9IcmDLine truncated
|
||||||
|
"text/plain": [
|
||||||
|
"<Figure size 1152x576 with 1 Axes>"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
"metadata": {
|
||||||
|
"needs_background": "light"
|
||||||
|
},
|
||||||
|
"output_type": "display_data"
|
||||||
|
}
|
||||||
|
],
|
||||||
|
"source": [
|
||||||
|
"torch.manual_seed(1)\n",
|
||||||
|
"\n",
|
||||||
|
"k = torch.tensor([2])\n",
|
||||||
|
"mu = dist.Uniform(x_i.min(), x_i.max()).sample(k) \n",
|
||||||
|
"std = torch.ones(k)\n",
|
||||||
|
"\n",
|
||||||
|
"def expectation(x, mu, std, k):\n",
|
||||||
|
" qz = torch.zeros((k, x.shape[0]))\n",
|
||||||
|
" \n",
|
||||||
|
" for z in range(k):\n",
|
||||||
|
" qz[z] = dist.Normal(mu[z], std[z]).log_prob(x).exp() \n",
|
||||||
|
" \n",
|
||||||
|
" # normalize zo that marginalizing z would lead to p = 1\n",
|
||||||
|
" qz = qz / qz.sum(0)\n",
|
||||||
|
" return qz\n",
|
||||||
|
"\n",
|
||||||
|
"qz = expectation(x_i, mu, std, k) \n",
|
||||||
|
"\n",
|
||||||
|
"plt.figure(figsize=(16, 8)) \n",
|
||||||
|
"plt.title(r\"$Q(z) = P(z|x; \\theta)$\")\n",
|
||||||
|
"plt.bar(x_i, qz[0], width=0.01)\n",
|
||||||
|
"plt.bar(x_i, qz[1], width=0.01)"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "markdown",
|
||||||
|
"metadata": {},
|
||||||
|
"source": [
|
||||||
|
"### Maximization\n",
|
||||||
|
"\n",
|
||||||
|
"Optimize $\\theta$ with the $Q(z)$ found in the expectation step.\n",
|
||||||
|
"\n",
|
||||||
|
"$$ \\sum_{z_i} Q(z_i) \\log ( \\frac{p(x_i| z_i; \\theta) p(z_i)} {Q(z_i)}) $$"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "code",
|
||||||
|
"execution_count": 370,
|
||||||
|
"metadata": {},
|
||||||
|
"outputs": [
|
||||||
|
{
|
||||||
|
"name": "stderr",
|
||||||
|
"output_type": "stream",
|
||||||
|
"text": [
|
||||||
|
"/opt/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:4: UserWarning: To copy construct from a tensor, it is recommended to use sourceTensor.clone().detach() or sourceTensor.clone().detach().requires_grad_(True), rather than torch.tensor(sourceTensor).\n",
|
||||||
|
" after removing the cwd from sys.path.\n",
|
||||||
|
"/opt/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:6: UserWarning: To copy construct from a tensor, it is recommended to use sourceTensor.clone().detach() or sourceTensor.clone().detach().requires_grad_(True), rather than torch.tensor(sourceTensor).\n",
|
||||||
|
" \n"
|
||||||
|
]
|
||||||
|
}
|
||||||
|
],
|
||||||
|
"source": [
|
||||||
|
"k = torch.tensor([2])\n",
|
||||||
|
"mu = dist.Uniform(x_i.min(), x_i.max()).sample(k) \n",
|
||||||
|
"std = torch.ones(k)\n",
|
||||||
|
"mu = torch.tensor(mu, dtype=torch.float, requires_grad=True)\n",
|
||||||
|
"std = torch.tensor([1., 1.])\n",
|
||||||
|
"std = torch.tensor(std, dtype=torch.float, requires_grad=True)\n",
|
||||||
|
"phi = torch.rand(2, requires_grad=True)\n"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "code",
|
||||||
|
"execution_count": 371,
|
||||||
|
"metadata": {},
|
||||||
|
"outputs": [
|
||||||
|
{
|
||||||
|
"name": "stderr",
|
||||||
|
"output_type": "stream",
|
||||||
|
"text": [
|
||||||
|
"/opt/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:4: UserWarning: Implicit dimension choice for softmax has been deprecated. Change the call to include dim=X as an argument.\n",
|
||||||
|
" after removing the cwd from sys.path.\n"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"data": {
|
||||||
|
"text/plain": [
|
||||||
|
"tensor(-512.5377, grad_fn=<AddBackward0>)"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
"execution_count": 371,
|
||||||
|
"metadata": {},
|
||||||
|
"output_type": "execute_result"
|
||||||
|
}
|
||||||
|
],
|
||||||
|
"source": [
|
||||||
|
"def get_elbo(mu, std, phi):\n",
|
||||||
|
" elbo = 0\n",
|
||||||
|
" for k in range(2):\n",
|
||||||
|
" z_i = F.softmax(phi)[k]\n",
|
||||||
|
" mu_i = mu[k]\n",
|
||||||
|
" std_i = std[k]\n",
|
||||||
|
" \n",
|
||||||
|
" for i in range(len(x_i)):\n",
|
||||||
|
" elbo += qz[k, i] * ((dist.Normal(mu_i, std_i).log_prob(x_i[i]) + torch.log(z_i) - torch.log(qz[k, i])))\n",
|
||||||
|
" return elbo\n",
|
||||||
|
"\n",
|
||||||
|
"elbo = get_elbo(mu, std, phi)\n",
|
||||||
|
"elbo"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "code",
|
||||||
|
"execution_count": 372,
|
||||||
|
"metadata": {},
|
||||||
|
"outputs": [
|
||||||
|
{
|
||||||
|
"name": "stderr",
|
||||||
|
"output_type": "stream",
|
||||||
|
"text": [
|
||||||
|
"/opt/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:4: UserWarning: Implicit dimension choice for softmax has been deprecated. Change the call to include dim=X as an argument.\n",
|
||||||
|
" after removing the cwd from sys.path.\n",
|
||||||
|
"/opt/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:21: UserWarning: Implicit dimension choice for softmax has been deprecated. Change the call to include dim=X as an argument.\n"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"name": "stdout",
|
||||||
|
"output_type": "stream",
|
||||||
|
"text": [
|
||||||
|
"215.0631103515625 tensor([5.8801, 2.9457]) tensor([1.5491, 1.5054]) tensor([0.3965, 0.6035])\n",
|
||||||
|
"214.99671936035156 tensor([5.6912, 2.8268]) tensor([1.5984, 1.4784]) tensor([0.4478, 0.5522])\n",
|
||||||
|
"214.94085693359375 tensor([5.5913, 2.7706]) tensor([1.6360, 1.4637]) tensor([0.4746, 0.5254])\n",
|
||||||
|
"214.90855407714844 tensor([5.5367, 2.7465]) tensor([1.6611, 1.4576]) tensor([0.4884, 0.5116])\n",
|
||||||
|
"214.89120483398438 tensor([5.5059, 2.7380]) tensor([1.6776, 1.4564]) tensor([0.4954, 0.5046])\n",
|
||||||
|
"214.88137817382812 tensor([5.4879, 2.7368]) tensor([1.6889, 1.4577]) tensor([0.4989, 0.5011])\n",
|
||||||
|
"214.87564086914062 tensor([5.4766, 2.7388]) tensor([1.6970, 1.4600]) tensor([0.5006, 0.4994])\n",
|
||||||
|
"214.8721160888672 tensor([5.4690, 2.7420]) tensor([1.7031, 1.4625]) tensor([0.5014, 0.4986])\n",
|
||||||
|
"214.8697509765625 tensor([5.4635, 2.7455]) tensor([1.7079, 1.4651]) tensor([0.5018, 0.4982])\n",
|
||||||
|
"214.8681182861328 tensor([5.4594, 2.7489]) tensor([1.7118, 1.4674]) tensor([0.5019, 0.4981])\n"
|
||||||
|
]
|
||||||
|
}
|
||||||
|
],
|
||||||
|
"source": [
|
||||||
|
"lr = 0.01\n",
|
||||||
|
"for i in range(10):\n",
|
||||||
|
" \n",
|
||||||
|
" # expectation step\n",
|
||||||
|
" with torch.no_grad():\n",
|
||||||
|
" qz = expectation(x_i, mu, std, k)\n",
|
||||||
|
" \n",
|
||||||
|
" # maximization step\n",
|
||||||
|
" for j in range(150):\n",
|
||||||
|
" negative_elbo = -get_elbo(mu, std, phi)\n",
|
||||||
|
" negative_elbo.backward()\n",
|
||||||
|
" \n",
|
||||||
|
" with torch.no_grad():\n",
|
||||||
|
" mu -= mu.grad * lr \n",
|
||||||
|
" std -= std.grad * lr\n",
|
||||||
|
" phi -= phi.grad * lr\n",
|
||||||
|
" std.grad.zero_()\n",
|
||||||
|
" mu.grad.zero_()\n",
|
||||||
|
" phi.grad.zero_()\n",
|
||||||
|
"\n",
|
||||||
|
" print(negative_elbo.item(), mu.data, std.data, F.softmax(phi).data)"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "code",
|
||||||
|
"execution_count": 384,
|
||||||
|
"metadata": {},
|
||||||
|
"outputs": [
|
||||||
|
{
|
||||||
|
"data": {
|
||||||
|
"text/plain": [
|
||||||
|
"[<matplotlib.lines.Line2D at 0x7f7ac0160470>]"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
"execution_count": 384,
|
||||||
|
"metadata": {},
|
||||||
|
"output_type": "execute_result"
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"data": {
|
||||||
|
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA6gAAAFlCAYAAADmqMVrAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdeZzO5f7H8dc1M8aIscQgxp49+9gisowtJSRbVDpFUTmd+knldFJJnTpp0SLFsZRUlN2MQgkxdlmHLMNgxpYZY9br98etjmowY+77/t4z834+Hh7MfX+/1+c95WHmM9f1vS5jrUVERERERETEaX5OBxAREREREREBNagiIiIiIiLiI9SgioiIiIiIiE9QgyoiIiIiIiI+QQ2qiIiIiIiI+AQ1qCIiIiIiIuITApwOkJlSpUrZypUrOx1DRERERERE3GzDhg3x1tqQzN7zyQa1cuXKREVFOR1DRERERERE3MwYc/By72mJr4iIiIiIiPgENagiIiIiIiLiE9SgioiIiIiIiE9QgyoiIiIiIiI+QQ2qiIiIiIiI+AQ1qCIiIiIiIuIT1KCKiIiIiIiIT1CDKiIiIiIiIj4hSw2qMaaLMWa3MSbaGPP0Fa5raoxJN8bcld17RUREREREJH+7aoNqjPEHJgJdgTpAf2NMnctc9yqwNLv3ioiIiIiIiGRlBrUZEG2t3W+tTQFmAT0yue5R4CvgxDXcKyIiIiIiIvlcVhrU8sDhSz6Oufja74wx5YGewAfZvfeSMR4yxkQZY6Li4uKyEEtERERERETykqw0qCaT1+yfPp4AjLLWpl/Dva4XrZ1krQ2z1oaFhIRkIZaIiPiqCxcu8Oabb5KRkcELL7xAvXr16Ny5M0OGDGHMmDHs3LnT6YgiIiLigwKycE0MUOGSj0OBo3+6JgyYZYwBKAV0M8akZfFeERHJI6y1HD58mNDQUBITE7lw4QKPP/44PXv25MiRIxw5coTo6GjatGnDlClT6N69u9ORRURExIcYazOd0PzfBcYEAHuADsARYD0wwFr782WunwossNZ+md17fxMWFmajoqKy+amIiIiTtm/fzvDhwwkODmbBggVXvHb9+vWkp6cTFhaGn58ffn469UxERCS/MMZssNaGZfbeVb8jsNamASNw7c67E5htrf3ZGDPMGDPsWu7N7icgIiK+beLEibRv356+ffvyzTffXPX6pk2b0qJFC9577z1eeeUVLyQUERGR3OCqM6hO0AyqiEjusXz5cgYMGMDatWupVKlStu5NSUnhwoULnD17lvPnz1OzZk0PpRQRERFfkaMZVBERkcs5cuQIAwcOZPr06dluTgECAwMpWrQoUVFRtGvXDu3iLiIikr+pQRURkWtWtGhR3nvvPTp27JijcXr27MnAgQN5+OGH8cWVPSIiIuIdalBFROSa/Pvf/+bUqVPceeedbhnvxRdfZMeOHcyaNcst44mIiEjuowZVRESuSZkyZShevLjbxgsKCmLatGmMHDmS2NhYt40rIiIiuYcaVBERyZatW7cye/ZsBg8eTLFixdw6dlhYGEOHDmX06NFuHVdERERyBzWoIiKSZcnJydx1112kpKR4rMZzzz3HG2+84bHxRURExHepQRURkSybNGkSNWrU4J577vFYjcDAQEqUKEHPnj05evSox+qIiIiI7wlwOoCIiOQOiYmJjBs3jsWLF3u8lp+fH8888wxly5b1eC0RERHxHZpBFRGRLHn77bdp06YNDRs29Eq9pk2bEhkZycKFC71ST0RERJynBlVERK7q3Llz/Oc//2Hs2LFerWutZdSoUWRkZHi1roiIiDhDDaqIiFxVkSJFWLlyJTVr1vRq3c6dOxMUFMS8efO8WldEREScoQZVRESu6MSJE7zzzjvUqVPH67WNMTzzzDOMGzcOa63X64uIiIh3qUEVEZErSk5OJjg42LH6d955JwkJCSxbtsyxDCIiIuIdalBFROSyjh8/TmBgIPfff79jGfz8/Bg9ejQvv/yyYxlERETEO3TMjIiIXNaYMWMICQlxvDns378/pUuXxlqLMcbRLCIiIuI5mkEVEZFM7d27lzlz5vCPf/zD6SgEBARw6623Mn/+fKejiIiIiAepQRURkUxNmDCBhx9+mOuvv97pKIDryJk5c+aQkpLidBQRERHxEC3xFRGRv0hISOCzzz5j69atTkf5XVBQEFOnTiUpKcnpKCIiIuIhmkEVEZG/mDVrFrfccguhoaFOR/mD9PR06tWrx+7du52OIiIiIh6gBlVERP4iPj6e4cOHOx3jL/z9/Rk8eDCvvvqq01FERETEA4wvHnweFhZmo6KinI4hIpIv+fpOuXFxcVSvXp2DBw9SrFgxp+OIiIhINhljNlhrwzJ7TzOoIiLyByNHjuSrr75yOsZlhYSE0L59e2bPnu10FBEREXEzNagiIvIHY8eOpWPHjk7HuKIhQ4bwySefOB1DRERE3EwNqoiI/C4yMpJffvnF55fOdunShZiYGA4fPux0FBEREXEjNagiIgK4nj196qmniIuLczrKVQUEBLBnzx4qVKjgdBQRERFxIzWoIiICwLp16zh37hwdOnRwOkqWFChQgCeffJK0tDSno4iIiIibqEEVEREAPvjgA4YOHYqfX+740hAQEEDt2rVJTk52OoqIiIi4Se74LkRERDzqzJkzzJ07l/vuu8/pKNnywAMP5IolySIiIpI1WWpQjTFdjDG7jTHRxpinM3m/hzFmqzFmszEmyhjT+pL3Dhhjtv32njvDi4iIe0yfPp0uXbpQunRpp6NkS0JCAg0bNlSTKiIikkdctUE1xvgDE4GuQB2gvzGmzp8u+xZoYK1tCAwBJv/p/XbW2oaXO4xVREScFRISwt///nenY2RbkSJF6NGjBzNmzHA6ioiIiLhBVmZQmwHR1tr91toUYBbQ49ILrLUJ1lp78cPCgEVERHKFCxcu0LdvX5o3b+50lGty//33M2XKFP73ZUhERERyq6w0qOWBSw+ai7n42h8YY3oaY3YBC3HNov7GAhHGmA3GmIdyElZERNxv4sSJjB071ukY16xNmzYkJiayYcMGp6OIiIhIDgVk4RqTyWt/+TG1tXYuMNcY0wZ4Eeh48a1W1tqjxpjSQKQxZpe19vu/FHE1rw8BVKxYMav5RUQkh5544olcvROun58f9913H1OmTCEsTE+SiIiI5GZZmUGNAS49CT0UOHq5iy82n9WMMaUufnz04u8ngLm4lgxndt8ka22YtTYsJCQki/FFRCQndu7cyccff0xQUJDTUXLk3nvv5ccff9QyXxERkVwuKw3qeqC6MaaKMSYQ6AfMu/QCY8yNxhhz8c+NgUDgpDGmsDEm+OLrhYFOwHZ3fgIiInLtpk2bxt69e52OkWMVK1Zk06ZNXPxSJCIiIrnUVRtUa20aMAJYCuwEZltrfzbGDDPGDLt4WW9guzFmM64df/te3DSpDLDKGLMFWAcstNYu8cQnIiIi2WOt5fPPP6dv375OR3GLn3/+mVdeecXpGCIiIpIDWXkGFWvtImDRn1774JI/vwq8msl9+4EGOcwoIiIesGHDBvz9/WnUqJHTUdyiTJkyeeZzERERya+yssRXRETyoN9mT/PKstiQkBBatWrF0aOX3SZBREREfJwaVBGRfMhay+zZs/PM8t7fTJkyhWeffdbpGCIiInKN1KCKiORDCQkJDBgwgJtuusnpKG7Vq1cv5s2bR2pqqtNRRERE5BqoQRURyYeCgoJ45ZVX8szy3t+EhoZSvXp1li9f7nQUERERuQZqUEVE8hlrLQ0bNuTIkSNOR/GI3r1789VXXzkdQ0RERK6BGlQRkXzGGMNPP/1E+fLlnY7iEb1792bu3Lmkp6c7HUVERESySQ2qiEg+M2PGjDy9023VqlWpXLkyW7dudTqKiIiIZJMaVBGRfCQ9PZ0nn3wSP7+8/c//jz/+qDNRRUREcqG8/R2KiIj8wQ8//ED58uW58cYbnY7iUdZaRo0aRUZGhtNRREREJBvUoIqI5CPffPMNPXv2dDqGxwUGBlK5cmWSk5OdjiIiIiLZoAZVRCSfsNYyb948br/9dqejeMWwYcNISEhwOoaIiIhkgxpUEZF8YufOnaSmplK/fn2no3jF7t27CQsLw1rrdBQRERHJIjWoIiL5xK+//soTTzyBMcbpKF5Rs2ZNChUqRFRUlNNRREREJIsCnA4gIiLe0aJFC1q0aOF0DK8xxtC7d2+++uormjZt6nQcERERyQLNoIqI5APx8fE0adIk3y137dGjB/Pnz3c6hoiIiGSRGlQRkXygZMmSfPXVV/lmee9vwsLCiI+P55dffnE6ioiIiGSBGlQRkXxg/vz5FCtWzOkYXufn50fXrl1ZsmSJ01FEREQkC4wvLvcKCwuz2tRCRMQ9kpOTKVOmDHv27KF06dJOx/G6M2fOEBwcjL+/v9NRREREBDDGbLDWhmX2nmZQRUTyuJUrV1KnTp182ZwCFC9enKlTp5KYmOh0FBEREbkKNagiInnc/Pnzuf32252O4ajjx49z9uxZp2OIiIjIVeiYGRGRPMxay/z581mwYIHTURz1zDPPkJGR4XQMERERuQrNoIqI5GHbt2/HGEPdunWdjuKoI0eOULt27Xx3zI6IiEhuowZVRCQP279/P4MGDcp3x8v8Wbly5UhLS2Pr1q1ORxEREZErUIMqIpKH9ejRg7Fjxzodw3HGGG677TYWLlzodBQRERG5AjWoIiJ51LFjxxg8eLDTMXzGbbfdlu+fxRUREfF1alBFRPKo4OBgHnzwQadj+Iy2bdvy888/Ex8f73QUERERuQw1qCIiedSePXto3ry50zF8RlBQEO3atWPJkiVORxEREZHL0DEzIiJ5UGpqKrfeeivR0dGEhIQ4HcdnvPrqq5QoUcLpGCIiInIZWZpBNcZ0McbsNsZEG2OezuT9HsaYrcaYzcaLine truncated
|
||||||
|
"text/plain": [
|
||||||
|
"<Figure size 1152x432 with 1 Axes>"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
"metadata": {
|
||||||
|
"needs_background": "light"
|
||||||
|
},
|
||||||
|
"output_type": "display_data"
|
||||||
|
}
|
||||||
|
],
|
||||||
|
"source": [
|
||||||
|
"fitted_m = [stats.norm(mu, std) for mu, std in zip(mu.data, std.data)]\n",
|
||||||
|
"\n",
|
||||||
|
"plt.figure(figsize=(16, 6))\n",
|
||||||
|
"plt.vlines(x_i, 0, 0.01, color=np.array(['C0', 'C1'])[z])\n",
|
||||||
|
"plt.plot(x, fitted_m[0].pdf(x))\n",
|
||||||
|
"plt.plot(x, fitted_m[1].pdf(x))\n",
|
||||||
|
"plt.plot(x, generative_m[0].pdf(x), color='black', lw=1, ls='-.')\n",
|
||||||
|
"plt.plot(x, generative_m[1].pdf(x), color='black', lw=1, ls='-.')\n"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"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