{ "metadata": { "name": "" }, "nbformat": 3, "nbformat_minor": 0, "worksheets": [ { "cells": [ { "cell_type": "code", "collapsed": false, "input": [ "import matplotlib.pyplot as plt\n", "import numpy as np\n", "\n", "from matplotlib import animation\n", "from JSAnimation import IPython_display" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 1 }, { "cell_type": "code", "collapsed": false, "input": [ "np.random.seed(13)\n", "data = np.random.random(100)\n", "plt.hist(data, bins=15, normed=True, color='black', alpha=0.5)\n", "plt.title('Histogram of $U(0,1)$ samples')" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "pyout", "prompt_number": 2, "text": [ "" ] }, { "metadata": {}, "output_type": "display_data", "png": "iVBORw0KGgoAAAANSUhEUgAAAs4AAAICCAYAAADfxkvjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xt0lNW9//HPTABJIOEihCQSSCItHCxQQqAQuaQtpo3I\nZS0uIlQMYiOpisixQnqsEC4Kpc1RqESP2ICoFUWPthKRcCTECKeEy7gARQMpgR4yUaACBQM/kv37\nw8UshoTMDplkArxfa80ys2fv5/k+yQ58eNyzx2GMMQIAAABQK2egCwAAAACuBQRnAAAAwALBGQAA\nALBAcAYAAAAsEJwBAAAACwRnAAAAwALBGQAAALBAcAYAAAAsEJwBVLNq1So5nU4dPnw40KVc93bu\n3KkhQ4YoLCxMTqdTBQUFfj/H9u3btXnzZr8f91JLly7VjfB5WvxuADc2gjNwA7j4l/327dtrfH36\n9OlyOr3/OHA4HHU+T2FhoTIzM3Xy5MmrqvNGU1lZqbvvvltlZWX6/e9/r1dffVU9evTw6zkOHTqk\nN954Qz/+8Y/9etzL/fjHP9acOXMa9BwAEGgEZwCSvIPylClT9O2336pLly51OgbBuW6OHDmikpIS\nzZgxQ2lpaZo0aZLCw8Ov2P9Pf/qThg8fLqfTqR49euiee+6R2+2WJH322WcaNmyYnE6nevXqpays\nLEnSb37zG/32t7/1Os7SpUv1yCOP6KmnntLkyZP11VdfWdf88ccfq2fPntXaExISZIzRxx9/bH0s\nALjWNAt0AQCahkv/N7vT6VSLFi38cqzGdPbsWYWEhATk3FfjYmBt06aNVf/7779fPXv2VGJiot54\n4w398Ic/9LzWs2dPbdmyRSkpKfrggw8kSUVFRQoLC1O7du08/V544QVt2LBB//M//yNJevHFFzV2\n7Fifgff111/Xpk2bdPr0ae3fv7/GPo888ogeeOABffjhh1bXAwDXGu44A6impnWc//rXv/T4448r\nNjZWwcHBioiI0E9/+lNt2bJFkjRv3jz95je/kSTFxsbK6XR6rdndu3evRo0apXbt2qlVq1ZKTEzU\nhg0bqp27sLBQP/rRjxQcHKyYmBgtXbpUOTk51eqZN2+enE6nPv/8c9133326+eab9YMf/ECSVFpa\nql/96lfq0aOHWrVqpXbt2mnUqFH67LPPvM518RhffPGF7rvvPrVr104dOnTQr3/9axljVFZWpnHj\nxqlt27YKDw/X/Pnzrb+Hvq43NTVVAwcOlCRNnTpVTqdTsbGxPo+bn5+vdu3aeYXmiw4cOKCkpCTP\n8+eff16pqalefRYvXuzV9otf/EL/+7//q08++aTW806aNEl/+tOfNGLEiCv2iY6O1oULF3TgwAGf\n1+FrPkmN93O8dC5NnjxZbdu2Vfv27fXLX/5Sp0+f9nktklRWVqa0tDRFRUWpZcuW+v73v6+lS5fW\n+ZoBNG3ccQZuIN98842OHTtWrb2iosLn2PT0dL311lt6+OGHddttt+nEiRPavn27XC6Xhg0bprFj\nx+qLL77Q2rVr9eyzz6pDhw6SpB49eujLL7/U7bffruDgYM2aNUutW7fWn/70J911111at26dxowZ\nI0n69NNPlZycrI4dO+qpp55S8+bN9dJLL6lVq1ZXXHM9YcIExcbGatGiRTp//ryk7+60FhQUaNy4\ncYqJidH//d//6YUXXtDQoUO1b98+derUyesY99xzj7p3767Fixdrw4YN+sMf/qA2bdrojTfe0KBB\ng7RkyRKtW7dO8+bNU58+fTR69Ohav1c21zt9+nTFxcVp3rx5evDBBzVkyBC1bt3a589hy5YtGjJk\nSI2vFRQUeK1lLiws1EsvveR5XlxcrMOHD+u2227ztLVq1Updu3bVRx99pNtvv93n+X1JTExUbm6u\nZsyYUWs/X/NJavyf48SJE9W5c2c988wz2r17t1auXKnDhw/7vIP+1VdfaeDAgaqsrNT06dMVGRmp\ngoICzZ49W0ePHtV//ud/Wl8zgCbOALju5eTkGIfDUevD6XRW619aWuppa9u2rXnkkUdqPc8zzzxT\nbZwxxowdO9a0aNHCfPHFF562U6dOma5du5quXbt62kaNGmWCg4PN4cOHPW3Hjh0z7dq1M06n0+u4\nc+fONQ6Hw4wbN65aHWfPnq3WdvDgQdOyZUuzaNGiasd44IEHvPreeuutxuFwmKeeesrTVlFRYdq3\nb29GjhxZ6/fA5nqrqqqMMcZs27bNOBwOs3r1ap/HNMaYCxcumLCwMJOVlVXj67/61a9MZWWlMcaY\nAwcOmH79+nm9npubaxwOhykuLvZq79u3b7XvwZVcnBtX8u6775p7773X53Fs5lNj/Rwvjr/zzju9\nxj/11FPG4XCYjRs3etpq+t345S9/aSIiIsxXX33lNf6JJ54wQUFBnr421wygaWOpBnADWb58uTZt\n2lTtcdddd/lcl9ymTRv97W9/09GjR+t0zsrKSm3YsEEjR47U97//fU97aGiopk+frsOHD2vPnj2q\nrKxUXl6eRo0apejoaE+/m2++WZMnT75ifenp6dXagoODPV+fPXtWx48fV1hYmL73ve9p586d1fo/\n8MADXs/79+8vSZo2bZqn7aabblLv3r118ODBel/v3r17az3GlezevVunT5/2Wo5xqQsXLnh2R/nH\nP/5R7Y7sP//5T0nf3WW+VEhIiOe1+rr55ptVXFzss5/NfGrsn+PDDz/s9fziXfOLa8ZrYozRunXr\ndNddd8nhcOjYsWOeR3JysqqqqjxLMa72dwhA00FwBm4g/fv3109+8pNqj6ioKJ9jlyxZoj179qhL\nly7q37+/fvvb317xTWKX+vrrr3X27Nkat1n7t3/7N0nfbZn21VdfqaKiQt/73veq9aup7aJbb721\nWltFRYWeeOIJRUVFqXXr1urYsaPCw8O1d+/eGnf8uHz3kItv1rs0wF9s/+abb65Yi2R/vVdjy5Yt\natu2bY3rm0tLSxUTE+NVR1hYmFefoKAgr/9edOHCBV24cOGqarpc+/btrXZVsZlPDflzrOkfCpfP\ns5tvvlnt2rWr9ef19ddf65tvvtHLL7+s8PBwr8cdd9whh8PheRPo1f4OAWg6CM4ArNx9990qKSnR\nihUr1LVrVy1btky9e/fWa6+9FtC6Lr0redGjjz6qrKwsTZw4UW+99ZY2btyovLw83XbbbaqqqqrW\n//IgeVFN66p93ZlvSB9//LEGDx5cY12bN2/2Wt9cVVVVbW/ujh07el671JkzZ6x39vDlSt/Ly9nM\np4b8Odry9fO+WMekSZNq/L85eXl5Gj9+vKSm+zsEwB5vDgRgLSIiQmlpaUpLS9PJkyc1cOBAZWZm\navLkyZJqDigdO3ZUq1at9Pnnn1d77WJbTEyMwsPDFRwcrC+//LJav5raarN27Vrdd999nr2MLzpx\n4oQnPDYU2+u9GsXFxbr77rtrfC0vL0/33Xef53mHDh2q3VW9uGtHeXm5IiIiPO0nTpxQXFzcVdV0\nuX/+85/W32Nf86mxf45ffvmlunXr5nl+7NgxffPNN7X+vDp27KiwsDD9v//3//STn/zE5zl8XTOA\npo07zgB8qqqqqva/xtu0aaOYmBivpQsX186eOHHC0xYUFKSf//znev/9973Wvp4+fVovvviiunbt\nql69eikoKEh33HGH/vKXv3htO3fs2DG9/vrrdbpr2KxZs2p3JP/85z+rrKzM+hhXy/Z6r0a7du1q\nvKu6bt06paSkeH2PIiMjvX4O0nfBuVu3bl7LA7766iuVlZVZhT4bJ06cUGRkZK19bOdTY/8c//jH\nP3o9X7ZsmSQpJSXlimOCgoI0btw4/fd//7d2795d7fWTJ0/qwoUL1tcMoGnjjjMAn06dOqVbbrlF\n48aNU+/evRUWFqZPPvlEH374odcbqi6+Ges3v/mNJk6cqBYtWuinP/2pFi1apLy8PA0ZMkQPPfSQ\nWrVqpZycHP3jH//QW2+95RmfmZmpDz/8UIMHD1Z6erqaNWumlStXKiYmRi6Xyzo8jxo1Sq+88orC\nwsJ02223yeVy6c0331RcXFy9l1rYjLe93rr65S9/qeXLl+vf//3f1bJlS0nffTDJF198oczMTK++\n3bt3V1lZWbUlG6mpqXrllVc8d65zcnI0aNAgry3uXn75ZS1YsEAul0tt27b1Ou7FIGuMqfHn4XK5\nfG5rZzufGvvnePToUd15550aMWKEPv30U61cuVLDhw/XHXfcUeuxFi9erC1btuj222/XtGnT9IMf\n/ECnTp3S3r179c477+jgwYNq0aKF1TUDaOJ8bbuxZcsWM3LkSHPLLbcYh8NhVq1a5XOrjvXr15sf\n/ehHJjQ01HTo0MGMHj3afPnll/XY/ANAfeTk5Bin02n+9re/1fj69OnTq21Hd+n2b+fPnzdPPPGE\n6du3r2nbtq1p1aqV6dWrl8nKyvJsf3bRokWLTJcuXUxQUJBxOp1my5Ytxhhj9u7da0aOHGnatGlj\ngoODTWJiovnggw+q1VJQUGAGDBhgbrrpJtOlSxfzzDPPmGXLlhmHw+G13de8efOM0+k05eXl1Y5x\n6tQpk5aWZjp16mRatWplkpKSzI4dO0xSUpL58Y9/7PMYl38/LhozZoyJjIys8Xt4OZvr3bZtm3E6\nndbb0RljTHZ2thk5cqSZOnWqmTJlinnllVeu2HfcuHFmx44dXm0XLlwwv/71r81DDz1k5s6dayZO\nnGjcbrdXn5deesmEh4d7bcO2fv16M3bsWBMeHm6cTqfp27evmTx5sjl58qTX2BEjRpg9e/bUeg22\n86mxfo4Xt6P77LPPzKRJk0ybNm1Mu3btzLRp08ypU6e8xl7+u3HR8ePHzWOPPWZiY2NNixYtTKdO\nncyQIUPMH/7wB3P+/Pk6/Q4BaLocxtT+z/YPPvhAn3zyifr27aspU6YoOztbU6ZMuWL/AwcOqGfP\nnpo1a5bS0tJ0+vRpzZ49WwcPHrTaoggALvfoo49q5cqV+te//lWvN3rdaN59910VFRVp0aJFdR57\n5swZrV+/XhMmTLAec+LECf385z/X9u3b63y+QJo3b57mz58vt9ut8PDwQJcDoAnzucY5JSVFCxcu\n1NixY6u9Q7smLpdLVVVVeuaZZxQXF6c+ffp4gvPl6+0A4HLffvut1/Ovv/5aa9as0dChQwnNdTR6\n9Ght27ZNZ8+erfNYl8vl9UY5G9nZ2Zo5c2adzwUA1wq/vznw9ttvV+vWrfXSSy+psrJSp0+f1qpV\nqzRgwAC1b9/e36cDcJ3p2rWrHn/8cf3Xf/2X5s2bp759++rs2bN66qmnAl3aNcfhcOh3v/ud5s+f\nX+exmzdvVt++fa37Hzp0SJ999pkmTZpU53MBwLXC728OjIyMVG5ursaMGaOHHnpIVVVV6tu3b62f\nvAQAF40cOVLvvPOO3G63mjVrpv79++vPf/6zBg0aFOjSrkkJCQn6+9//rg8//FA/+9nPrMbs2bNH\nPXv2rNMd/meffVbZ2dlXW2ZAORwO/m8GACs+1zhfKjQ0VM8//3yta5xLSko0cOBATZ06VZMmTdKp\nU6c8d4o++ugj/nACAADANcnvd5xffPFFRUdHa8mSJZ62V199VdHR0dq2bZsSExO9+nfr1k0HDx70\ndxkAAACAx6233qoDBw7U6xh+D87GmGpvIrz4vKaPSD148GBAP8IWTdO8efM0b968QJeBJoZ5gcsx\nJ1AT5gVq4o9VDz7fHHjmzBm5XC7PbhmlpaVyuVw6cuSIJCkjI0PDhw/39B81apR27dqlBQsWqLi4\nWLt27dLUqVPVpUsX9evXr94FAwAAAIHgMzgXFRUpPj5e8fHxqqio0Ny5cxUfH6+5c+dKktxut0pK\nSjz9Bw8erLVr1+q9995TfHy8UlJS1LJlS23YsEHBwcENdyUAAABAA6rTmwMbpACHg6UaqCY/P19J\nSUmBLgNNDPMCl2NOoCbMC9TEH5mT4AwAAIDrnj8yp98/AAUAAAC4HhGcAQAAAAsEZwAAAMACwRkA\nAACwQHAGAAAALBCcAQAAAAsEZwAAAMACwRkAAACwQHAGAAAALBCcAQAAAAsEZwAAAMACwRkAAACw\nQHAGAAAALBCcAQAAAAsEZwAAAMACwRkAAACwQHAGAAAALBCcAQAAAAsEZwAAAMACwRkAAACwQHAG\nAAAALBCcAQAAAAsEZwAAAMACwRkAAACwQHAGAAAALBCcAQAAAAvNAl0AcKk5c+bI7XYHugxrERER\nWrx4caDLAAAAjYDgjCbF7XYrJiYm0GVYO3ToUKBLAAAAjYSlGgAAAIAFgjMAAABggeAMAAAAWCA4\nAwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWPAZnAsK\nCjRq1Ch17txZTqdTq1evtjrws88+qx49eqhly5aKiopSRkZGvYsFAAAAAqWZrw5nzpxR7969dd99\n92nKlClyOBw+Dzpr1iytX79ev//979WrVy+dPHlSZWVlfikYAAAACASfwTklJUUpKSmSpNTUVJ8H\n/OKLL/THP/5Re/bsUffu3T3tffr0ufoqAQAAgADz+xrn9957T3FxccrNzVVcXJxiY2OVmpqqr7/+\n2t+nAgAAABqN34NzSUmJSktL9eabb+qVV17RmjVrtH//fo0cOVLGGH+fDgAAAGgUPpdq1FVVVZXO\nnTunNWvWqFu3bpKkNWvWqHv37tqxY4f69+/v71MCAAAADc7vwTkyMlLNmjXzhGZJ6tatm4KCgnT4\n8OEag/O8efM8XyclJSkpKcnfZQEAAOAGkp+fr/z8fL8e0+/BefDgwbpw4YJKSkoUFxcn6bvlG5WV\nleratWuNYy4NzgAAAEB9XX4zNjMzs97H9LnG+cyZM3K5XHK5XKqqqlJpaalcLpeOHDkiScrIyNDw\n4cM9/YcPH674+Hjdf//9crlc2r17t+6//34NHDhQCQkJ9S4YAAAACASfwbmoqEjx8fGKj49XRUWF\n5s6dq/j4eM2dO1eS5Ha7VVJS4unvcDj0/vvvKzw8XEOHDtXPf/5zdenSRe+9917DXQUAAADQwBwm\nwFtdOBwOdtuAR2pqqmJiYgJdhrVDhw5p1apVgS4DAAD44I/M6fft6AAAAIDrEcEZAAAAsEBwBgAA\nACwQnAEAAAALBGcAAADAAsEZAAAAsEBwBgAAACwQnAEAAAALBGcAAADAAsEZAAAAsEBwBgAAACwQ\nnAEAAAALBGcAAADAAsEZAAAAsEBwBgAAACwQnAEAAAALBGcAAADAAsEZAAAAsEBwBgAAACwQnAEA\nAAALBGcAAADAAsEZAAAAsEBwBgAAACwQnAEAAAALBGcAAADAAsEZAAAAsEBwBgAAACwQnAEAAAAL\nBGcAAADAAsEZAAAAsEBwBgAAACwQnAEAAAALBGcAAADAAsEZAAAAsEBwBgAAACwQnAEAAAALBGcA\nAADAAsEZAAAAsEBwBgAAACwQnAEAAAALBGcAAADAgs/gXFBQoFGjRqlz585yOp1avXq19cGLi4sV\nGhqq0NDQehUJAAAABJrP4HzmzBn17t1bzz33nIKDg+VwOKwOfP78eU2cOFHDhg2zHgMAAAA0Vc18\ndUhJSVFKSookKTU11frAs2fP1g9/+EMNHTpUW7ZsueoCAQAAgKagQdY4r1+/XuvXr9fy5ctljGmI\nUwAAAACNyucd57o6evSo0tLS9O677yokJMTfhwcAAAACwu/B+d5771V6err69+9vPWbevHmer5OS\nkpSUlOTvsgAAAHADyc/PV35+vl+P6TB1WEsRGhqq559/XlOmTLliH6fTqaCgIM9zY4yqqqoUFBSk\n7OxsPfDAA94FOBws54BHamqqYmJiAl2GtUOHDmnVqlWBLgMAAPjgj8zp9zvOe/fu9Xr+7rvvatGi\nRSoqKlJUVJS/TwcAAAA0Cp/B+cyZMyouLpYkVVVVqbS0VC6XSzfffLOio6OVkZGhoqIibdq0SZLU\ns2dPr/Hbt2+X0+ms1g4AAABcS3zuqlFUVKT4+HjFx8eroqJCc+fOVXx8vObOnStJcrvdKikpqfUY\n7OMMAACAa12d1jg3SAGsccYlWOMMAAAagj8yZ4Ps4wwAAABcbwjOAAAAgAW/76pxvZszZ47cbneg\ny7AWERGhxYsXB7oMAABwjSDrXBnBuY7cbvc1twYXAADAFlnnyliqAQAAAFggOAMAAAAWCM4AAACA\nBYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIz\nAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAA\nYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAAYIHg\nDAAAAFggOAMAAAAWCM4AAACABYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAAYMFncC4o\nKNCoUaPUuXNnOZ1OrV69utb++fn5Gj16tKKiotSqVSv16dNHOTk5fisYAAAACIRmvjqcOXNGvXv3\n1n333acpU6bI4XDU2n/btm3q06eP5syZo8jISG3YsEFpaWlq2bKl7rnnHr8VDjs7d+5UampqoMuw\ntnPnTsXExAS6DAColzlz5sjtdge6DGsRERFavHhxoMsAmjyfwTklJUUpKSmSZBXAMjIyvJ5Pnz5d\nmzdv1ttvv01wDoBvv/32mgqihYWFgS4BAOrN7XZfU3/2Hjp0KNAlANeERlnjfPLkSbVv374xTgUA\nAAA0CJ93nOvr/fff10cffaStW7c29KkAAACABtOgd5w/+eQTTZ48WcuXL1dCQkJDngoAAABoUA12\nx7mwsFAjRozQggUL9OCDD9bad968eZ6vk5KSlJSU1FBlAQAA4AaQn5+v/Px8vx6zQYJzQUGB7rrr\nLs2fP18zZszw2f/S4AwAAADU1+U3YzMzM+t9TKvt6IqLiyVJVVVVKi0tlcvl0s0336zo6GhlZGSo\nqKhImzZtkvRduh8xYoQefvhh3XPPPZ7teIKCgtSxY8d6FwwAAAAEgs81zkVFRYqPj1d8fLwqKio0\nd+5cxcfHa+7cuZK+23KnpKTE03/16tWqqKjQ0qVLFRkZqaioKEVFRelHP/pRw10FAAAA0MB83nFO\nSkpSVVXVFV+//FMBc3Jy+KRAAAAAXHcaZR9nAAAA4FpHcAYAAAAsEJwBAAAACwRnAAAAwALBGQAA\nALBAcAYAAAAsEJwBAAAACwRnAAAAwALBGQAAALBAcAYAAAAsEJwBAAAACwRnAAAAwALBGQAAALBA\ncAYAAAAsEJwBAAAACwRnAAAAwALBGQAAALBAcAYAAAAsEJwBAAAACwRnAAAAwALBGQAAALBAcAYA\nAAAsEJwBAAAACwRnAAAAwALBGQAAALBAcAYAAAAsEJwBAAAACwRnAAAAwEKzQBcgSa+//nqgS7DS\nokULXbhwIdBlAAAAIACaRHDesmVLoEuwcvbsWZ07dy7QZQDADWfOnDlyu92BLsPazp07FRMTE+gy\n0EQwf68fTSI4R0ZGBroEK0eOHFFFRUWgywCAG47b7b6m/iIvLCwMdAloQpi/1w/WOAMAAAAWCM4A\nAACABYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAAYIHgDAAAAFggOAMAAAAWfAbngoIC\njRo1Sp07d5bT6dTq1at9HnTPnj0aNmyYQkJC1LlzZy1YsMAvxQIAAACB4jM4nzlzRr1799Zzzz2n\n4OBgORyOWvufOnVKd9xxhyIjI7Vjxw4999xzWrp0qbKysvxWNAAAANDYmvnqkJKSopSUFElSamqq\nzwO+9tprqqio0OrVq3XTTTepZ8+e2r9/v7KysjRr1qx6FwwAAAAEgt/XOG/btk1DhgzRTTfd5GlL\nTk7W0aNHVVpa6u/TAQAAAI3C78HZ7XarU6dOXm0Xn7vdbn+fDgAAAGgUfg/OvtZAAwAAANcin2uc\n6yoiIqLaneXy8nLPazXJz8/3fB0TE6OYmBh/lwUAAIAbSH5+vlfG9Ae/B+dBgwZp9uzZOnfunGed\nc15enm655RZ17dq1xjFJSUn+LgMAAAA3sKSkJK+MmZmZWe9jWm1H53K55HK5VFVVpdLSUrlcLh05\nckSSlJGRoeHDh3v6T5o0SSEhIUpNTdW+ffv0zjvvaMmSJeyoAQAAgGuaz+BcVFSk+Ph4xcfHq6Ki\nQnPnzlV8fLzmzp0r6bs3/JWUlHj6h4WFKS8vT0ePHlVCQoIeeeQRPf7443rsscca7ioAAACABuZz\nqUZSUpKqqqqu+HpOTk61th/84AfasmVL/SoDAAAAmhC/76oBAAAAXI8IzgAAAIAFgjMAAABggeAM\nAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAA\nWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABaaBboAAI1nzpw5crvdgS7DWkREhBYv\nXhzoMgAAkERwBm4obrdbMTExgS7D2qFDhwJdAgAAHizVAAAAACwQnAEAAAALBGcAAADAAsEZAAAA\nsEBwBgAAACwQnAEAAAALBGcAAADAAsEZAAAAsEBwBgAAACwQnAEAAAALBGcAAADAAsEZAAAAsEBw\nBgAAACwQnAEAAAALBGcAAADAAsEZAAAAsEBwBgAAACwQnAEAAAALBGcAAADAAsEZAAAAsEBwBgAA\nACxYBecVK1YoNjZWwcHBSkhIUGFhYa39c3NzNXDgQIWFhaljx44aM2aMiouL/VIwAAAAEAg+g/Pa\ntWs1c+ZMPfnkk3K5XEpMTFRKSoqOHDlSY/8DBw5ozJgxSkpKksvl0qZNm1RRUaE777zT78UDAAAA\njcVncM7KytLUqVM1bdo0de/eXcuWLVNkZKSys7Nr7O9yuVRVVaVnnnlGcXFx6tOnj2bPnq2DBw/q\nxIkTfr/xgGs7AAAWiklEQVQAAAAAoDHUGpzPnz+vXbt2KTk52as9OTlZW7durXHM7bffrtatW+ul\nl15SZWWlTp8+rVWrVmnAgAFq3769/yoHAAAAGlGtwfnYsWOqrKxUp06dvNrDw8PldrtrHBMZGanc\n3Fw9+eSTatmypdq2bat9+/bpr3/9q/+qBgAAABqZ33fVKCkp0ZgxYzR16lTt2LFD+fn5Cg0N1YQJ\nE2SM8ffpAAAAgEbRrLYXO3TooKCgIJWXl3u1l5eXKzIyssYxL774oqKjo7VkyRJP26uvvqro6Ght\n27ZNiYmJ1cbk5+d7vo6JiVFMTEwdLgEAAADwlp+f75Ux/aHW4NyiRQv169dPGzdu1NixYz3teXl5\nGj9+fI1jjDFyOr1vZF98XlVVVeOYpKSkutQMAAAA1CopKckrY2ZmZtb7mD6XasyaNUurVq3Syy+/\nrM8//1yPPvqo3G63pk+fLknKyMjQ8OHDPf1HjRqlXbt2acGCBSouLtauXbs0depUdenSRf369at3\nwQAAAEAg1HrHWZImTJig48ePa+HChSorK1OvXr2Um5ur6OhoSZLb7VZJSYmn/+DBg7V27VotXrxY\nv/vd7xQSEqJBgwZpw4YNCg4ObrgrAQAAABqQz+AsSenp6UpPT6/xtZycnGpt48aN07hx4+pXGQAA\nANCE+H1XDQAAAOB6RHAGAAAALBCcAQAAAAsEZwAAAMACwRkAAACwYLWrBoCa7dy5U6mpqYEuw9rO\nnTuvqU/mvNa+vxEREVq8eHGgy7AyZ84cud3uQJdh7Vqbu9caftcAOwRnoB6+/fbba+ov88LCwkCX\nUCfX2vf30KFDgS7Bmtvtvqa+t9fa3L3W8LsG2GGpBgAAAGCB4AwAAABYIDgDAAAAFgjOAAAAgAWC\nMwAAAGCB4AwAAABYIDgDAAAAFgjOAAAAgAWCMwAAAGCB4AwAAABYIDgDAAAAFgjOAAAAgAWCMwAA\nAGCB4AwAAABYIDgDAAAAFgjOAAAAgAWCMwAAAGCB4AwAAABYIDgDAAAAFgjOAAAAgAWCMwAAAGCB\n4AwAAABYIDgDAAAAFgjOAAAAgAWCMwAAAGCB4AwAAABYIDgDAAAAFgjOAAAAgAWCMwAAAGCB4AwA\nAABYIDgDAAAAFgjOAAAAgAWCMwAAAGDBKjivWLFCsbGxCg4OVkJCggoLC32OefbZZ9WjRw+1bNlS\nUVFRysjIqHexAAAAQKA089Vh7dq1mjlzprKzszV48GA9//zzSklJ0Weffabo6Ogax8yaNUvr16/X\n73//e/Xq1UsnT55UWVmZ34sHAAAAGovP4JyVlaWpU6dq2rRpkqRly5Zpw4YNys7O1tNPP12t/xdf\nfKE//vGP2rNnj7p37+5p79Onjx/LBgAAABpXrUs1zp8/r127dik5OdmrPTk5WVu3bq1xzHvvvae4\nuDjl5uYqLi5OsbGxSk1N1ddff+2/qgEAAIBGVusd52PHjqmyslKdOnXyag8PD5fb7a5xTElJiUpL\nS/Xmm2/qlVdekSQ9/vjjGjlypLZt2yaHw+Gn0gGgadm5c6dSU1MDXYaVnTt3KiYmJtBlAFflWvpd\nk/h9u574XKpRV1VVVTp37pzWrFmjbt26SZLWrFmj7t27a8eOHerfv7+/TwkATcK33357zfzlaPMm\nb6CpupZ+1yR+364ntQbnDh06KCgoSOXl5V7t5eXlioyMrHFMZGSkmjVr5gnNktStWzcFBQXp8OHD\nNQbn/Px8z9cxMTHX1C8DAAAAmp78/HyvjOkPtQbnFi1aqF+/ftq4caPGjh3rac/Ly9P48eNrHDN4\n8GBduHBBJSUliouLk/Td8o3Kykp17dq1xjFJSUlXWT4AAABQXVJSklfGzMzMrPcxfe7jPGvWLK1a\ntUovv/yyPv/8cz366KNyu92aPn26JCkjI0PDhw/39B8+fLji4+N1//33y+Vyaffu3br//vs1cOBA\nJSQk1LtgAAAAIBB8rnGeMGGCjh8/roULF6qsrEy9evVSbm6uZw9nt9utkpIST3+Hw6H3339fM2bM\n0NChQxUcHKzk5GRlZWU13FUAAAAADczqzYHp6elKT0+v8bWcnJxqbREREXrzzTfrVxkAAADQhFh9\n5DYAAABwoyM4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAM\nAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAA\nWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4\nAwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAA\nABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWLAKzitWrFBs\nbKyCg4OVkJCgwsJCq4MXFxcrNDRUoaGh9SoSAAAACDSfwXnt2rWaOXOmnnzySblcLiUmJiolJUVH\njhypddz58+c1ceJEDRs2TA6Hw28FAwAAAIHgMzhnZWVp6tSpmjZtmrp3765ly5YpMjJS2dnZtY6b\nPXu2fvjDH2r8+PEyxvitYAAAACAQag3O58+f165du5ScnOzVnpycrK1bt15x3Pr167V+/XotX76c\n0AwAAIDrQrPaXjx27JgqKyvVqVMnr/bw8HC53e4axxw9elRpaWl69913FRIS4r9KAQAAgACqNThf\njXvvvVfp6enq37+/9Zj8/HzP1zExMYqJifF3WQAAALiB5Ofne2VMf6g1OHfo0EFBQUEqLy/3ai8v\nL1dkZGSNYzZv3qyCggJlZmZKkowxqqqqUvPmzZWdna0HHnig2pikpKSrLB8AAACoLikpyStjXsym\n9VFrcG7RooX69eunjRs3auzYsZ72vLw8jR8/vsYxe/fu9Xr+7rvvatGiRSoqKlJUVFS9CwYAAAAC\nwedSjVmzZunee+/VgAEDlJiYqBdeeEFut1vTp0+XJGVkZKioqEibNm2SJPXs2dNr/Pbt2+V0Oqu1\nAwAAANcSn8F5woQJOn78uBYuXKiysjL16tVLubm5io6OliS53W6VlJTUegz2cQYAAMC1zurNgenp\n6UpPT6/xtZycnFrHpqamKjU1tc6FAQAAAE2J1UduAwAAADc6gjMAAABggeAMAAAAWCA4AwAAABYI\nzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAA\nAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAF\ngjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMA\nAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABggeAMAAAAWCA4AwAAABYIzgAAAIAFgjMAAABg\ngeAMAAAAWCA4AwAAABasg/OKFSsUGxur4OBgJSQkqLCw8Ip98/PzNXr0aEVFRalVq1bq06ePcnJy\n/FIwAAAAEAhWwXnt2rWaOXOmnnzySblcLiUmJiolJUVHjhypsf+2bdvUp08fvf3229q3b5/S09OV\nlpamP//5z34tHgAAAGgszWw6ZWVlaerUqZo2bZokadmyZdqwYYOys7P19NNPV+ufkZHh9Xz69Ona\nvHmz3n77bd1zzz1+KBsAAABoXD7vOJ8/f167du1ScnKyV3tycrK2bt1qfaKTJ0+qffv2da8QAAAA\naAJ83nE+duyYKisr1alTJ6/28PBwud1uq5O8//77+uijj+oUtAEAAICmpMF31fjkk080efJkLV++\nXAkJCQ19OgAAAKBB+Lzj3KFDBwUFBam8vNyrvby8XJGRkbWOLSws1IgRI7RgwQI9+OCDV+yXn5/v\n+TomJkYxMTG+ygIAAACuKD8/3ytj+oPP4NyiRQv169dPGzdu1NixYz3teXl5Gj9+/BXHFRQU6K67\n7tL8+fM1Y8aMWs+RlJRkXzEAAADgQ1JSklfGzMzMrPcxrXbVmDVrlu69914NGDBAiYmJeuGFF+R2\nuzV9+nRJ3+2iUVRUpE2bNkn6LuGPGDFCDz/8sO655x7PWuigoCB17Nix3kUDAAAAjc0qOE+YMEHH\njx/XwoULVVZWpl69eik3N1fR0dGSJLfbrZKSEk//1atXq6KiQkuXLtXSpUs97TExMV79AAAAgGuF\nVXCWpPT0dKWnp9f42uWfCpiTk8MnBQIAAOC60uC7agAAAADXA4IzAAAAYIHgDAAAAFggOAMAAAAW\nCM4AAACABYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAAYIHgDAAAAFggOAMAAAAWCM4A\nAACABYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACA\nBYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIz\nAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAAYIHgDAAAAFggOAMAAAAWCM4AAACABYIzAAAA\nYIHgDAAAAFggOAMAAAAWrILzihUrFBsbq+DgYCUkJKiwsLDW/nv27NGwYcMUEhKizp07a8GCBX4p\nFgAAAAgUn8F57dq1mjlzpp588km5XC4lJiYqJSVFR44cqbH/qVOndMcddygyMlI7duzQc889p6VL\nlyorK8vvxQMAAACNxWdwzsrK0tSpUzVt2jR1795dy5YtU2RkpLKzs2vs/9prr6miokKrV69Wz549\nNXbsWM2ePZvgjDo5dOhQoEtAE8S8wOWYE6gJ8wINpdbgfP78ee3atUvJycle7cnJydq6dWuNY7Zt\n26YhQ4bopptu8up/9OhRlZaW+qFk3Aj4Qw81YV7gcswJ1IR5gYZSa3A+duyYKisr1alTJ6/28PBw\nud3uGse43e5q/S8+v9IYAAAAoKlr5u8DOhyOOo+50nrppqaqquqqrg8AAADXAVOLc+fOmWbNmpl1\n69Z5tf/qV78ySUlJNY6ZMmWKGTFihFfb9u3bjcPhMIcOHarW/9ZbbzWSePDgwYMHDx48ePBosMet\nt95aW+y1Uusd5xYtWqhfv37auHGjxo4d62nPy8vT+PHjaxwzaNAgzZ49W+fOnfOsc87Ly9Mtt9yi\nrl27Vut/4MCB2koAAAAAmgSfu2rMmjVLq1at0ssvv6zPP/9cjz76qNxut6ZPny5JysjI0PDhwz39\nJ02apJCQEKWmpmrfvn165513tGTJEs2aNavhrgIAAABoYD7XOE+YMEHHjx/XwoULVVZWpl69eik3\nN1fR0dGSvnvDX0lJiad/WFiY8vLy9NBDDykhIUHt27fX448/rscee6zhrgIAAABoYA5jjAl0EQAA\nAEBTZ/WR2/XBx3WjJnWZF/n5+Ro9erSioqLUqlUr9enTRzk5OY1YLRpDXf+suKi4uFihoaEKDQ1t\n4AoRCFczL5599ln16NFDLVu2VFRUlDIyMhqhUjSmus6L3NxcDRw4UGFhYerYsaPGjBmj4uLiRqoW\nDa2goECjRo1S586d5XQ6tXr1ap9jrjpv1vvthbV44403TPPmzc3KlSvN/v37zSOPPGJat25tDh8+\nXGP/kydPmk6dOpm7777b7Nu3z6xbt86EhoaaP/zhDw1ZJhpZXefF008/bX7729+arVu3mr///e8m\nOzvbNGvWzLz++uuNXDkaSl3nxEXnzp0z8fHxZsSIESY0NLSRqkVjuZp58dhjj5nvf//75i9/+Yv5\n+9//blwul/nggw8asWo0tLrOi+LiYtO8eXMze/Zsc/DgQeNyuczPfvYz061bt0auHA0lNzfX/Md/\n/IdZt26dCQkJMatXr661f33yZoMG5wEDBpi0tDSvtu9973smIyOjxv4rVqwwbdq0MRUVFZ62hQsX\nmltuuaUhy0Qjq+u8qMmECRPM2LFj/V0aAuRq58TMmTPN/fffb1atWmVat27dkCUiAOo6L/bv32+a\nN29u9u/f3xjlIUDqOi/eeustExQUZKqqqjxtH330kXE4HOb48eMNWisaX+vWrX0G5/rkzQZbqsHH\ndaMmVzMvanLy5Em1b9/e3+UhAK52Tqxfv17r16/X8uXLZXirxnXnaubFe++9p7i4OOXm5iouLk6x\nsbFKTU3V119/3RgloxFczby4/fbb1bp1a7300kuqrKzU6dOntWrVKg0YMIC/R25Q9cmbDRac+bhu\n1ORq5sXl3n//fX300UdKS0triBLRyK5mThw9elRpaWl67bXXFBIS0hhlopFdzbwoKSlRaWmp3nzz\nTb3yyitas2aN9u/fr5EjR/KPq+vE1cyLyMhI5ebm6sknn1TLli3Vtm1b7du3T3/9618bo2Q0QfXJ\nmw3+5sC64OOs4csnn3yiyZMna/ny5UpISAh0OQiQe++9V+np6erfv3+gS0ETUlVVpXPnzmnNmjUa\nPHiwBg8erDVr1mj79u3asWNHoMtDgJSUlGjMmDGaOnWqduzYofz8fIWGhmrChAn8g+oGVZ+82WDB\nuUOHDgoKClJ5eblXe3l5uSIjI2scExERUS3pXxwfERHRMIWiUV3NvLiosLBQd955pxYsWKAHH3yw\nIctEI7qaObF582ZlZmaqefPmat68uR544AGdOXNGzZs318qVKxujbDSwq5kXkZGRatasmbp16+Zp\n69atm4KCgnT48OEGrReN42rmxYsvvqjo6GgtWbJEffr00ZAhQ/Tqq69qy5Yt2rZtW2OUjSamPnmz\nwYLzpR/Xfam8vDwlJibWOGbQoEH6+OOPde7cOa/+V/q4blx7rmZeSN9tNXPnnXcqMzNTM2bMaOgy\n0YiuZk7s3btXn376qecxf/58BQcH69NPP9W4ceMao2w0sKuZF4MHD9aFCxe8PpSrpKRElZWV/B1y\nnbiaeWGMkdPpHXcuPq+qqmqYQtGk1Stv1uutiz6sXbvWtGjRwqxcudJ89tlnZsaMGSY0NNSzZcyc\nOXPMT3/6U0//kydPmoiICDNx4kSzd+9e8/bbb5uwsDCTlZXVkGWikdV1XmzevNmEhISYJ554wrjd\nblNWVmbKysrMV199FahLgJ/VdU5cLicnh101rkN1nRdVVVWmX79+ZtiwYWb37t1m165dZujQoWbQ\noEGBugQ0gLrOi48//tg4nU4zf/588+WXX5qdO3ean/3sZ6Zr167m7NmzgboM+NG//vUvs3v3brN7\n924TEhJi5s+fb3bv3t0gebNBg7Mx3235ERMTY2666SaTkJBgPv74Y89rqampJjY21qv/nj17zNCh\nQ03Lli1NVFSUmT9/fkOXiACoy7xITU01TqfTOBwOr8flcwfXtrr+WXGpnJwc9nG+TtV1XpSVlZnx\n48eb0NBQEx4ebn7xi1/wj+zrUF3nxVtvvWX69etnWrdubcLDw83o0aPN559/3thlo4Fs3rzZkw0u\nzQtTp041xvg3b/KR2wAAAICFJrWrBgAAANBUEZwBAAAACwRnAAAAwALBGQAAALBAcAYAAAAsEJwB\nAAAACwRnAAAAwALBGQAAALBAcAYAAAAs/H8ckUzP0yiInwAAAABJRU5ErkJggg==\n", "text": [ "" ] } ], "prompt_number": 2 }, { "cell_type": "code", "collapsed": false, "input": [ "def normal(x, mu, sigma):\n", " \"\"\"Normal distribution PDF.\"\"\"\n", " return np.exp(-0.5 * ((x - mu) / sigma)**2) / (np.sqrt(2 * np.pi) * sigma)" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 3 }, { "cell_type": "code", "collapsed": false, "input": [ "def _fit_gmm(data, num_components, num_iters=100):\n", " \"\"\"Fit a single GMM with EM with one random initialization.\"\"\"\n", "\n", " # Random initialization.\n", " mu = np.random.choice(data, num_components, replace=False)\n", " sigma = (np.random.random(num_components) * 0.15) + 0.1\n", " prob_population = np.ones(num_components) / num_components\n", "\n", " # Keep track of results after each iteration.\n", " results = []\n", " results.append((mu.copy(), sigma.copy(), prob_population.copy(), float('-inf')))\n", " \n", " last_log_likelihood = None\n", " for i in xrange(num_iters):\n", " # E-step.\n", " probs = [normal(x, mu, sigma) for x in data]\n", " probs = np.array([p / p.sum() for p in probs])\n", "\n", " # M-step.\n", " for k in xrange(num_components):\n", " k_probs = probs[:, k]\n", " mu[k] = (k_probs * data).sum() / k_probs.sum()\n", " sigma[k] = np.sqrt((k_probs * (data - mu[k])**2).sum() / k_probs.sum())\n", " prob_population[k] = k_probs.sum() / len(data)\n", "\n", " # Bookkeeping.\n", " log_likelihood = np.log(np.product([(normal(data[n], mu, sigma) * probs[n, :]).sum()\n", " for n in xrange(len(data))]))\n", " results.append((mu.copy(), sigma.copy(), prob_population.copy(), log_likelihood))\n", " if last_log_likelihood is not None and log_likelihood <= (last_log_likelihood + 0.01):\n", " break\n", " last_log_likelihood = log_likelihood\n", "\n", " return results\n", "\n", "\n", "def fit_gmm(data, num_components, num_iters=10, num_random_inits=10):\n", " \"\"\"Find the maximum likelihood GMM over several random initializations.\"\"\"\n", " best_results = None\n", " best_results_iters = None\n", " best_ll = float('-inf')\n", "\n", " # Try several random initializations and keep the best.\n", " for attempt in xrange(num_random_inits):\n", " results_iters = _fit_gmm(data, num_components, num_iters=num_iters)\n", " final_log_likelihood = results_iters[-1][3]\n", " if final_log_likelihood > best_ll:\n", " best_results = results_iters[-1]\n", " best_results_iters = results_iters\n", " best_ll = final_log_likelihood\n", "\n", " return best_results, best_results_iters" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 4 }, { "cell_type": "code", "collapsed": false, "input": [ "colors = 'bgrcmy'\n", "\n", "def gmm_fit_and_animate(data, num_components, interval=200):\n", " _, best_results_iters = fit_gmm(data, num_components, num_iters=200, num_random_inits=10)\n", "\n", " # Remove initial random guess (before doing a single iteration).\n", " best_results_iters = best_results_iters[1:]\n", " \n", " fig = plt.figure(figsize=(12, 8))\n", " ax = plt.axes(xlim=(0, 1), ylim=(0, 2))\n", " \n", " line, = ax.plot([], [], label='GMM Fit', color='black', alpha=0.7, linewidth=3)\n", " \n", " ax.hist(data, normed=True, bins=15, color='lightgray', alpha=0.2, label='Real Data')\n", " ax.legend()\n", " ax.set_title('{0} Components'.format(num_components))\n", "\n", " def animate(i):\n", " mu, sigma, prob_population, _ = best_results_iters[i]\n", " xs = np.linspace(0, 1, 1000)\n", " ys = [(normal(x, mu, sigma) * prob_population).sum() for x in xs]\n", " line.set_data(xs, ys)\n", " \n", " for k in xrange(num_components):\n", " ys = [normal(x, mu[k], sigma[k]) * prob_population[k] for x in xs]\n", " ax.plot(xs, ys, alpha=0.2, color=colors[k % len(colors)])\n", "\n", " # Things like to crash if I try to do too many frames, I guess, so limit\n", " # the number of frames.\n", " num_iters = len(best_results_iters)\n", " frames = np.arange(0, num_iters, max(1, num_iters // 20), dtype=int)\n", "\n", " return animation.FuncAnimation(fig, animate, frames=frames, interval=interval)" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 5 }, { "cell_type": "code", "collapsed": false, "input": [ "gmm_fit_and_animate(data, 3)" ], "language": "python", "metadata": {}, "outputs": [ { "html": [ "\n", "\n", "\n", "
\n", " \n", "
\n", " \n", "
\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
\n", " Once \n", " Loop \n", " Reflect \n", "
\n", "
\n", "\n", "\n", "\n" ], "metadata": {}, "output_type": "pyout", "prompt_number": 6, "text": [ "" ] } ], "prompt_number": 6 }, { "cell_type": "code", "collapsed": false, "input": [ "gmm_fit_and_animate(data, 6)" ], "language": "python", "metadata": {}, "outputs": [ { "html": [ "\n", "\n", "\n", "
\n", " \n", "
\n", " \n", "
\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
\n", " Once \n", " Loop \n", " Reflect \n", "
\n", "
\n", "\n", "\n", "\n" ], "metadata": {}, "output_type": "pyout", "prompt_number": 7, "text": [ "" ] } ], "prompt_number": 7 } ], "metadata": {} } ] }