gps_notebook.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "\n",
    "import h5py\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "\n",
    "from sklearn import mixture\n",
    "from mvn import MVN"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# constants\n",
    "RESAMPLE = 6\n",
    "GMM_DYN_COMPONENTS = 10\n",
    "INITIAL_VARIANCE_RATIO = 0.1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# load data\n",
    "f = h5py.File( 'data.h5', 'r' )\n",
    "x = f[ 'x' ]\n",
    "x = x[ ::RESAMPLE ]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.text.Text at 0xaf6af7f0>"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAa0AAAEACAYAAAAA1CkHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl8VPX1//HXZF9YEgJJWEIIO0EQEJBFJBBRUFxwqUtV\n1KpVSr/aWqvirxWXr1q1alsrikgF9au1olYtFESJoiAEZQ0g+xKBBNmXLJPk/v64MzBMZjLbvXM/\nM3Oej8c8SGYmMyeXyT33876fey8IIYQQQgghhBBCCCGEEEIIIYQQQgghhBBCCCGEEKIJO4A1wEpg\nubWlCCFE1JoJVABrXe5rBXwGbAIWABkujz0EbAY2AheGqcaIsB19wQkhhDDPCKA/ZzatZ4DfO75+\nAHja8XUhsApIBDoBW4C4sFQZAbYDWVYXIYQQMaATZzatjUCO4+tcx/egj7IecHnef4EhZhfniypd\nUwMWAiuAOyyuRQghYkkOemSI419nA2sHlLs8rxxoH8a6PEqwugCH4cBeoA16troRWGxpRUIIEXs0\nx62pxy2lykhrr+Pf/cCHwGDXB7t06eJckIrf2mpQ4fi3rQYLNLhUgbrkJje5WXc7rkGuVe+/Bd8q\n0GNBgLZApePrH4E8l+d1cNxnKRWaVhrQ3PF1OvoMFde8la1bt6JpmvK3rVv30KlTNpq2B03bw1VX\njeG99z429D0eeeQRy39PVW6yLGRZqL4sGho04uLS2b59ryXvD3TxYx38MTDR8fVE4COX+68DkoAC\noBsKzO5WIR7MQR9dgV7P2+jTLiNOdTWkpJz+PiVFv08IEZvq6qChQan1wDvASKA1sBv4I/pswfeA\nX6AffvQzx3PXO+5fD9QBk9BHb5ZSoWltB/pZXYQRqqshOfn098nJSn1YhRBh5vz7V2g9cL2X+y/w\ncv+TjpsyVIgHo0ZNTeORVk2Nse9RVFRk7AtGMFkWp8myOE2lZeH8+zd6PRDLVBhpRQ2j48HD1Yd5\nufRlEuISaJPWhuz0bEYOGxl6oVFCpZWT1WRZnKbSslBwpBXxZKRlICOb1sJtC+k7rS8/HPiByhOV\nlOws4ZklzzDotUH88NMPxhQshDCVNC3jyUjLQDU1jfdpBRoLVNmreGDhA3y48UNev+x1Luxy5um+\nXvvuNUb8YwQzL5/J+O7jDahaCGEWiQeNJ03LQJ5GWseOBfYa9y24j91Hd7PmrjVkpmY2evyOc+7g\nrOyzuOZf13Df0Pv4zdDfhFi1EMIsMtIynjQtA4UaDy7avohPNn3C2rvXkpGS4fV5Q/OGsuz2ZZwz\n/RyGdBjC0LyhIVQthDCLNC3jyT4tA4USD56oPcHtn9zOK5e80mTDcmrfoj0vXfwSEz+ayEn7ySAr\nFkKYSeJB40nTMlAoI62Hv3iY4XnDuaT7JX6/39WFVzOo/SAeWvhQgJUKIcJBRlrGk6ZloGCb1je7\nvuG9svd4ceyLAb/n38b9jTkb5rBo+6KAf1YIYS5pWsaTpmWgYOPB3y74LS9c9AKtUgO/Dmar1FZM\nv3Q6t/77VokJhVCMxIPGk6ZloGBGWqU/llJxvIKrC68O+n0v7nYx/XL78fr3rwf9GkII48lIy3jS\ntAwUTNOatmIadw28i/i4+JDee8qIKTy75Flq62tDeh0hhHGkaRlPmpaBAo0HD1Yd5IMNH3Bb/9tC\nfu/B7QfTo3UP3lrzVsivJYQwRk0NxMVJPGgkaVoGCnSk9caqNxjffTzZ6dmGvP/DIx7m6a+fpr6h\n3pDXE0KEproaMjJkpGUkaVoGCqRpNWgNTFsxjUmDJhn2/iPzR9I6rTVzNswx7DWFEMGrroaWLaVp\nGUmaloECiQcXbltIs6RmDO1g3NksbDYbU0ZM4cnFTzqvWiqEsFBNjTQto0nTMlAgI62XS19m0sBJ\n2Gw2Q2u4pNslaGjM2zLP0NcVQgTOGQ/KPi3jSNMykL9Na/+J/ZTsKOGGPjcYXoPNZmPyoMn8Y9U/\nDH9tIURgJB40njQtA/kbD87fOp+iTkWkJ6WbUsfVhVezYOsCjlQfMeX1hRD+kXjQeNK0DOQ+0kpO\n1u9z3700b8s8Lu52sWl1ZKZmMrpgtEzIEMJiEg8aT5qWgdybVkICxMeD3X76vvqGehZsXcC4ruNM\nreXnfX7O22vfNvU9hBBNkynvxpOmZaCamjObFujfu25lrdizgpz0HPJa5play/ju41m5dyU/Hv3R\n1PcRQngn8aDxpGkZqLr6zH1acDoidJq7ea7poyyAlIQUJvScwLvr3jX9vYQQnjknYkg8aBxpWgZy\njweh8QxCs/dnufp5X4kIhbCSxIPGk6ZlIF/x4P4T+9l0YBPDOw4PSz0j80dScaKC9fvXh+X9hBBn\nknjQeNK0DOQrHpy/dT6jCkaRFJ8Ulnri4+K5/qzreXuNjLaEsILEg8aTpmWQhgaorW3ctFzjwXDt\nz3J13VnX8f6G98P6nkIInWvTkjOrGUOalkGcBxa7n5XJ2bTCNdXd3YC2Azhac5StB7eG9X2FEPrf\nfmqqfvhLrVzqzhAqNa14YCXwidWFBMP9bBhOzrNilO4ppW3ztqZPdXcXZ4tjbNexci5CISzg3M/t\nfuiLCJ5KTeseYD0QkYNoTzMH4fRIq2RHCcUFxeEvDLi468XStISwgHO94M9VzIV/VGlaHYCLgRmA\nsac9DxNfTWvFnhUMbj84/IUBY7qMYfHOxVTZqyx5fyFilXNylvvxmiJ4qjStF4D7gQarCwmWP/Hg\noHaDwl8YkJGSQb/cfpTsKLHk/YWIRfX1+i0xUUZaRlKhaY0HKtH3Z0XkKAuaHmlVHjnKkeojdGnV\nJfyFOVzcTSJCIcLJuT/LZpN9WkZKsLoAYBhwGXo8mAK0AGYDN7s+aerUqae+LioqoqioKGwF+qOp\nprVp3y4GdhlInM26bYSLu13MhH9O4C9j/2L4hSeFEI25Hrcp8aBxVGhaUxw3gJHA73BrWHBm01JR\nU/Hghp/KGTTCmmjQqU92H2rqath8cDPds7pbWosQscD1DDkSDxpHhXjQXdTNHtzx0z4Gtbe2adls\nNsZ1Hce8zRIRChEOrusEiQeNo1rT+hI9Kow43ppWcrJG+cGfLJuE4eribhczd8tcq8sQIiZIPGgO\n1ZpWxPJ0slyAag5DfTIdWnQIf1FuijsXs2T3Ek7aT1pdihBRT+JBc0jTMoink+UCVFTvIie5oxKT\nH1okt6Bn656s3rfa6lKEiHoSD5pDmpZBvMWDe6u3kZXYPvwFedEvpx+r9q2yugwhop7Eg+aQpmUQ\nb/HgrhObyUjMDX9BXvTLlaYlRDhIPGgOaVoG8RQPNmgN7Dy+ieZxra0pyoN+uf1YVSFNSwizSTxo\nDmlaBvEUD24+sJlmafFodR6GYBbpm9OXdZXrqG+ot7oUIaKaxIPmkKZlEE/xYOmeUrrn5Cu1hdU8\nuTltm7Vl88HNVpciRFSTeNAc0rQM4ikeLP2xlMLczsp9WGW/lhDmk3jQHNK0DOIpHly+Zzl92nWX\npiVEDJJ40BzStAziHg/W1NWwpmIN/fJ6KreFdXbO2dK0hDCZxIPmUOGEuVHBPR5cuW8l3bO6k9ks\nTbkPq4y0hDCfezyo0HpgB3AUqAfswGCgFfBPIN/x+M+Aw9aU1zQZaRnEPR5cunspQzsMVe3DCkCH\nFh2wN9jZd3yf1aUIEbUU3qelAUVAf/SGBfAg8BnQHfjc8b2SpGkZxD0eXFq+lCEdhqj2YQX0M773\ny+0np3MSwkSulytScJ+W+3nlLgNmOb6eBVwR3nL8J03LIO7x4NJyfaSl4IcVkP1aQphN4XhQAxYC\nK4A7HPflABWOrysc3ytJ9mkZxPUDWn60nOq6arq26srx40p9WE/pl9uPeVvk2lpCmMWKeLCkpISS\nkhJfTxsO7AXaoEeCG90e11D4uobStAziGg9+W/4tQzoMwWazKRkPgt60nv76aavLECJqWREPFhUV\nUVRUdOr7Rx991NPT9jr+3Q98iL5fqwLIBfYBbYFKM+sMhcSDBnGNB52TMAASEqC+HurqLCzOg56t\ne7Lj8A65tpYQJlE0HkwDmju+TgcuBNYCHwMTHfdPBD4Kf2n+kaZlENcPqHN/FoDNptzMIQCS4pPo\n0boH6yrXWV2KEFFJ0dmDOcBiYBWwDPgUWAA8DYwBNgGjHd8rSeJBgzg/oDV1NayuWM2g9oNOPebc\nykpPt7BAD5zHaw1uP9j3k4UQAVH0jBjbgX4e7j8IXBDmWoIiIy2DOPdprdq3iu5Z3WmW1OzUYwpt\nZZ2hR1YPNh+QE+cKYQY5I4Y5pGkZxLlVtbR8KUPaDznjMYW2ss7QObMz2w9vt7oMIaKSovFgxJOm\nZRDnB3Rp+VKG5g094zFVt7IKMgrYdmib1WUIEZUUjQcjnjQtA2ja6emtrjMHnVTdypKRlhDmkXjQ\nHNK0DGC3Q3w8bD28idr6Wrq26nrG46puZbVKbUV9Qz2Hqg5ZXYoQUcc1HkxI0GcSq3boSySSpmUA\n54dzxvczmHj2RGy2M0/rpepWls1mk9GWECZxP7WbqhuvkUaalgH0GEBj1upZ3D7g9kaPqxoPAhRk\nyn4tIczgfhJtVTdeI400LQNUV0NDfBW92/SmW1a3Ro+rvIXVOaMz2w/JSEsIo7lfrkialjEi5uDi\nITOG+HxOx5YdGZk/kpGdRlLYppA4W3h6cnU1VGmHuPOcOz0+rvKHtSCzgLLKMqvLECKqOCdnJSWd\nvi85Wd3EJZJETNN6ceyLTT6uaRqbDmziy51f8vy3zwMwffx0ijsXm17blsrd1NqOM6HnBI+PqxwP\nds7szKebPrW6DCGiSm2t3rDiXLabVd54jSQqNK0U4EsgGUgC/g085P6kIR18j7SG5g1lYj/9nI/z\nNs/jln/fwvhu43lmzDM0T27u46eD989V/yar+RUkJyR7fFzleFCO1RLCeO7RIEjTMooK+7SqgVHo\n58Pq6/j6vFBfdFy3cay9ey019TX0faUvOw7vCPUlPaqtr+Xj9Qton5nl9Tkqf1g7ZXRi15Fd1DfU\nW12KEFHDfeYgSDxoFBWaFoDz+hhJQDz6yRtDlpGSwczLZzJp4CRumHMD9nq7ES97hg83fEiH9C5k\nNk/1+hyV48HUxFRapbZiz7E9VpciRNRwnzkIam+8RhJVmlYc+qnyK4BFwHojX/y+YfeRkZLBIyWP\nGPmyHKs5xv2f3c91PW5ptFXlSuV4EOTMGEIYTeJB86iwTwugAT0ebAnMB4qAEtcnTJ069dTX7lfn\n9CXOFscbV7zBgFcHMKrTKMZ0GRNywQB/WPQHijsX06PubFakeH9eSgqcOGHIW5rCeazW+fnnW12K\nEFFB4kHzqNK0nI4A/wEG0kTTCkZ2ejazJ8zmxg9uZOUvV5LTLCek11v+43LeXfcuZZPKmP9R460q\nVykpcNCQwNMccqyWEMaSeNA8KsSDrYEMx9ep6FfPXGnGG40uGM1t/W/j1n/fiqZpQb+Ovd7OnZ/c\nyZ8v/DNZaVkeowBXqn9YCzIL2HZYZhAKYRSJB82jQtNqC3zB6cs/fwJ8btabPTLyESpPVDL9u+lB\nv8afl/6ZnGY53NDnBsBzFOAqIvZpyUhLCMNIPGgeFeLBtcCAcL1ZYnwib054k/PfOJ/izsWNzsju\ny+zVs/nb8r/x9a1fnzoxrqcowJXKswdBP1ZLJmIIYRyJB82jwkgr7Hq16cUfzv8DN394M3UN/l8r\nYObKmUz5fAqf3/w5BZkFp+6P9HiwXfN2HDh5gCp7ldWlCBEVJB40T0w2LYDJgyeTlpjGk4uf9Ov5\n07+bziMlj/DFxC/o2brnGY9FejwYHxdPx5Yd2Xlkp9WlCBEVJB40T8w2Lec0+LfWvMUNc27gwMkD\nHp+39eBWrn3/Wp76+ikWTVxE96zujZ4T6fEg6Pu15HROQhhD4kHzxGzTAujQogOr7lpFbrNc+kzr\nwwcbPmDn4Z2UVZbxbfm33DPvHs6dcS59s/uy7u51Xvd/RXo8CI79WjIZQwhDSDxoHhUmYlgqLTGN\n5y96ngk9JzB53mQOVx+mWVIzmic159z257L+V+vJTs9u8jUiPR4EGWkJYSRpWuaJ+ablNCJ/BKvv\nWh3Uz0ZDPFiQWcDS8qVWlyFEVKipkX1aZonpeNAo0RAPykhLCOPISMs80rQMEC3x4PbD20M6U4gQ\nQidNyzzStAwQDSOtjJQM4mxxHKxS+CSJQkQImfJuHmlaBoiGfVogEaEQRpEp7+aRpmWAaIgHQa6r\nJYRRJB40jzQtA/iKB5OTobYWGhrCV1MwCjIKZKQlhAEkHjSPNC0D+IoHbbbTjUtlEg8KYQyJB80j\nTcsAvkZaEBkfWIkHhTCGxIPmkaZlAF/7tCAy9mtJPCiEMSQeNI80LQP4igchMmYQ5mfkU360PKDL\ntQghGpN40DzStAwQLfFgUnwSOek5lB8tt7oUISKaxIPmkaZlgGiJB0EmYwhhBIkHzSNNK0R1daBp\nkODj1MOREA+CNC0hjOApHkxKArsd6uutqSlaSNMKkfPDabM1/bxIiQbkulpChM5TPGizRc7Gq8qk\naYXIn2gQIiwePCwjLSFC4W29ECnrAZVJ0wqRPzMHIXK2sCQeFCJ03tYLkbIeUJk0rRD5M3MQIige\nzJR4UIhQaJr3kVakrAdUJk0rRNEWD+ak53DCfoJjNcesLkWIiFRXp++/8jQ5K1LWAyqTphWiaIsH\nbTabPhlDTuckRFCaWidEynpAZdK0QhRt8SBIRChEKJpaJ0TSekBV0rRCFG3xIEDnDJmMIUSwmlon\nRNJ6QFXStEIUbfEgyAxCIUIh8aC5VGhaecAioAxYB/yPteUEJlrjQTlWS4jgREA8OBbYCGwGHrC4\nloCp0LTswG+A3sAQ4FdAL0srCkBUxoOZnWWflhBBUjwejAdeQm9chcD1RND6FtRoWvuAVY6vjwMb\ngHbWlROYaIwHnbMHG7QGq0sRIuIoHg8OBrYAO9AHDO8Cl4fpvX8PTAQuDeVFfJzmFYAuwNVAS2A3\nsJTTTcZonYD+wLIz7x7PJ5+Y9I4hKi31v2lt2oQlv0eXLlBY6P/z05PSSd0/nNn/OkRWapYpNfXr\nB3l5wf/8li2wYYNx9QQrPx/69m18f0UFLF8e/nqaYrPBqFGQnt74sW++gYMHw1+Tu8GDISen6ecc\nOABLloSnnmCsWdN001q+HDIzzXnv9HQYPbrJp7RHX487lQPnmlNNI32AZwEtlBfxp2mNBuYA+cDd\nwCVAW+B54O1Q3txNM+B94B70EZeLO/nd73449V1WVhZZWa0NfOvQ3HST7+cMHKg3rOnTza/H1cGD\nEB8PX30V2M9VvfMGf10eR3sTetbOnXDeefDyy8G/xv33Q3k55OYaV1egjh7VbytXNn7s73+H996D\nbt3CX5c3330Hf/kLXHPNmfdrGhQVwUUX+T7xs5k2b4af/Qwee6zp573+uv531EvhUOuyyzzfP2YM\nvPmmeesBTSvhq69KmnyKOe/sly8d758BXAQcBkoBwzeXfgn0c3w90fFvMnoWerdB75EIzAfu9fK4\nJoKzbJmmDRwY+M+ltTqkPfbxTOML0jTt9dc17ZZbQnuNiy7StLlzjaknWGVlmtajh+fHfvc7TfvT\nn8Jbjy833qhps2Y1vr+2VtPi48Nfj7s//Ulfbr5Mnappf/iD+fVEAxo3qSHAf12+f4jwTca4w+Vr\nG/Aa+ugrIP7s05oOFAELgSuA8eiR4XL00VGobMDrwHrgRQNeT7gINkO31Sez7dh64wvCmFzf332J\nZmrq91ChPnfe6lWlVn8/F6rUG6FWAN3Qd8UkAdcCH4fpvf8X+BB4BL2XbATWOh7zO6L0Jx7U0JvJ\n39GjwqHoO9IOAO/4X69Xw4EbgTWAM2h5iDO3BkSQgp1iW29PZPPRtb6fGAQjpv36e6iBmZr6PVSo\nz523elWp1d/PRXU1tFZn70CkqQMmoydb8egDhnDtHf5f9PkKQ9DnSZyLvjvoSyAduNKfF/GnaTnZ\n0X/R+QGV6dvXqDGLMSoF2yDsNfFsPGzOfBtpWtaIpqalQr0RbJ7jFm5/cfz7rct9rdFnNE7290UC\naVoiAgVzXEhdnb5zvsFWS+WJSrLTsy2vyZ2/x8eZKTm56XjQ6vrceatXlVqbWp6uVKlXGOInYC4B\nTMaQEU6UC2b/kb7PwEbv7N6s32/8fq1o2aflbL6ah/lYKo4GZKQlFPat76fopGlFuWCiOOdKobB1\nIWWVZUrU5E6FFVdCgn44gd3e+DEV6nMnTUtEA2laUS4pCWproSGAk1s4ozezRlrREg+C999Flfpc\nqV6rv58LVeoV1pCmFeVsNv0PvLbW/59xRm+FbQop22/OSCsa4kFQfxq5K9VrlSnvwh/StGJAoHGc\nM37p3ca8fVqhjLQ0TZ2tbdUjN1eq1yrxoPCHNK0YEGgc52wIuc1ysTfY2X9iv6X1uKur00eQCQrM\nfVU9cnOleq0SDwp/SNOKAYHGcc74xWazUdim0PDRVqjxoErxkOqRmyvVa5V4UPhDmlYMCDYeBD0i\nNHq/VqjxoErxkOqRmyvVa5V4UPhDmlYMCDYeBEwZaSUkQH29fguGSvGQ6pGbK9VrlXhQ+EOaVgwI\nNh4Ec0ZaNltoEaFK8ZDqkZsr1WuVeFD4Q5pWDAglHjRjpBVMTa5Uioc8/R4qzW50pXo82NQZRlyp\nUq+whjStGBBKPNiueTtq6mr46eRPltbkSqWG4On3sNv1M2XEx1tTkzeqx4POZebpDCOuVKlXWEOa\nVgwIJR5UcQahSvGQp99DpfpcqR4Pgu/PRUODfshDUlL4ahJqkaYVA0KJB8GxX8vgcxBGczyoUn2u\nVI8HwffnwnmGd5stfDUJtUjTigGhxIMAfXP6srpitaU1uVIpHvLWtFSpz5Xq8SD4/lyoVKuwhjSt\nGBBKPAgwsN1ASveUWlqTK5XiLE/XgFKpPlfREA+qVKuwhjStGBBqPHh27tls2L+B6roQT80eQk2u\nVI+zVKrPlbPBus/OU6leX58LlWoV1pCmFQNCjQfTEtPoltWNtRVrLavJlUoRUSTFg3Fx+oHd7mf8\nV6leiQeFL9K0YkCoIy2AQe0GsWLPCstqcqXS1nYkxYOg/mxHiQeFL9K0YkCo+7RA369ldNOKhn1a\nkRQPgvr1SjwofJGmFQNCjQfB0bT2Gte0JB60hqflrlK9Eg8KX6RpxQAj4sE+2X3YfGAzJ+0nLanJ\nlUpb2xIPGkviQeGLNK0YYEQ8mJyQTGGbQlbvM+Z4LYkHraF6vRIPCl+kacUAI+JBMPZ4LYkHrSHx\noIh00rRigBHxIBg7GUPiQWtIPCginTStGGBEPAjGNy2JB8NP9XolHhS+SNOKAUbFg73b9GbnkZ0c\nqzkW9ppcqRQRSTxoLIkHhS+qNK2ZQAVg3CkXxClGxYOJ8Yn0zenLyn0rw16TK5W2tiUeNJbEg8IX\nVZrWP4CxVhcRrYyKBwEGtjUmIpR40Bru9dbX69enSky0riZXEg8KX1RpWouBQ1YXEa2MigfBuP1a\nEg9aw325q3Z9KokHhS+qNC1hIqPiQdCb1vIfl4e9JlcqbW1720ekSn3u3Je7arVKPCh8SbC6AH9N\nnTr11NdFRUUUFRVZVkukMTIe7NWmFwerDrLv+D5ym+WGrSZ/6ws31fcRuXOvV7VaJR4UvkRk0xKB\nCWRU09AAdjskJXl+PM4Wx7C8YXyz6xuuKrwqLDW5U2nFFen7tFSrVZqW8EXiwRgQyP6jmhq9YTW1\nj+O8jufx9a6vw1aTO5X2a0T6Pi3VavU0G9OVcx+ciF2qNK13gCVAd2A3cKu15USXQEY1/mzJntfx\nPL7eHVrTipaRVkKCPjqtqzt9n2qRmyuJB0WkUyUevN7qAqKZ62XWfc0S82clNrDdQNbvX8/x2uM0\nS2oWVE3Rsk/LZjv9uyQ4/ppUXrFKPCginSojLWGi+Hj9Zrf7fq4/cVFKQgr9c/uHNIswWuJB8NwI\nVKrPVSTEgzLlXTRFmlaM8DeO83dLNtT9WtESD0Lj/TAqjQTdRUI8KFPeRVOkacUIf+M4f1cKRjSt\nYOLB+nr9psoZHED9yM2V6rVKPCh8kaYVI/yN4/yNX4blDePb8m+pa6jz/eQQ6nGn2hkcQOJBI0k8\nKHyRphUjjI4HW6W2omPLjqypWGNqPe5U3NKWeNA4Eg8KX6RpxQij40GA4XnDg44Ik5P1iSGaFtjP\nqbjSUj1yc6V6rRIPCl+kacUIo+NBCG2/ls2m75cKdL+WivGQxIPGkXhQ+CJNK0YYHQ+C3rS+2f0N\nWqDDpQBrcqXilrbEg8aReFD4Ik0rRpgRD3bK6ATA9sPbTa3JlYorLdfmW1ennyEjQZXD9t2oHg96\nOsOIK9XqFeEnTStGmBEP2mw2RuaP5MsdX5pakysV4yHXRqDi7EZXqseDNpv38w9qmpx7UEjTihlm\nxIMARZ2KKNlZYmpNrlTc0nZdyao4EnSlejwI3kfgdrs+EouTtZZZpgLlwErHbZzLYw8Bm4GNwIVh\nr8yF/PfHCDPiQXA0rR0lQe3XisZ4UMWm6kr1eBC8b8yoWGuU0YDngf6O2zzH/YXAtY5/xwIvY2Hv\nkKYVI8yIBwG6teqGvd4e1H6taIwHVazPlerxIHj/XKhYaxTyFGxfjn4lDjuwA9gCDA5jTWeQphUj\nzIoHbTbbqdGWWTW5UnFr23XEqOJI0FUkx4Mq1hqFfg2sBl4HMhz3tUOPDZ3KgfZhrusURec4CaMF\nEg+2bBnYazub1m39bzOlJlcqrrhcRwYqNlVXzlmNdXX61yrWK/Fg8EpKSigpKWnqKZ8BuR7ufxiY\nBjzm+P5x4M/AL7y8TnDHuRhAmlaMMCseBBjVaRRPfPUEmqZhC2DaXDTFgydO6F+rWJ8753Jv1kzN\neiUeDF5RURFFRUWnvn/00UfdnzLGz5eaAXzi+PpHIM/lsQ6O+ywh8WCMMCseBOjaqiv1Wj3bDm0z\npSZXKm5tR1I8COrXK/GgZdq6fD0BWOv4+mPgOiAJKAC6AcFfTC9EMtKKESkpcOiQ7+cF0xRc92t1\nadUloJo51e/UAAAWaklEQVSioWlFUjwI6s92lHjQMn8C+qFHf9uBXzruXw+85/i3DpiEhfGgjLRi\nhLcDNt0Fe/BmUX7gx2v5W5MrFQ8ujaTZg9C4yapWr8SDlrkZ6AucDVwBVLg89iTQFegJzA9/aadJ\n04oRZsaDENzxWtEy0lI9bnOner0SD4qmSNOKEWY3ra6tuqJpWkD7taKlaUk8aCyJB0VTpGnFCLPj\nQed+rUU7FhlekyuJB0Mn8aCIZNK0YoTZIy3QI8Ivtn9heE2uVNzaVn3k4k71eiUeFE2RphUjwtG0\niguK+WL7F37v14qWphVJJ8yFyNinJfGg8EaaVowwOx4EKMgsIC0xjbL9ZYbW5EriwdBJPCgimTSt\nGBGOkRboo63Pt31uaE2uVNzaVj1uc+esV9XrU0k8KJoiTStGhK1pdS7m8+2x1bQiNR6srYWkJPWu\nTxVsPHi89jjTSqfx2nevsWT3Eo5UHzGvSGEZOSNGjAhHPAgwumA0d316F3UNdSTENf3xknjQGs74\nTdVam4oHMzMb33+o6hB/W/43Xlr+EiPyR9A8qTnTv5/Ohv0baJnSkl6te9GzdU/6ZPdhQq8JZKdn\nm/9LCNOo0rTGAi8C8egnavyTteVEn3CNtLLTs+nYsiMr9qxgSIchhtRkZH1m8BUP1tbX8vm2zynd\nU8p3e7+jrLKM6866jikjppCWmGZZvSouS/A/HtQ0jdmrZ3Pfgvu4tMelLL51MT1a9zj1eIPWwO4j\nu/nhwA9s/Gkji3ct5oGFD1DUqYhb+t3C+O7jfW5YCfWo8D8WD7wEXIB+5uBS9BM0brCyqGjjT4Mw\nah+Hcxahyk3rp5M/8XLpy/Ru05sxXcbQIrlF0K/lLR6sPFHJqyteZdqKaXTO7MzwvOHc2OdGOmd2\n5tklz9L75d68eNGLXNbjMp9nx7fX22nQGkhOCH1o5GwKqkaZ/sSDR2uOMuk/k1i5byWLJi6iT06f\nRs+Ps8WRn5FPfkY+F3bRrxB/rOYY769/n2e+eYYHFz7I46Me56rCq4izNZ2RaprGziM7OVpzlLOy\nz/L5fGEeFZrWYPQrYe5wfP8u+pUypWkZyJ9rV9nt+jWW4uNDe6/izsW88O0LTBkxJeSa3IW6otU0\njffK3uPe+fcyrus4luxewi3/voUhHYZw/7D7T63cAuE+0kpIquPBhf+PV797lWsKr2HBTQs4K/us\nM37m3avf5fNtnzN53mQe/fJRiguKKepUxIC2AzhUfYi9x/ay59geVu5bybIfl7F632psNhvD84Yz\ntutYLul2Cd2yugW1DCJhpNVU01pTsYYr/3klxQXFlN5RGtBotXlyc27tfyu39LuFz7Z9xpTPp/D0\nN08zedBkerTuQUFGAW3S27Dj8A7KKsso219G6Z5Svi3/Vv/5pOYcrj5McedixnQew8XdLia3mafL\nU+ns9XZ2HN5B11Zd/b5sT4PWwOKdizlUfYjuWd3pktnFkI2VaKFC02oP7Hb5vhw416JaopY/164y\nah/H+fnnc/2c66myV5GamBpSTe5CqfGnkz/xi49/wdaDW/no2o84t4P+MTtee5z/bPoPt/37Nm7q\nexOPjXqMxPhEv183OVmf1NDQAIeOVzF18R/pMnQ9G3+1kZxmOV5/rrhzMWvuWsOyH5exaPsinv/2\nedZUrCErNYu2zdvStllb+mT34fFRjzOw3UAatAa+2P4F/93yX576+iku6XYJTxY/2eRK01u9x46p\nvU/LWzx4qHYvF711Ec+NeY6f9/150O9hs9m4sMuFjOk8hjkb5jBnwxxe/e5Vth3axoGqA+S1yKN3\ndm8KWxdybe9refGiF+nYsiM2m41dR3axcNtC5m+dz30L7qN7Vncu7X4peS3yOGk/yUn7SfYc28Oy\nH5exat8qmic3p2VyS3416FdM7DfR46i+QWtg/f71vLP2Hd5a+xYZKRnktchj88HN7Dy8k7yWeZzT\n9hwGtRvEwHYD6dWmF23S2gR0/bpoocJvfBX6Pq07HN/fiN60fu3yHC2QE7GKxux2SE3Vr1jrTWUl\n9O4N+/eH/n7DXh/G46Mep7hzsdfnrFwJt94Kq1b5/7o5ObB6NeQGtp7mYNVBRs8azeiC0Tx9wdMk\nxSc1ek7liUomfjSRozVHmX3FbDpndva5UqhvqNdXclmtmbb4Xe76eXtuvnsf0++9xtQI6WjNUZ74\n6glmrpzJ/cPu5zdDf+Pxd/Lk+eehvByuvx4mTYLSUtPKDMpnn8Ezz+j/uhp3aRXfZd/Lk78axO0D\nbjft/f2ZRORUW1/L4p2L+XTTpxysPkhaQhppiWm0TmvN4PaDGdR+EM2TmvP1rq95qfQlPtv6Gf3b\n9icnPYec9BxsNhvf7/2elftWkpWaxZW9ruSmvjdxdu7Zp97DXm9n88HNrNiz4tTthwM/UNdQR7dW\n3RjXdRz3D7/fYzN0fH5VWM8bRoWRlvtVMfPQR1tnmDp16qmv3a/OKXxLSND3WTkvs+6JkXFRcYE+\n9b2pphWufVpHa44y9q2xFBcU89yFz3ltRNnp2fznhv/w7DfPMui1QdTU19CxZUfaNW9HTV0Nh6sP\nc6j6EFX2Kuoa6qhrqMPeYCcjJYM62zZmf/cevTJmMXFgS+JMXk20SG7BM2Oe4Y4BdzBp7iS+2vUV\n71/zfpMjWyez48E9x/Ywc+VM1lauZevBrWw5uIW0xDQGtB3AOW3P4fz88xldMNrr/4Onz8WhqkN8\nva2Mq84fze0DrjW+aBeBTM5Iik+iuHNxk59zgBH5IxiRP4K9x/aytnItlScqqTheQV1DHVNGTOGc\ntueQlZbl8WcT4xMpbFNIYZtCbj775lP3H6w6yA8//cAr371C9791548j/8gdA+4IKCUQwUkAtgKd\n0K+MuQro5fYcTYQuNVXTjh/3/vimTZrWpYsx77Vo+yJt0PRBTT5n61ZN69QpsNdNTta0kyf9f/7x\nmuPaeTPP0+7+9G6toaEhoPc6Wn1UK6ss0xZsWaB9teMrbW3FWm3X4V3agZMHtCPVR7STtSe1uvo6\nTdM0LTtb0/bu1bSBAzVt2bKA3iZktXW12nXvX6eNnjVaO1ZzzOfzZ8zQtFtv1bQFCzStuDi492xo\naNDKKsu0xTsXaxv3b9QOnDygra1Yq93y0S1a5tOZ2qRPJ2lvr3lbW1a+TPvpxE/ajkM7tDnr52gP\nf/6wdtbLZ2mDpg/S5m6a6/H/ZNkyfTk6lVWWaQNeHaB16LtZW7gwsP/DWPH9nu+14lnF2tAZQ89Y\nplh4sUazqDDSqgMmo19YLB54HZmEYQrnFmx6uufHjdzyHpY3jJ1HdrL5wGavEwYCHWkFOrvxRO0J\nLn3nUrq26spLF78UcP7fPLn5qS1cX6yckZcYn8hbE97izk/u5KK3LmLuDXNpmdLS6/M91bpk9xK+\n3vW1Pi1/z3ckJyTTJ7sPfbL70LVVV1okt6BZUjNsNhvzt8zn/Q3vU2WvokOLDuw/uZ/9J/aTmpjK\n5EGT2fzrzY1GDVlpWeRn5HNlryt5bNRjzFk/h9999jse++oxfj3410zoOeHUKPFUfXU1PLn4SV5e\n8TKPFj3KW+ldSE2NqqTLMP3b9uezmz5j99HdUb+fS4WmBTDPcRMm8tUkjGxaSfFJ3Nz3ZmZ8P4M/\njfF82F2gTaumxv8zOByvPc74/xtPQWYBMy6dYfoUZdcDdq2YkRcfF89rl73Gb+f/lnNnnMusK2ad\nmmjizjUejEus4fo5t7D8x+WM7zaey3tczuOjHqe2vpZ1letYW7GWj374iGM1xzhWe4za+lqK8ot4\nc8KbDGo3KKgVZJwtjmt6X8OVva7kw40f8tr3rzF57mSu7X0tfXL6ULbBzs4D19Fn2gjOyj6LVb9c\nRfsW7Zmh6GxHVdhsNjq27Gh1GaZTpWmJMPB1BgqjzzZx+4DbGfnGSB4f/bjHSQKBnhHD3/qO1Rzj\nkv+7hG6tuvHaZa+F5Zga10Zg1Yy8OFscL459kX+V/YvL372c2/rfxiMjH2k0XdrZYJduX8nCndv5\nZXou6+5e12h/WN+cvtD48CfDxMfFc3Xh1VxdeDW7j+xm9urZrN63mtSknsTXpzHjshmcn3/+qeer\nOttRhJc0rRgSzpEWQI/WPeie1Z1PfviEqwqvavS4c+WpaeDPBrs/9Z20n2Tc2+MobFPIK+NfCdtB\noCodsHtN72s4P/98fvnpL+kzrQ/DOw6nc0ZnCjILqDhewdzSGr7dUsTyFh9T1O0eXhh7pbUFA3kt\n83j4/IcBqKiAdzTOaFigxrIV1pOmFUPC3bQA7hhwBzNWzvDYtJwHMtvteuzni6/66hvqufGDGynI\nLAhrwwLr40F3Oc1y+PDaD1lavpQN+zew9dBWPt30qT4VO38cBzPO5tbh5/LDhhCPJDeBXE9LNEWa\nVgwJdzwIcHXh1dw7/152Ht5Jfka+15r8aVq+6rv/s/s5VH2Id69+N+yn2VEhHnRns9kYljeMYXnD\nzri/tBQ+qwe7gicfBrmelmianEArhlgx0kpNTOWGs25g5sqZQdXkqqn6Xlr+EvO2zOODn33g90G2\nRnL+Hna7+itWlaJMT5KT9eXofj4BVesV4SVNK4ZY0bRAn5Axc9VM6hvqA67Jlbf6vtzxJU8ufpK5\nN8wlM9XDtSvCIDkZjhzRR4yqzzhW/dyDNpu+HN1TAVXrFeElTSuGWBEPApydezbtm7fngw0fBFyT\nK0/1NWgN3LfgPl4c+yIFmQUGVBuclBS9aak+ygL1r6cFjSPCujp95OXtbC4idkjTiiFWjbQAphZN\n5Q+L/kBdw5knPwx1pDVn/RxA33dmJWfTioSRgOrxIDS+AoDKtYrwkqYVQ6xsWhd1uYi2zdsya9Ws\ngGpy5V6fvd7Ow188zFPFT1l+faPkZDh8ODJWrKrHg9D4c6FyrSK8pGnFEKviQdBnsj1V/BRTv5xK\nlb3K75qaqu+NVW+Q1zKPCzpfYHC1gZN40Fju8aDKtYrwkqYVQ6wcaQEM6TCEge0G8vfSv/tdkyvX\n+qrsVTz65aM8VfyUEudai6R4MClJ30dUVaVuvRIPCm+kacUQq5sWwBOjnuCZb57hcPVhv2py5Vrf\n30v/zrkdzmVw+8EmVRqYlJTIiQdtttOzHVWtV+JB4Y00rRhiZTzo1Du7N5d0v4TnljznV02unPVV\n2at4bslzPFb0mImVBsa5TytSIizV65V4UHgjTSuGqDDSApg6cirTVkyj8kRlUCOtmStnMqTDEHpn\n9za30ABEUjwI6tcr8aDwRppWDFGlaeVn5HPDWTfw1OKnAm5aiUn1PLvkWR4870FziwyQ6k3Aner1\nSjwovJGmFUO8ndPNKZwRzMPnP8zsNbOx244H1LQ2H11HQWYBQzoMMbfAAKket7lTvV6JB4U30rRi\niHvk4i6cEUxus1zuHHAnq3/61u99WtXVGt/s+ZwHh6s1yoLTI4NIGQ2oXq/Eg8IbaVoxRJV40On+\n4fez5VgZPx484NfzN+7bTmJSPRd2udDkygLnXG6RsmJVvV6JB4U30rRiiErxIECr1FYMye/P/B9K\nfD5X0zSW71rLuJ6jlTguy51zuUVKhKV6vRIPCm+kacUQleJBp5FdzuXHQwd4YekLXp+jaRr3/Pce\nNHsiQzv1D2N1/lN95OJO9XolHhTeSNOKIarFgwAtmyUzoet1/HX5X3l1xauNHtc0jfsW3MfS8qWc\n1+4C0tLU/Miq3gTcqV6vxIPCGzXXAMIUqsWDoL9fotaChTct5InFT/Dm6jdPPaZpGg8sfICSHSUs\nuHEB9fYkZSMi1eM2d6rXK/Gg8EauThNDVIwHnTV1adWFBTcuoHh2MVO/nEpNXQ1VdVV0bNmRhTct\nJDM1U+mISPWRi7uUFIiPV/f6VCkpcPz46e9V/r8X4aXoR1aYQcV40LWmXm16UTapjINVB0lJSCE5\nIZmMlAwS4hIsq89fkdi0VK7VUzzYrJl19Qh1SNOKIarGg641ZaZmkpma6fG5KkdEqsdt7pKT1a7V\nUzyYlWVdPUIdsk8rhqgcD/pD5YhIRlrGktmDwhtpWjGkqXiwrg40Lfz7OEK5crFKpGkZS2YPCm+k\nacWQphqEc6UQ7uN2o6VpOaM2VetzJ01LRCqrm9Y1QBlQDwywuJao19S1q8JxLS1PgrmeloqcM/FU\nrc9dJOzTco8HVa43SjS1Pn4I2AxsBFzPo3YOsNbx2F/CUKPlTWstMAH4yuI6IkZJSUnQP5ucDLW1\n0NDQ+DGrtmRDGWmFsizMYOXoJdBlEc0jLdU+FxHE2/q4ELjW8e9Y4GXAmclMA34BdHPcxppdpNVN\nayOwyeIaIkoof5DOy6x7GtlI0wqdNC3jSNOyhLf18eXAO4Ad2AFsAc4F2gLNgeWO580GrjC7SKub\nlggzb9PerZpO7msavpOmqT3lHdSP3FypXqucEUMp7YByl+/LgfYe7v/Rcb+pwjFX7DMg18P9U4BP\nwvD+wkVaGlx7beMVwJEjkJpqTT3798Ollzb9PE1T+wwOoP8uaWlWV+Ef1WtNS4OtW09/Ltats+bz\nGYVkfWyQRTQ9EWMLoMlNbnKTm9wCum0hcIs4c338oOPm9F/0eDAX2OBy//XAK0G8X0RahD4LRQgh\nhLXc18eFwCogCSgAtnJ6IsYy9AZmA+YShokYVpsA7AaqgH3APGvLEUKImNXU+ngK+qhtI3CRy/3O\nKe9bgL+Gp0whhBBCGGIsemffDDxgcS3hloc+TC8D1gH/47i/FfrO1E3AAiDDkuqsEQ+s5PQO41hd\nFhnA++j7E9ajxzOxuiweQv8bWQv8H5BM7CyLmUAF+u/u1NTv7u0AYWGQePQhZycgET1T7WVlQWGW\nC/RzfN0M+AH9938G+L3j/geAp8NfmmV+C7wNfOz4PlaXxSzgNsfXCUBLYnNZdAK2oTcqgH8CE4md\nZTEC6M+ZTcvb7+7cL5WIvty2IIc8GW4o+iwVJ/cZLLHmI+AC9K2kHMd9uY7vY0EHYCEwitMjrVhc\nFi3RV9TuYnFZtELfmMtEb96fAGOIrWXRiTOblrff/SHOTKv+CwwxuzgzqNxp26PvFHRyHtAWizqh\nb1EtQ/9AVjjur+D0BzTavQDcD7iehCoWl0UBsB/4B/A98BqQTmwui4PAn4FdwB7gMHo0FovLwsnb\n7+7tAOGIo3LT0qwuQBHNgDnAPcAxt8ecx2JEu/FAJfr+LG/noY+VZZGAfgzNy45/T9A4gYiVZdEF\nuBd9o64d+t/KjW7PiZVl4Ymv3z0il4vKTetH9MkITnmcuaUQCxLRG9ab6PEg6FtPziPa26KvzKPd\nMOAyYDv6OdBGoy+TWFwW5Y5bqeP799Gb1z5ib1kMBJYAB4A64AP03QqxuCycvP1NuK9POzjuizgq\nN60V6GcN7oR+UNu1nN4BHwtswOvos8NedLn/Y/SdzTj+/YjoNwX9D64AuA74AriJ2FwW+9Bj8+6O\n7y9Anz33CbG3LDai75dJRf97uQD97yUWl4WTt7+Jj9H/dpwHCHfj9IluhYHGoe9o3YK+IzGWnIe+\n/2YVeiy2Ev0QgFboExKifTqvNyM5vfESq8vibPSR1mr00UVLYndZ/J7TU95noacTsbIs3kHfl1eL\nviFzK03/7t4OEBZCCCGEEEIIIYQQQgghhBBCCCGEEEIIIYQQQgghhBBCCCGEEEII9fx/YRnEyVu/\nl0gAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0xaf6cc6f0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot data\n",
    "fig, ax1 = plt.subplots()\n",
    "ax1.plot( x[ :,1 ], 'g-' )\n",
    "ax1.set_ylabel( r'$\\theta$' )\n",
    "\n",
    "ax2 = ax1.twinx()\n",
    "ax2.plot( x[ :, 5 ], 'b' )\n",
    "ax2.set_ylabel( '$F$' )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "g1 = MVN( x[ 10, 1:], np.array( [ 0.5, 0.5, 0.5, 0.5, 0.5  ] ) )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "x_dyn = np.column_stack( ( x[ 1:, 1:5 ], x[ :-1, 1:6 ] ) )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "GMM(covariance_type='full', init_params='wmc', min_covar=0.001,\n",
       "  n_components=10, n_init=1, n_iter=100, params='wmc', random_state=None,\n",
       "  thresh=None, tol=0.001)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "gmm = mixture.GMM( n_components=GMM_DYN_COMPONENTS, covariance_type='full' )\n",
    "gmm.fit( x_dyn )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "g1 = MVN( x_dyn[ 13, :], np.diag( x_dyn[ 13, :] * INITIAL_VARIANCE_RATIO ) )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "ValueError",
     "evalue": "operands could not be broadcast together with shapes (5,) (4,) ",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mValueError\u001b[0m                                Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-10-a62bc3925394>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m()\u001b[0m\n\u001b[1;32m----> 1\u001b[1;33m \u001b[0mg\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mg1\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mcond\u001b[0m\u001b[1;33m(\u001b[0m \u001b[0mx_dyn\u001b[0m\u001b[1;33m[\u001b[0m \u001b[1;36m13\u001b[0m\u001b[1;33m,\u001b[0m\u001b[1;36m4\u001b[0m\u001b[1;33m:\u001b[0m \u001b[1;33m]\u001b[0m \u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[1;32m/home/user/work/gps/mvn.py\u001b[0m in \u001b[0;36mcond\u001b[1;34m(self, a)\u001b[0m\n\u001b[0;32m     49\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     50\u001b[0m                 \u001b[0mt\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mdot\u001b[0m\u001b[1;33m(\u001b[0m \u001b[0mcov12\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mnp\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mlinalg\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0minv\u001b[0m\u001b[1;33m(\u001b[0m \u001b[0mcov22\u001b[0m \u001b[1;33m)\u001b[0m \u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m---> 51\u001b[1;33m                 \u001b[0mmi\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mmi1\u001b[0m \u001b[1;33m+\u001b[0m \u001b[0mnp\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mdot\u001b[0m\u001b[1;33m(\u001b[0m \u001b[0mt\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0ma\u001b[0m \u001b[1;33m-\u001b[0m \u001b[0mmi2\u001b[0m \u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m     52\u001b[0m                 \u001b[0mcov\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mcov11\u001b[0m \u001b[1;33m-\u001b[0m \u001b[0mnp\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mdot\u001b[0m\u001b[1;33m(\u001b[0m \u001b[0mt\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mcov21\u001b[0m \u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m     53\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mValueError\u001b[0m: operands could not be broadcast together with shapes (5,) (4,) "
     ]
    }
   ],
   "source": [
    "g = g1.cond( x_dyn[ 13,4: ] )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}