{
 "metadata": {
  "name": "",
  "signature": "sha256:cca7df3778414be6ed71f2d5da083a906dea73e58bd2125ed98a63bf0b3a1da0"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "# Low Rank Approximations\n",
      "\n",
      "# Gaussian Process Summer School, Melbourne, Australia\n",
      "### 25th-27th February 2015\n",
      "### written by Neil D. Lawrence and James Hensman\n",
      "\n",
      "In this lab we are going to consider low rank approximations to Gaussian process models."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%matplotlib inline\n",
      "import numpy as np\n",
      "import matplotlib.pyplot as plt\n",
      "import GPy\n",
      "import pods"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Populating the interactive namespace from numpy and matplotlib\n"
       ]
      }
     ],
     "prompt_number": 1
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## Low Rank Approximations\n",
      "\n",
      "In the worst case, inference in a Gaussian process is $\\mathcal{O}(n^3)$ computational complexity and $\\mathcal{O}(n^2)$ storage. For efficient inference in larger data sets we need to consider approximations. One approach is low rank approximation of the covariance matrix (also known as sparse approximations or perhaps more accurately parsimonious approximations). We'll study these approximations by first creating a simple data set by sampling from a GP."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "X = np.sort(np.random.rand(50,1)*12,0)\n",
      "k = GPy.kern.RBF(1)\n",
      "K = k.K(X)\n",
      "K+= np.eye(50)*0.01 # add some independence (noise) to K\n",
      "y = np.random.multivariate_normal(np.zeros(50), K).reshape(50,1)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 21
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Build a straightforward GP model of our simulation. We\u2019ll also plot the posterior of $f$."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "m = GPy.models.GPRegression(X,y)\n",
      "m.optimize()\n",
      "fig = plt.figure()\n",
      "ax = fig.add_subplot(111)\n",
      "m.plot_f(ax=ax)\n",
      "m._raw_predict?\n",
      "mu, var = m._raw_predict(X) # this fetches the posterior of f\n",
      "\n",
      "plt.vlines(X[:,0], mu[:,0]-2.*np.sqrt(var[:,0]), mu[:,0]+2.*np.sqrt(var[:,0]),color='r',lw=2)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 22,
       "text": [
        "<matplotlib.collections.LineCollection at 0x1115ec990>"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAW0AAAD7CAYAAAChScXIAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3Xd4XNWZ+PHvnaZR7733Llly78ZgY1NiagKB9GSzqZuQ\nyqawyWYDvySbTU/YJGTTCCQYCLbBYGzcsNyLeu+9d2navb8/RhY2lmSV0RTpfJ7HD7Z0594jY706\n855z3ldSFAVBEATBNagcPQBBEARh9kTQFgRBcCEiaAuCILgQEbQFQRBciAjagiAILkQEbUEQBBei\nWewHSJIk9hQKgiDMg6Io0rs/tqCZtiRJekmSzkiSdFmSpFJJkp6c5sGTv5544onr/uzMv8RYl/dY\nXWWcYqxLc6zTWdBMW1GUcUmSblEUZVSSJA1wUpKkTYqinFzIfQVBEISpLTinrSjK6MRvdYAa6F3o\nPQVBEISpLThoS5KkkiTpMtABvKUoSulM12/btm2hj7QbMdbF4SpjdZVxghjrYnHGsUoz5U7mdCNJ\n8gVeB76uKMrRaz6uPPHEE5PXbdu2zSn/IgRBEBzp6NGjHD16dPLP3/nOd1CmWIi0WdAGkCTpW8CY\noig/uuZjii2fIQiCsBxIkjRl0F7o7pEgSZL8Jn7vDuwALi3knoIgCML0FrpPOxz4oyRJKqw/AP6s\nKMrhhQ9LEARBmIpN0yNTPkCkRwRBEOZsUdIjgiAIgn2JoC0IguBCRNAWBEFwISJoC4IguBARtAVB\nEFyICNqCIAguRARtQRAEFyKCtiAIgh3JioK8gLMri965RhAEYbm7VNlBU9cgBpM88REJjVpCq5KI\nCfUmKzEElXTDOZopiRORgiAIi6R3cIxD5+sJDfIl2N9zymu6+0bo6Bki0MeNLbnRaDVqYPoTkSJo\nC4IgLIKmjkHeLm4hKzkcaRaz6DGDiZrGbvw8tWxfGY9OqxZBWxAEwR76hsZ4/VwD2cnhc36twWSm\nqbWPezYni9ojgiAIi81ikTl4po6spLB5vd5NqyEuMnDaz4ugLQiCYENvnq8jKTZ4VimR+RBBWxAE\nwUa6+0cYMch46HWL9gwRtAVBEGzk6OVmUuKCF/UZImgLgiDYQG1LH16e+kVLi1wlgrYgCIINXKnt\nIjrMb9GfI4K2IAjCAnX0DaNW2+eAuQjagiAIC3S2pI2EqAC7PEsEbUEQhAUwW2TGTPKi57KvEkFb\nEARhAc6WthIb6W+354mgLQiCsAAdfaN4ubvZ7XkiaAuCIMxT39AYSPYNoyJoC4IgzNO5sjYSo6ev\nE7IYRNAWBEGYp5FxM2r1Epxpt/UO2+MxgiAIdjM8ZkTGPjtGrrWgoC1JUrQkSW9JklQiSVKxJEmf\nn+q6U8VtHDhVjcksT/VpQRAEl3O2tJUEO6dGYOEzbRPwRUVRMoF1wGckSUp/90Vp8SGEhfrz/Ftl\n1LUNLPCRgiAIjjc0akQ30RrMnhYUtBVFaVcU5fLE74eBMiBiqmv1Og15aVFcqe3m2KUGRDcbQRBc\nldFkwWB2TAyzWU5bkqQ4IA84M9N1yTFBSBote49WiHSJIAguqbC6g5hw+x2ouZZNgrYkSV7AC8C/\nTcy4ZxTo50lSbAjPHyljcMRgiyEIgiDYTUffKN6e9jtQc60Fl6WSJEkL7AX+oijKy1Nd88sfPzn5\n+9XrN7Fm/Wa0WjW5qRHsL6hh24oYIoK8FjoUYY5MZgujBjM+Hjq71U0QhKXAYLJ9luBswQnOFZwE\nQJanT70sqBu7ZP1O/yPQoyjKF6e5Rilu7J/xPiU1HaxOCSE+YvFr0S53iqJwtqyNxs4hADRqNWaz\nBb1Oxe1rEtBp7b+wIgiuZGB4nMOXm0mNXbwONWazzIqEgCm7sS90pr0ReBQolCTp0sTHHlcU5eBc\nbpKZGMqFqg4klURcmO8ChyRMZ2jEwP6CGnbfsYq1pjFKrvlhajRZ+PuRcnaviyfQ18OBoxQE53ap\nsp24CMfks2GBQVtRlJPYKC+enhDK2fJ2JCBWBG6ba+oY5GRRC5nJ4fiZxlCAzp4hhscMRIX6odNq\nWJEeycGzddy9IQkfB+XrBMHZDY6ZCA62T8ODqTjuyVPITAzjdGk7GrWKyGBvRw9nyahv62cwPYuH\nB1o4WdTK99N28/eoVQx84y8AqFQSeelRfPT+DeSkRPBqQQ3vuzVd5LkFYQomB231u8rpao9kJYdx\norCZnoFRRw9lSWjtGuJsRSc5Ay28GJHHv3z7WX6bsIUBnQe+3noiQ31BgQslTXzue39n31vFRIUH\ncPxyk6OHLghOp3dwDI3ases+TjXTviorOZzXztTxwLZU9DqnHKJL6B8epyU1m/v6m/j3rHv4W8xa\nMJjY3lnOv1W9ibqoBIDB4XH+uu8sB46V8Nt/vI2iKKQmhNEzMCry24JwjcLqDqIdvGHC6WbaAJIk\nkZ0Swd6j5VgsrnMAR1EUimo6ef1MLQcKajhQUMPrZ2opre+2+wnQcaOZ/W/XkD3Qwqfz3s/fYtai\ntxj5/Ae28bvzfyRnoGXyWh8vPZ96eAuffWQrAL974RS9/cOcKGy265gFwdkNj5lx0zp2IumUQRtA\no1aRmhDGyycqnf7Iu6IonCxs4rJ/LJHp8USGBxAXFURcVBCR4QF09Bt44Wglxy41YpEX/4eQxSLz\n0rEKUhNC+HT+IxwKy8TXOMpfz/yenRvTkZEo8Qmnsr6TqvpOSqrbGRgaY9fmDD5y3zoAfvnscYZG\njbT1iAqNgnCV0QlOcTt17sHdTUtosB+Hz9dz2+p4Rw9nSt39I9QkZJEz1I6PeRyAtmsW8CRJIiTQ\ni5BAL4ZGDPzjSAUrU0NJjl68zs0vn6gkPiaIX/zlGMdD0vA3jvDsmd+SOtTB3pp2mk9XsTk3hky9\nFgBZUThV1ExJdT97tudwsbSJK+UtHDhWjJe7jnu3pCzaWAXBVYwZzDg+ZDvxTPsqfx93FLWGc6Wt\njh7KDSoaezhyqZm1ffWTAXsm3p5uZKdGUNkyyOtnahflHcRrBTWEBvvxx5fOcPx8NV6mcf7v3B+I\nGu3jTEA8t66M4fa1iXhMBGwAlSSxKSea21fHUVzVxuce3Ya7XsuZK/WY3v9+sSgsCMCVKsfVG7nW\ngk5EzuoBszgRORu1TT2kRfuS7ID6tVO5UtVBY9cwCdFBJMcG8lp4Fm8HJlHqE85gahZ6Nw2xkYHk\npESwbkU8Hnrdda8fHjFQ29TFezYl4+mum+Ypc3PsUiMyKg4cL2bvG5fRadX88cSvSRzu4uxbl7hj\nfeJNt/GNjBt55WQ1pTXt/OWVs+T0N/GJH3yBuzYm22SMguCqXi2oITYqyC7PWswTkXaTEB3IldoO\nfL30hPh7OnQsV6o68MrL4q6Rbn7ytwI+se3LtLlfs6Lc3gdAVUMXb54qx02r4bYNaTxwex7BAdYa\nK16ebmSmRLDvVA1r08MXfIS/xiuERIuJ3/50L3vfuIxapeLrn9hJ/MH/pMvNizs3JM3qPp56Hesz\nI9Bo1AQZhij0i+ZMcRO3r01Eq3H6N2aCsGicIZ8NLpAeuVZafAiHLzQwZjA7bAzl9T00dA4TZBjm\nU/mP8LM/H6XN3Y+UoXb+vewAL5z6Nb/89vt46kt7+MSDG8hMCsdgMnPgWDH/8u1neXb/OYwm6/jV\nKhXZKREU1/dyurjlJk+e3umSFrSymTdD0/nTP88iSfDYR7YTExlIWUExGYNtc7pfXLgfwb7ufLz2\nBAAFl2o5eaVx3uMTBFcnKwpGJ9nJ5jLpkavMFpmymjbee0saKpV9f+bUtw8wmpJO2PgAj679GGU+\nEXjodXzv9J/Z03oFFda/y5J3fb0Nrb08/+oFjp+vBiA82JdPP7yZvIzoyWu6+kYYHBzmjvVJaObQ\nKPToxXoMFomhe+/jK7kPAvC5R7eSlx6DpFjYlhczr69VVhTOBCXymfxH6NN58tlHtvCl92+Y170E\nwdW1dg1ysaaH+MjF20BwrZnSIy410wbrVsCk6BD2n6q263M7+0Y4W9pG5Fgf71/7ccp8Iogf7uJn\n33iQe1svTwbsqcRGBPDVj+/gqS/tISbcn7auAb71s/386Jk3GRgeAyDY35M7tqRzPiRpVo2QjSYL\ne98qQ6V1o7Smja/l3A/Axx/YQH5GDLLFNO+ADdbFyYjxAfa0XAbg7Ut1NHWKVnHC8lTW2Et0mHNU\nIXW5oA3g7q7F39eboxcb7PK8oREDh87Xkxofwmfz3k+FTziJw508d+a3hAR6Ue0ZzGXfKPadrOKV\nk5VU1ndSXteB0Wy57j5ZyRH89BsP8uF71+Gm1XD0bBWf/s5zHD9fjaIoaBSZdb11nClrZ//bVYyM\nmW4Yi6wonCttpTA0kffcmsXhgjJ+9uejyJKKxyreYNPKJEwGI7eujFvw1x070kPkU99FrVJxpaSB\n8qx1C76nILiicZNlTu+AF5PLLES+W4CfBy0dJi6Ut7EyLXzRnjNuNPPPt6vJTY3g/146zcngZAIN\nw/zx7DMoQGVdB1ua6kjy92TFNa8bGTNyqqiZ3mEDsRGBk10utBo1D9yex4a8BH7+l6MUVbbyg98d\n4sjpCp7yDCJxpJvkmGBMZgv1YbGYJDVtB48iIWFRFEbGzdy+cwUSCp/Pe5iDEznsJ4pf4Y62Ii6N\nG9i5xnZ72tdmRpCfGc25ogbeCk3jVllBpRKFpITlxbgITQ/my+Vy2u9W09RNapQfqTG23wpoNFn4\nx1vlZKdEUFTZyjd/ug+NbOGvZ35HiGEIi6QicbhzxntYLDIFxc209oySHBdyXZMBRVF4/WQZz+wt\nYHTciEa28EDzBW75w/8QGepHZsyNb8dMkopXIlbwVNouut288dDreOzD20mMCcZkMLLDhgH76hi/\n/8cT/O4fp4gc7eMnv/o3VqVP2btZEJYks0XmxRPVZCaG2u+ZS2HL33QSo4MoqetCo1aRGGm7je8W\ni8y58FQe6qnlQkU7P/3TWwB8oepNwsYH6LpUysrUsJveR61WsSk3BqPJwqFzdRjliebGkoQkSeza\nnMHanDj+/MpZDp0s4bmYNTz3xN9ISwhld9J2UoY6cLeY6Ne5c8U3mgPh2XTqfQBY2VvPJ3/x71hk\nBWSLzQM2WE90bsiJYf+fXqXFw58TlxpE0BaWlbq2fgL9nKdwmsvPtK8qr+sgJz7IJsfDTWaZF94q\n475dubjJFj715D85eKKUlLgQnvvtZ2lx9ydrYH5b9Dp6RzhR2Iy3tztRIdc3e9CnJ/G/CVv4Z8J6\nDKbptzUmDnfyqZpj3NtyiRePlpMc4Utu8uLNAswWma/nvY+90StZlRXD377/kNPk9wRhsb15ro7Q\nEH+7pgVnmmkvmaAN1sMsSRE+ZCWEzPseo+MmXjpWQUZyBHmJgRT6RnLPps+iklT89Bv3Mzhs4CEb\nNAgoqe2iuK6byDB//H3cASbTIecru7hQ0kjP49+kxisEk6TCw2IkY7CNdT21rO6rp1fnSY1nMAm1\nxQT5Lf5ho1fDs/lM/iO4adX89okH2ZwXt+jPFARncKCghjg7nYS8akmnR66VHBtMQ2sf3QMNbMuL\nnfPrW7qGOHaliey0SNQqFTIS/5HxHhQF9tyWjcEos3N1nE06umQmBJMeH8S50lbKqtvw8NCTOfE5\nd72WTSsTySx/7YbXtbv5cNkvGl/TGGt66xY8jtla3VtP2mAr5T4RvHWhTgRtYdlwpkVIWGJBGyA2\nwp/egTH+fqSMHavj8Pd2v+lrFEXh+OUmeoeN5KZGTn58X0QOl/xj8Pfx4IGdefT0Ddm0KYBKklib\nGcnaTGjsGODVgmqMJhljTQegYPF9ZywKEs2vHSUp0o+8RVh0vZlg4zBreuoo94ngQmkTZossUiTC\nkifLCiaL5eYX2tGSC9oAAb7u+Hm7cfhiEx46NVtXxODprr3hOllRKKrupLKpj/BQP1Li3skxm8wW\nfpy8A4AP7FlDY3s/ezbNrn7HfMSE+hITan2+oiiYzDJyZz0qlQqNWkIlSeQu2tNnQVFIOVgIP3+V\nitoOThc3syl3/od3BMEVtPUM4eV584mfPS3JoA2gUqlITwjFbJY5eLYetQp0GhVajYSsSJjMMmNG\nM8EB3mQk37jP+9DbZTR6BpI43MnmVUn09AzYrfWZJEnXbQ10FndtSuHpFwpoauuj6L4PsKnmmKOH\nJAiLqqKxl8h3bRhwtCUbtK/SaFSkXbMweXXhdaa8tNFk5rlXLwDwWOUh6pv/lbs2JCzuQF2Ar5ee\n3NQImtr6KAhM5JPioI2wxI0ZzU6XBnSu0djB1f3RMzlcUEHvwCjpg63sbi9Gr5VEg+EJ67KtC7wX\n/aIpre9y8GgEYXGZzM7X6nDZBe2bsVhk9r5hLZL06Zqj1HkGkb+Ix+Rdze4NycSMdDOidefQ6RpH\nD0cQFs3VtSVnI4L2u5y8WEN79yDhwT7sbitmQONORKCXo4flNAJ83Mnvs9bWvlTR4vRNlwVhvobH\nTKidLDUCImhfR1EUXnj9EgD378zDrFLhJjuu4YKzWt1bD0BFbfusysgKgisqr+8iwskWIUEE7etc\nKGmkrrkHfx8Ptq9Lodw7nLShdkcPy+nc2lWGt3GMzt5hkSIRlqyeQQNeHrbp32pLCw7akiQ9I0lS\nhyRJRbYYkCNdnWXfc1sOOq2Ghn1HUCnOl9NytFDDMHn91hTJhbJmB49GEBaHyWKxyelnW7PFTPsP\nwC4b3MehymrbKa5qw9NDx+7NmfT0jxAT6u3oYTknRSHxXz4IQEl1O8NjRgcPSBBszxl3joANgrai\nKCeAPhuMxaFeOVwIwO7NmXi46+joHiInyX71c13Nlvw4VJJEXXM3Ry/WO3o4gmBTFlnB4pwxW+S0\nAbr7hnn7Ui0qlcRd27IAcNOqxMGRGWzMjSU+KhBFgXPFIkUiLC09A6O4OenZDLuM6pc/fnLy96vX\nb2LN+s32eOysvXqsBFlW2LwykSB/L4ZGDIT4OVe9AWej1ajITgmnpqmboqpWZHE6UlhCyht6iAqz\n786RswUnOFdwErAWqpqOXYL2Zx573B6PmReD0czBE6UAvGd7DgCNrX3cuyXZkcNyCasyo3n5cBGV\n9Z0U1nSwIvnmnXwEwRUMj5sIDLTvTHvN+s2TE1qzWebpn/1gyuuWfXrk2LkqBkfGSYoJJi3BmsN2\n06rQapb9X81N7ViTSLC/FyNjRk5canD0cATBZpzxJORVttjy9zfgFJAiSVKTJEkfWfiw7ENRFPa9\nZd2pePf2bCRJwmi24OVxYxlX4UYh/p6TP+iulM+v/ZogOCOzs65CYpvdIw8rihKhKIqboijRiqL8\nwRYDs4eS6jbqmnvw83Zny0prrez6ll5WzaJhr2CVnWxt8lte30nP4JiDRyMIC2e2yFicuDzDss4B\nvHLEOsvevSUD7UT9asUi4+uld+SwXMrmvFi0GjUtHf28dd5+7c8EYbG0dA3h48QxYNkG7c7eIU5f\nrkOtUrF7i7U7o6Io6HXO13zAma1KjyAxxtr09GKZSJEIrq+6pY/wIB9HD2NayzZov3qsBFlR2LQy\ngQBfazfzvsExIoNFRb+50KhVZCZaS9cWV7dhsTjvAo4gzIbBZHG6xgfXct6RLSKD0cwbJ8sAuOuW\n7MmPt3UNkBkf7Khhuaz8dGteu7qxiwsVbQ4ejSAsjDPvHIFlGrSvbvNLjg0mLf6do+puWrVT/4R1\nVlvz4wgJ8GZs3MTpoiZHD0cQFsSZd47AMgzaiqKw/2gxYJ1lX63iJSsKeidspusKIoK8SY239uG8\nUtHq4NEIwvwZTRanb+yx7IJ2SXUbtU3d+HrrJ7f5AbR3DZEc7e/AkbkuSZLITrLmtSsbOmnrHnLw\niARhfurb+wn083T0MGa07IL2/ress+xdm97Z5gcwMDRKXLifo4bl8lZmRKLTqmlu7+dkoUiRCK6p\nqXOIABG0nUd33zCnLlur+V3d5neVVqtC5YQFz13F+uxoEqKti7iXxNY/wUU5+84RWGZB+2o1v415\nCQT5v7O1z2KRcRf57AVxd9OQkWhd1C2uacdosjh4RIIwd+YF7BwZN5qob+mhuqGL3oERG47qes5Z\nMHYRGE1mDp60VvO7dpsfQEvnIFmxgY4Y1pKSmxLJi4euUF3fwbnSFjbmxjh6SIIwJ6Z57Bwpr23n\n+dcucqW8+brJSrC/F5tWJrJrcwaRobZLvS6boH3ifDWDw+MkRgeRkXh9bZHhkXEiQ5z3BJSrWJcd\nSViQD+3dg5wrE0FbcC1mi8xcNo6YTBb+9x8nee146eTHIkN9cdNpae8apKtvmJfevMI/DxeyY2Ma\nj969Bn9fjwWPc1kEbUVR+OdEnZFrt/ldpRP5bJtIigwgOTaY9u5BCivFIRvBtXT1j6DXz67Cp8Fo\n5vtPH+RCSRMatYp7d+Ry9y3Zk6erZVmhqqGT106U8tbpSl4/WcbJizV8+uEtbFmVtKCGwcsip32x\ntInapm78fNzZujrpus+J/dm2o1JJZE1U/atq6KKmudfBIxKE2atq6iUi+ObvuC2yzPeffp0LJU34\neOn50dfu40P3rJsM2GD9XkiND+ULH7yFX377vazMjGZk1MgPf/8m//PHIxiM5nmPc1kE7b1vXAJg\nz/YcdNrr31z0DowSFSy6rttKbkrYxNa/Pi6UiYM2gusYGTfPqi/kPw5e4kJJIz5eep56bA9JMTOX\nvogK8+c/Pnsnn31kK246DUdOV/LVH71EZ8/8zjMs+aBdWd9JYUUr7nrtDdv8ADq6B0kTi5A2sy4r\niqRY6z/ii+J0pOBCZnN8vbiqlWf3nUOS4MsfvY2YiIBZ3VuSJHZtzuBHX72PsCAfahq7+cL3X+By\n+dybYi/5oP3C69ZZ9u7NmXh5uN3weZ1GjdrJ92W6El9PN9LirAu9ZbUdDI8aHDwiQZgd000qVJot\nFn717HFkReHB2/PJz4ie8zPiowL5n8fvZ2VmNIMj43z7Z/vZf7R4Tkfnl3S0aunop+ByLRqNij23\n5tzweUVRcNMu6b8Ch8hNm8hr13dytkQctBGcnywrmGfogA7w2vFSGtv6CA/24eE7V837Wd6eer79\nmTt4cFcesqzwm+dO8Mtnj2O2zO5sw5KOWC8euoyiwPa1qVPWExgdN+PvdePsW1iYFSmhRIT4MGYw\ncamy3dHDEYSbGhozolZNHw6HRsb5y76zAHzs/g3XlcCYD7VKxYfuWceXP3orWo2agydK+dZP9zM4\nPH7T1y7ZoN3dN8zh0xVIEty/c8WU1zS395KdGDrl54T5y4gLJjHGWvVPNEYQXEFFQw+RIb7Tfn7/\n0WJGRo3kpEawNjfOZs/dtiaFp760B38fD4oqW3nsqb00tM6862rJBu2/v3YRs1lmY17iDKeRJDzd\nRed1W9OoVWQnW6v+VdV3Ulzb5eARCcLM+obG8fLQTfm5cYOJfRPnPB66Y9WC9lhPJTU+lP95/H6S\nYqxnHL78gxc5X9ww7fVLMmi3dw/yxttlqCSJR+5ePe11bhpxoGax5KeG46HX0tI5wMVysYtEcG5m\nizxtMH79ZBmDI+OkxoeSnRKxKM8P8vfiqS/vYfOqJMbGTfzX0wenvXZJBu3nDpzHbJHZuiaZ6PCp\na2SbLTLuoonvolmfHU1agnUXyfnSZqcvLC8sb9PVHLHIMi8fvgLAg7fn2XyWfS29TstXP3Ybj9y9\nesbj9EsuaNc2d3P4dAVqlWrGFd6mtj4yE0PsOLLlxUOvJSfVOispq2mnurnPwSMShKkpijLtHu1L\npU109Q4TFuTDmpy4RR+LJEk8fOcqvvHJXdNes6SCtqIo/O4fp1AUuHNbJhEzLCyMGUyE+C28eIsw\nvbzUCLQaNXUtPZwRW/8EJzVutDDdUvnrEw3Ad25KR6WyXzp1dXbstJ9bUkH7bFEDhRUteHm43XQf\npU6jWtS3OgKsy4wiLcG6O6fhi18D8fctOKH69n6CppjA9Q6McKawHpVK4rb1qQ4Y2dSWTNAeN5h4\n+rmTADx85yq8PfXTXquIIlF2ERLgOblwc8Y/nja9KH8rOJ+WziGCpjjHcbigAllWWJMdd10xKEdb\nMkH72f3n6OwdIjE6iLu2Zc147cCwgdAAkRqxh1WZUWg1aor8oyn0iXT0cAThBkaLPGUpi2PnqgDY\nsTHN3kOa0YKDtiRJuyRJKpckqUqSpK/ZYlBzVVHXwcuHC1FJEp99dOtNa4m0dvaTETdzZS7BNtZl\nRpGdYt2z/XZQkjhoIzgd0xQtxhpbe6lv6cXTQzevGiOLaUFBW5IkNfALYBeQATwsSVK6LQY2W2Pj\nJn70zJvIssKeW3NIjr35jhC1SkIn0iN2ER7oRU6qdYZ9JiCe8+WiOYLgXKbaOXLiQjUAG/MS0Gqc\nK1YsdKa9BqhWFKVeURQT8BywZ+HDmh1FsRZbaesaJC4ykA/uWTur1+k0SyYr5BJWpkfhbRqj0iec\niyJoC07EbJF5d50oRVE4ft4atLesSnbAqGa20OgVCTRd8+fmiY/ZxStHCjl8ugKdVs1XPnbbrIq4\niEM19rcxJ4otXZUAnC6sY3BElGsVnENn3wju72oxVtfcQ0vHAL7e+kU7AbkQC+0ROatjbr/88ZOT\nv1+9fhNr1m9e4GPhbGE9v3+hAIAvfmg7sbMsRt7WNURGzOyuFWzD39udWzorOBCRy8XSJo5ebOQ9\nm51vBiMsPzUtfTe0GDt9pQ6A9SsS7Fpr/2zBCc4VWHfAyTOUiV1o0G4Brs3SR2OdbV/nM489vsDH\nXO9CSSPf/9/XkRWFh+5YyeZVSTd/0YShkXGiQpxrYWE5CNz3InG/eYP6ll5OFzVw54ZE0XxCcLiR\nMRMBAdeHwTOF9QCstcMJyGutWb95ckJrNss8/bMfTHndQr9rzgPJkiTFSZKkA94HvLLAe87o2Lkq\nvvfrg5jNMnduzZqxINRUtGpJdF53gPVZUaxIiwKgsKKZguK5t1kSBFt7d+OD7r5hahq7cdNpyE1z\nzi2qCwraiqKYgc8CrwOlwPOKopTZYmDvZjSZeWZvAT/8/ZuYzBbu3JrFvz60ac6nGsUipGO46TSs\ny41D76bjngkGAAAgAElEQVShpLqd9jvuEUWkBId7d6GosxOz7PyM6BuagDuLBY9KUZTXgNdmuuZb\nP9tPXnoU+RnRxEYEzCnQKorC2aIGntl7ipaOAVQqiU88uJG7tmXNOWCPjpvw85y6Zq6w+NZmRLAu\nN4GjZys5GpxKZEkL67OiHD0sYQrNXYMU1XRhMMkk3rUNWVLR8cYxNudG4+62NGrQy4pyQ19IR6VG\n5sIuP0oulTZxqbSJZ/YWEODrQVZyBEkxwSTFBpMQFYSX5/Utv0wmCw1tvVwsaeTouSoaW60V4qLC\n/Pi3D9xCemLYvMbR2N7PjnyRz3aUuHA/1q+I4/iZcl4PzyK7upM16REit+1EGtoHuFDRQf4dG9g9\n2kNJYz+Zg9Z66PLaJC77RaM+e5bcZNfv+DQ8ZkRzTYuxcYOJKxUtSBKsyopx4MhmZpeg/eWP3sql\n0mYulTXROzDK8fPVk/sgAdz1Wrw99Gg0KsbGTQwOj2OR3/kJGODryX07crlza9aCerMpsoyPp+gJ\n6UhZCcHc1lHKG2FZNH/jO5R8+AI5/SK/7Wgms8zrZ2pQVBpSE0KJG+1BRqK7b5gRtQ5PixEVCvn9\njbzaN465op2VqfObPDmLyoYewq+pBFpU2YrZLJMcG4yfj/OWubBL0N62JoVta1JQFIWG1l4q6jqo\naeymurGLhtZexsZNjI2bJq+XJIgK9SM1PpQN+Qnkp0cvuJEmgJvIZzvc2sxI2tpLOBSawfMxq9nZ\nUUpTxyDRoaKYlKM0dQzSmp7Lzv4mqhp6KKtp5+f5j1AQmMjg43+G279D6PgA7288ywcaThMb7kdd\ncy9l9d2kxwU5evjz1jc8TljoO9t/L5Zaj5w427H1d7Nrpl2SJOIiA4mLDJz8mKIojIwaGR4zYDZb\n0Ltp8fHS23wRQFEUsQjpBNQqFSnDndzVWsi+yBW8GJXP6sJmHrwlDa34/2N3jR4BDGv0rB1qo0fn\nyfd+fdC6TznMWnTNx1OPcXCQDr0v/5Oyg7/GrOXr9Z2kxIVwubKV+Ag/9DrnXLC7GaNZuW4n2aWJ\noJ3n5EHb4d8lkiTh5elGWJAPUWH+BPl7Lcqqbe/AKJFBXja/rzB3aT1NpH/7K2hlM/+MyCX5fXdw\nNixFFJOyI0VReOtCA7KkImOojbP+ceze9HlOX6nDXa/lM9VHOHnkKZ79749Q+voT/PXMb8nva6BT\n78PX//tlzhc3kJ4QysEztY7+UubNfM2/t86eIZo7+nHXaydrwDsrhwdte+noGSY1NvDmFwqLTqtR\nkRIdwIfqT6FIKh7Pvo+c/maeP1LG6DVpMmFxyLLCKyerQKMhbrSHA2FZfGDNx+jS+5CVHM6vvv0Q\nX648ROT4AAASsKGnlr+d/i3vazyL0WTh//3uEO1dg+jddFQ19Tr2C5ons/md7X5XUyO5qVFo1M5d\n5mLZBG21SnK6al3L2bb8WLZ3lpM43EmNVwg/T95OVkoE+07VcLqkRcy6F4nBaOb5I2WEhvgREuDF\nSxEr+FzewxjVGh5tKOC/vvgefL3d+fuRMvadrKKhpZtC30gqvELRKRaeLH6JTfmJjI2b+M9fv0ag\nvydXarsc/WXNmcFowXJNFY7LZdbF8PxM506NgJ1z2o7kphGnIJ2JWq1COnacRxp7+d4v9vPbhC34\nnKvmlrUpxCSEUugZQuurR/Dz1JEaF4y/tx61HXv0LUV9Q2McOFVLVnI4Wq2aU5dq+X85D6BIKr5Y\neYjPVR/htfYBkC3s2ZiE29VcdX8zXX0jnEnMImOojS986BaaO/qpb+nh2f3nuWNLFhWNPaTGuM47\n2fr2foInutXIskJhpbWHaW6qc56CvNayCNoms+yyiyVL2drMSBo6hvhm6QG+m3k3P/3zW/j7eJBp\nNpA30IQuOphxo4lTJW0YjGZUKgmNyppe8fXUkRUfjK/X9G3lhHc0dgxwqriVFemRSJLE+eIGfvC7\nQ1hUaj5XdZjPVx+h1Duc2GAPMhNurEkf7O9JUE8t/zxRRajJwhc+dAuPPbmXV44Ucuu6VIrrul0q\naDd3DhIU5AdAQ2svg8PjBPl7ztgM3Fksi/RIc8cAGfGu8w9quZAkiQ1ZkWzuruJD9acwm2We+PkB\nXojMn7xGr9OSGB1IRmIoafEhJMWGEBsZhJvenSOXmnnhaAUHTlVzrrRlyg4kAhTXdHKhqovslAgk\nSaKwooXvP/06ZovMx+pO8sWqN9l7rAKPitIpA/ZVkiSxZ3MyTW29xEYEcOe2LGTZWtNe76aluWvQ\njl/VwoybZDQTh7oKK6ypkZzUSJdo9r0sgvbouIFQf7FzxBlFh/owrHHjq+UHufe2XCyyzFdyH+ST\n+Y/S1jUw7evc9VpS4oJJTwzjzo3JZObEUxYUwxW/aA6cruHw+TpqW/qXfSA/cqGepp4xUmKt7fXK\na9v57q9exWiysGtzBt8oO0C1Vwg58UHEhfvd9H6SJHHP5hRKa9p59D2r8fHUU1LdRv/gKJcqOhb7\ny7GZa9dMrlS4TmoElkl6RKdWoRL5UKeV3dvA82+W8dG0SCJD/fjDH9/gjbBMDn37WXJTo1idHUNK\nfCghAd74+bijVt041/CwmMgZsH7zaSKDyIjxo03vy0uvnkatUuGmVeHhpmFtRgSe7ku//szIuJFX\nT9UQGR5A6EQKqbqxiyd+cYBxg5lta5L51MObObyzET8PDauiZ19jXqdVsz4jgpLGPu7ZkcufXj7D\nc69e4AN71mIwmt/JhTuxq4WiLBaZ4kprN6VsEbSdh5t2WbyhcFlqlYpd6xI4dL6RXZszePSjt/Oj\n1J38M24tl8ubuVz+zjF3lUoi0NeTIH8vggK8CPb3Yk3UKvL6G0ke7uTqj2YJiBgfIDMpfPK1RrOF\ng+fqUQFJkb5kJzr3ftz5Kq/v4VJ1J5lJYZN1XWoau/jmT/YxMmpk/Yp4vvih7YyOmbCYjKxKm/uO\nifgIP4pru9i1MYMX37hMcVUbBoOZk4XN3LoqzsZfkW1ZZAXLREnW6sYuRseNhAf7EhLg7eCRzc6S\nD9rDYwYCvMVilbML8HEnPzmYkoYuMg1D/LBwL/c//0sKLtdRWt1GXXMP3f3DDAyN09U3TFffMEyc\n63gx534AIsb62d1WxNrOATKneIZOoyYt3hqoO3qHeeGtCrLig0iLWxrrHUMjBg5fqMfdXU9O6jtt\nsmqbuvnmT/cxPGpgbU4cX/34DhQFGlq6ee/2+ffhvm11HPtO1fKe7dk8u/88+44W8b7d+Td/oYO1\ndg3i4+UOQHGVtRiWM7YVm86SD9pNbQPsXhvn6GEIs5AcHYDZIlPuHUraUAfennp2bkxn58Z3AovJ\nZKGnf4SuvmG6+4bp7Bmi5+nfcd4/jlZ3P36fsJlnnniWXXnv56sVr0/7rNAAL0IDvKhr66Okvotd\naxPxdHfNkqN9Q2O8XdhM/B1b2TPYRmlj/+Tnapu6+eZP9jE0YmBNTixf/5edaDVqLpc3c9+WlAUt\nvLm7aQn21bN1TQp/P3iR88UN3LEli7aeYcIDnXcNqba1n9BA66y6tLodgMzk8Jle4lSWfNAGBU+9\na34zLkfpcUE0VFbwbHELGWbLDQeitFo1YcE+hF3T1y/zX29FRuKyXzTPR6/ixfi1vBaezZuh6dz3\nzzM8dOeqaQ9WxYT7Y7bI7CuoIT7Mm7UZrpHXHBwxcKmynf4RIxYFkmKCyBy8vtP9xdImnvzf1xkb\nN7E6O5bHP3E7Wo2a8rpOttqoLvaW3Gj2Hq9iU34iR89WUVLdSoi/O+HrZ98C0N5GDWaCtGpkWaG0\nxvp3dm0azdkt+aCtE4dqXE5smC+hAZ4cvlDPyJiFhOjAGzpmv9vVsqH5/Y3c+defs/+Rz/NCVD7P\nv3aRM4UNfO0TO4gO85/ytRq1iuzkcLp6h9n7Vjm71ifgqbfPYqUsK/QOjdPeM8TA0Dgmi4xFYTI3\nryggKzIWi4IFsFgUTGYZlVpFTLjf5F7jaymKwsuHC/m/F09jkWU2r0zksQ/filarprljgPhQb6JC\nbFNVUa1WEeClZceGdI6ereLQqTLy02NQFMVpt89dbXzQ1N7H0IiBQD/PyZm3K1jSQVtRFNzE0XWX\npNdpuHN9EuNGM+fL2+jo6p/o5yehKNbyvVdZfCNRkEgdakcvmwn08+QHRXt5sPk8n7/7q9S39PCl\np17kyx+9lTUzdCQJDvAiwNeDAwV1pEb7kZu0OAuVfUNjnC1tZcRgwWyR0bvp2Lo1kzSLkdK6HnLi\nrv/hUtzQN+sA2Kz347u/eo1zRQ0A3L9zBR+6Zx0qlUT/4BjIFvJtXAd7y4pYXjxRRXxUIHXNPdS3\n9lDd1Euykx62MU/sHCmtts6yM5LCnfYHzFSWdNDu6BkiIcL5TzgJ09PrNGzKucnuhv5mzBaZgqJm\nWntGSDKaAVjd18DPv/lefvKnI7x9sZb//PVrPHL3Gt67K3/aLaBqtYrMpDBaOwd4+UQld6xLRGeD\nWu5grVt9vqIdRVKRFB14XcceX/O49flTjGs2AaV/cJQfpezgmbhNjBU14Omh4wsf3M76FfEAmMwW\nWtp7eXABC4/T0WnVeOs1bF+Xwu9fKOByWRM5yWFOGbRlWZls5ltS7XqpEVjiQbunf5StOa6RoxQW\nRqNWsXlFDBaLzKuna4h08yHMMIi7XsvXP7GTfxy8yJ9fOctfXjlLfXMPX/zw9hn3E0eE+GIye/PS\n8Upyk0JIW0CFyOFRI+VRKQQahkkuq54MwkaTmeaOfto6B7gcs4YRtRs9RwopjF6Nu8WEt3mcAOMI\n3Z0D+Hjq8XDXTf6wGRs30d49SFVDJ+eKGjhf3IgpaTsAm1cl8bH71xM0caDMIssUV7by3u3pizaj\nXJcZQc/gOCqVxIWSJnZvyXTKFEl77zDeHtbuVSJoOyGNWjV5VFVYHtRqFXdvTOb4mWKKLBCBdab6\n3t0riY8K4oe/f5OTF2vo6hvmW5/aNWNbKa1GRVZKBPXtA5Q39rJjVeycDubIssyxy010D47zQF8D\nEnC4Z4iTF2o4V9RARX0H5qsnNrPutf73729D9n3X3+jbzwLWPepajXUBzWS23PC87R1lfKr2GO6/\nuXLdGArLW7l3a4rN3jFMJdDXA28PHXnp0VwoaaSstoOGjhjiwm5+ytKeapr7CA30obN3iK7eYTzd\ndcRGzP5gkTNY0kHbTScC9nK1JTea45cb6egdJjTAOuNcnR3LD796L//xiwNU1HXwpR+8yH985k6i\nw6deoLwqOszXOoM/U4enm4atK6JvGryvVHdgXLOWdQOt1NR1cjoggd/Fb+LIN/86eY0kQWSoL5Gh\nfiTu/zseZiMdH/0UHn/+A+NqLYMad/p0HnQkZTI4PM7ouBHDROrHTashwM+TxOggMpPDWb8inq25\n1s72JRP3t8gyhRWt7NmUZJeF1bhQH1ZnxXKhpJHLZU2U1ic4XdAeNZjxD1BN5rPTE8Nc7rT0kg3a\nBqMZL7cl++UJs7BlRQz7T1UzotfiOfGWODYigB9/7X6++6tXqWro4ss/fJFvfHIXOTc5wqxWq8hI\nDMNslnntbD1alYSPh5bU2ECC/TxRFIWWrkGqmnoZGjez+fbVBBuGOBSawQ+efIGadZ8AQKtRsyEv\nno35ieSkROI10Wg686n3A1Dy0LNkfvXu655d8px137XJbLHWzJCsQXum1MPImIGq+k7u25KKh522\nvOamhFLa2IPeTUNlfSctndPXjnEUk0VGkiSXTY3AEi4Y1djWR/Yirf4LruOO9YnUNnVjkd8pEOTv\n68GTX9rDutx4RkaNfOun+3nj7bJZ3U+jUZGeEEpSXAh+/j6YYuMoD4zm5ber0SfGsWNdEmkJYTS7\n+/PedZ/kkys/QE1jN4GGYb5Q+SZ/+P6jfOVjO9iQlzAZsGdLq1Gjd9Oi12lnDNiNbX10dg/y0G0Z\ndgvYACpJwt9Tz6qsWACKq9rpGRyz2/Nn42rNEVcO2kt2KmoyW/D3mts3hbD0qCSJuzcmsa+gluxr\nTr3pdVoe/+RO/vDiaV5+8wo/+/NRymra+eRDm9DrZhfodFo1saPWVluqxDDCxwep9Qzi+08f5NSG\nTwEQYBjmwQ/u4osfuQU32UyxtztDIwaa2vtRqyR0Ggm1SqLlTC2SpKC09dJYUI3ZohB75zYAurqH\nCAn0uq4J7VRa9b60uvuRHOFLUtTMKZ/FkpsSwvnyEE5eqKGitp2L5W3sWJPgkLG8mywrmC0yQyPj\nNLb2odWoSY6dvhSts1qyQVunUTndyrXgGJ7uOlYmh1De3EtC1DuLTmqVio8/sIGYcH9+89wJDp0q\np7iqjX/74DaykudWi6K1c4A/Zd/H3sh8LJfq0FuMfLz2JP9Sd5yGZ77MoEZPi7sfrc09hAZ6sGfj\nLLYSTmxlrGvtp769F4NJxmyZeIsPyChIE1+HTi2hq6ydU7W+xRAZ5E1SdBA6rZqK+k5auoccOp5r\ndfSO4OXhRmmN9eh6SlwI2kVcnF0sSzJoi0M1wrslRwfQ0jVEb/8oAX7X7xjZuTGd5NgQfvTMmzS0\n9vL1//4nG/MTePjOVcRFTr/VT1EUin0i+N+ELRx44m/I0atRKTI7N6bz3e9+gGDDMCU+4dQ3daGt\nrGdFbCAr5jhujVpFcnQAydcEY0VRkBXrf1Uq6aYzcHsL9NGTlx7NmcJ6LpY288jOLKfoz1rT3Eto\nkA9vnioHXDM1Aks0aPcOjBIR7OnoYQhOZmteDP84Uo63txvad3Xcjo8K5CePP8DfD15g7xuXefti\nLW9frCU1PoRVWbEkxQTj46XHYpHp7h+hsr6TC8WNNG36HAAalcQD9Wf515pjjPymmr6nPGhx9ye1\ntZpsGxeikiQJtQTvHHZ3LusyIzl+sZ4zhfXUNHZxpaqDVemOr6I3bDDhp1FN5rMzkmx7MtRe5h20\nJUl6EPgPIA1YrSjKRVsNaqHau4dYv9l5C9YIjiFJEndvSmbv8Qry0qJu+LxWq+aRu9dw+6YM/nHw\nIkfOVFJR10lFXee09/Q3jnBPy2W2/O03bFvxNUbUOvZVtnJrQw1p/stz4uDrpSc9MRSVJFFR10lt\na79TBG2TWcFgMlPT2IUkWbf7uaKFzLSLgHuBp200FpvRqCWneDsmOB93Nw0bMiMpqusmKSZoymuC\n/L341MNb+Mh96zlf3EhpTRtNbdbiQmq1hJ+3BwnRgWQkhvO+3bloFZmSAC8OXWjEYjTy3rUJTpey\nsLcQX3dSE0Ipq2mnqLqdBxfxNOZsmSwylXWdmC0yCdFBeLq75kaFeQdtRVHKYXZ1EezNTbNkdzIK\nNpAQ4Udr9xCdvcOEBExf91nvpmXTykQ2rUyc9hqtYt1KWNfcQ4ivnjX501+7nKxMCycpJpiymnYq\n6zpp6x0mwoGV9Kzdama31U9RFKoau1GhoNOokBWF4TEzcZEBeHk4PtAvuZz2yJiRQB/RqUaY2aac\naPYeLcfXW4+bdmHfBqXeYSSEeZMRH2yj0bm+QF8P0hNC2fdWEeW1bVyp6nBo0O7oHcbTw+2dyn7T\npEYGhsZobu/jlrwYgq9Jb1ksMieLmmnpHCQ1zrH/n2f81ypJ0iFgqq/u3xVF2Tfbh/zyx09O/n71\n+k2sWb951gOcq4a2Pu5YG79o9xeWjvdsSuH5w6XkpEVO2Sz4ZhRF4ZJfNNGjvQSJgH2DhEh/QgO9\n6egZoqK+m93rHLfOVNXUS4i/J+W11o7xU3WqGRoZp6t7kAdvSbshg6BWq9i6IoaWriGOFzaTnWz7\ncq5nC05wruAkYN1TPp0Zg7aiKDtsMZjPPPa4LW4zKxKITjXCrGg1Ku7flsreY5XkpEbMKXCbzBZK\nqtrY1VBBoI/7Io7SdWUnBJOWEEZHzxBFVW0YTZZFLVo1k1GDmf6RfsYMJsKDfQjwvX6R2Giy0NjW\ny4PbbgzY14oM9mb3mjhePVNHTkqETQP3mvWbJye0ZrPM0z/7wZTX2Sr56zSJbdGpRpgLdzct92xK\npqi8BYPJPKvX9A2OUVnXwQO3pImAPYOIYB9S4qwnDlveOEp9gON2kBjNMiXV1ia+GYk3zrLLatu5\nZ9Psemb6ebuzc1UspTUdNh/nbMw7aEuSdK8kSU3AOuCAJEmv2W5Y82M0W/CeQ+lMQQDw8tDx3lsz\naGjuprVzcNrrFEWhvK4Ts9HIg7ek4eaCp+nsLSspFL2bhgqfcJrcHXO0HsBsViiZpolvXUsva9PD\n5/QuIMjPk9zEIOpaem06ztmYd9BWFOUlRVGiFUVxVxQlTFGU3bYc2Hw0tPSSk+R6tQQEx9NqVNy7\nJZVwfzfKatqpberBYLJgkRV6+kcormqjtrGLrTkRbMuLccpdU84oIzaItHjrstgF/1iHjGHUYEaR\n3mkvdu3OEbNZRjGbSYyc+w+U1JhAvNzU9PSP2myss7Gkdo9YzDIB3mLniDB/WQkhZCWE0Nk7QlVT\nLwazhTB/D9alJTssH+vKEqP8SY4N4XJ5Mxf9YujuHyXIb/rGE4uhrK4LlSTRPzSGn7c7ESHvtCCs\nqO/k7g3zL2i1dUU0/zhSjr+PHtU8FrPnY0ltaNZpRZEowTZCAjzZmBvN9pVxZCSEiIA9T5IkkZtm\nndle8o/hbFmr3cfQMzhGfUs3YD26fjVGjIwZCPLV4+42/40LkiRx+9p4ymY4NWtrSyZoWywy7qJT\njSA4nbSYQDIGWjCqtRRVtdv9+UazTNlEZb9rUyO1zb1szb1J0+hZ8PXSExPkSU//yILvNRtLJso1\ndw6SHjf1sWRBEBwnNzmMzAHrDLuwstXafceOrDtHruazrTtYxo1mAr10qG3UQ3ZdVhQdXQMoyvT7\nq21lyQTtkdFxwoMcd+JKEISpadQq8vsbAWg6eY5mH/sdRDKZLfQOjtHWNYi7m5b4KGup3ZrGLrbl\nx9n0WZtzo6hu7LHpPaeyZIK2Xqte9kV6BMFZpQ+24WscpckzkFJv+9Wxrmnpp7vH2oghIykMtVqF\nxSLjpdeitXGNorBAb9y1MDZusul9321JBG2LRcZDNPEVBKcV31pDfE4KABf8Yuz23KbOAepbrbPf\n7BRr8+aapm425txYmtcWblsVT01j16Lc+6olEbQb2vrIShC1HwTBWfl4upGWYG20fdkvhuFRo12e\nazDJFFdZ8+nZKdZ8tgT4zrGp8myp1Sqy4gNnPKS1UEsiaBuMZkIDlmfBeUFwFbkTQbPQL5LTxc12\neWZn3ygtHQO4u2lJjAmivXuI1OjFPZmZmRDCwNDooi1KLomgrRf1swXB6a3JiCBtsA2TWsv58pZF\nf57RZKG60bp/OiMpDI1aTf/gCGmx0/f9tJVb8qKprF+cNInLRzuD0Yy3h6jqJwjOLi7cj5x+6wy7\nqLJt0bfHldV309Y1AEBWcgQWi4y3u9YuB/ACfT3wdFMzZrD9oqTLB+26ll5WprlmV2VBWE5UKomV\n/Q0AVNZ3Utvat6jPa+8doWpitpudEkFDay+r7Bgrbl0VR1W97U9KunzQlhQFL1HZTxBcQtSrL+Hl\n4UZ33zBvX2lc1Gd19o3Q0tGP3k1DUmwwsiwTYMdSuhq1ivSYQNq7h2x6X5cO2oqioNeJmhCC4CpW\npUWQFGPd6XVpEfPaiqJQMVEPJCMxHLNFwd/T/pO7FSmh9PYP2zQV5NJBu6t3mPgIP0cPQxCEWdJp\n1WQmWUu1ltR0YDDOrvHEXHX2jV6zPzuCxpY+Vmc4pgnD5pxoKmy4KOnSQbu7f4TU6ABHD0MQhDnI\nS7MebKlt6uJUYdOiPKMjIZXW46cBa9DWalhQNb+FCAv0xE0j2ewHlEsHbTeNCpVKHF0XBFeybWUs\nMeH+WCwKBcWLE7R7tZ7UeYWgd9MQHuJHiK9j28LtXB1PeZ1t2pO5bNAeGzfh64AclSAICxPgrScl\n3no6sqiqDXkRtv4V+1pTIRmJ4bR1DZKfGmbzZ8yFTqsmNcqfzt7hBd/LZYN2XXMPazMjHT0MQRDm\nSJIkcib6NFbXd1Je323T+w+OGCj0taZgclIjcNeq0Ggcv2FhZVo4XT2DC16UdNmgrVFL6HWiSJQg\nuKId65Lw8tDR3T9i87z2xYo2rvhZg3ZybCiRQc5T4uKWvJgFn5R0yaBtMsv4eIjUiCC4quSoANIS\nJhr+lto2aA/ctYdOvS/B44PoNGpykkJtev+FCPb3JMBLy8DQ2Lzv4ZJBu66lh5Vpjs1RCYIwf2qV\nxIpUa3qztKaDVhsdQJEVhUsTpV+3dlfi6a5xur6xW/JiaVjAaVCXDNqKLOPnJbquC4Ir25QXj06r\nprGtl2MXG2xyz5rmPkomFiHz+hpJCPe9ySvsTyVJbF0RTeU8j7i7XNA2mi34itSIILi8dVmRJMWG\nAHC+xDYpkpK6Tgp9o5EUmejRPlJjnbNvbESQFz7uGgaHxuf8WpcL2jWN3awTu0YEweW5adWTNbaL\nq9ttkiK5XNGGUa1hRX8z3ua5B0R72r4yjvrWHmR5bo2OXS5oa1Tg6S5KsQrCUrA2Jxq1SqK6oZMT\nlxc22+4bGqO4ytp1fV1PDeHjA7YY4qKRJIk71iVQXN0+p9e5VNAeGTMS4C1y2YKwVNyyMp6U+FBk\nReF0UcOC9jAXFDZTVmMN2lkDLYSOL17LL1vx9dKTkxBMQ2vvrF8z76AtSdIPJUkqkyTpiiRJL0qS\ntOgZ/7qWHtZnLU5DTkEQ7M/Hw40VE7VISqtaKVvAQZsr1e2MjptIHuogwsln2dfKjA9CK8HA8Oy2\nAS5kpv0GkKkoSi5QCTy+gHvdlKIoeGjVNm97LwiCY63KjEKrUVHV2MW5srZ53aOxY4DyWmttjzW9\ntcSP2PaU5WK7bXUcja29mMyWm1477wioKMohRVGuZtDPAIs6BW5s62NFivNskhcEwTY250STkRiO\nosDZooZ5VcO7XNlBYYW1Pnfwt/4dX9P8D684giRJ3LclleKqm7dhs9W09aPAqza615RGx4xEh/gs\n5vRTSMEAAAkcSURBVCMEQXCAYH9PVmZGA3CprIkTV+a2IGkyy1wsb6V/aIyoUD/S4pxzm9/N6LRq\n7liXQFHlzO82ZgzakiQdkiSpaIpfd19zzTcAo6Ioz9pm6DfqGxwjJsR7sW4vCIKDrc6Mwt/Hnbau\nQVofeBiLZfbb4E4UNlExUfY0Jy3KpQvJBfi4sy0/esYyrjNWXFIUZcdMn5ck6cPAHcCtM133yx8/\nOfn71es3sWb95pkuv0FrRx8P3JI2p9cIguA6NuXGkJ8Zw+GCCs4GxJNQ0sKmnOibvk5WFJo6Bjlf\n0oAkwarMaJc9LX306FGOHj0KWA8RTmfeZfIkSdoFfAXYqijKjLvYP/PY/NcoR8YMhAR4onKy+gGC\nINiOl7uOdbnxHC6o4FBoBvmNPWzIjrrp9/3Z0lZqm7oxm2WyUyOIC3fd9oPbtm1j27Ztk39+8r++\nN+V1C8lp/xzwAg5JknRJkqRfLeBe06pr7mVzttjmJwhL3YqkUDZ3VmBUaymqbOXIhfoZrx8zmKhr\nG+CNt8sAyEmJZEOW66ZGZmshu0eSFUWJVRQlb+LXp205MID+wTGigjxRq8U2P0FY6vLTwrijvQSA\nN94up7NvjLYZOr28fqaW3v4R2rsHCQ/2ZWVmNG7LoMa+U0fD5o4+NohZtvD/27v32KzuOo7j708L\npbR0MC4C45IC2woWRHAhQEASYZFNhBkT3NQ4NRpv0Wm8oonGfzTRmM3M6B+4LZsTlgWYzmQxMBDR\nuZmxjcu4jNkNuXQtUMbdrn3o1z+eZ5GNtrDN9ncO/bySJuecnvR80j79PL9ze471CWUSdWeamH2s\ngXOtbexpeIVNW//Nudb2i9Z9YuchqqoH8sjG7QAsnFPHlPF94yHfmS3txqOnmFo7LHOfhWtmPafu\ndBPLGrcBsHb9NkaNqGHdln20nCxed10438GGp1/mbFsH+15upuHAMYbUDGTyhJHUjR+WMnqvyWRp\ntxfOc+b0Oeonvit1FDPrRVWFNgauXsXcGRP5z2vt/Gb13/n4B6dybOwE1v71BXYOq2XunMkMGljB\nyoefAOC2JTdwzfBBfWaAl8nS3tvQzE1zJqWOYWYJzK4fzeL311NdVcEzuw5wz4R51J1pZvLEUcw4\neZCKjgJ3P7iZU2dbmT55DONHD2X+9EtfHnilyFxpNxxs4X11I/3QXrM+auTVg6ip7M+Xbi3ez/GT\nKR/izusWca61jSMDavjizE/y1Pb9DKzszxeWz2NIdQX9+tDFCnqnj3O/5AakeP7Aictat+noaaoH\nyA85MOvjWk6eY/O2wxx4pYW7H9hEqIx+/cooFIp3StZUD+DHX13Ca23nuWX+tVfkIE8SEXHRMZ/M\nvD01Hj1FOQUXtpkxbHAVIwYPYNa0WlY+8ztmtbxEodBB5fk2Fjbv5mff+ghXD65m4uirrsjC7k4m\nRtovHWphWE0Fc/1Z2WZWEhE8tHEPy2+cSjnBU3uaeG/9WCo7CuzYf5y9DU0s/8CU1DF7TCZH2qfP\ntrLjhUam1Q51YZvZG0jiowvqeHLYRM4jaqorqewoUFAZO/ceZsnca1NHTKLXR9odERxofJVzrW2M\nHlrFnKljKC/LzFEaM8uYtvbzPL51P2daC5SXiYpysXj2pCv+gShdjbR7pbQf2bKP8jJRVib6l4v6\niSMYM9wftWpmly8i+sy12JC4tHt6G2ZmV5pMHtM2M7O3xqVtZpYjLm0zsxxxaZuZ5YhL28wsR1za\nZmY54tI2M8sRl7aZWY64tM3McsSlbWaWIy5tM7MccWmbmeWIS9vMLEdc2mZmOdLrpb158+be3uTb\n5qw9Iy9Z85ITnLWnZDGrS7sbztoz8pI1LznBWXtKFrP68IiZWY64tM3McqRXHjfWoxswM7tCJXlG\npJmZ/f/48IiZWY64tM3MciRJaUv6uaQ9krZLWidpcIoc3ZG0WNJeSS9K+m7qPF2RNE7SXyTtkvS8\npK+lztQdSeWSnpP0p9RZuiNpiKQ1pdfpbkmzU2fqiqQVpb//TkmrJA1InQlA0r2SmiXtvGDZUEkb\nJO2TtF7SkJQZX9dF1kz2VKqR9nqgPiKmA/uAFYlydEpSOfArYDHwbuA2SVPSpupSO/CNiKgHZgNf\nyXBWgDuA3UDWT6b8EngsIqYA7wH2JM7TKUm1wOeBmRExDSgHbk2Z6QL3UfwfutD3gA0RcT2wsTSf\nBZ1lzWRPJSntiNgQER2l2X8CY1Pk6MYs4F8RsT8i2oGHgGWJM3UqIpoiYltp+gzFcrkmbarOSRoL\n3Az8FrjorHhWlEZU8yPiXoCIKETEycSxunKK4ht3laR+QBVwOG2kooj4G/DqmxYvBe4vTd8P3NKr\nobrQWdas9lQWjml/FngsdYg3GQMcvGD+UGlZppVGXTMovsCy6E7g20DHpVZMbAJwVNJ9kp6VtFJS\nVepQnYmI48AvgANAI3AiIh5Pm6pbIyOiuTTdDIxMGeYtyExP9Vhpl45b7ezk68MXrPMDoC0iVvVU\njrcp67vuF5E0CFgD3FEacWeKpCXAkYh4jgyPskv6ATOBX0fETOAs2dmNfwNJk4CvA7UU97AGSfpE\n0lCXKYrXG2f+fy1rPdWvp35wRNzY3fclfZrirvLCnsrwDhwGxl0wP47iaDuTJPUH1gIPRsQfUufp\nwlxgqaSbgUrgKkkPRMSnEufqzCHgUEQ8XZpfQ0ZLG7gB+EdEtABIWkfxd/37pKm61ixpVEQ0SRoN\nHEkdqDtZ7KlUV48spribvCwiWlNkuIStwHWSaiVVAB8DHk2cqVOSBNwD7I6Iu1Ln6UpEfD8ixkXE\nBIonyjZltLCJiCbgoKTrS4sWAbsSRurOXmC2pIGl18Iiiid6s+pR4PbS9O1AVgcZme2pJHdESnoR\nqACOlxY9GRFf7vUg3ZB0E3AXxbPx90TETxNH6pSkecAWYAf/29VcERF/Tpeqe5IWAN+MiKWps3RF\n0nSKJ0wrgAbgM1k9GSnpOxQLsAN4Fvhc6QR6UpJWAwuA4RSPX/8Q+CPwMDAe2A8sj4gTqTK+rpOs\nP6J4tUjmesq3sZuZ5UgWrh4xM7PL5NI2M8sRl7aZWY64tM3McsSlbWaWIy5tM7MccWmbmeWIS9vM\nLEf+C6QEhxTLLrqmAAAAAElFTkSuQmCC\n",
       "text": [
        "<matplotlib.figure.Figure at 0x111579910>"
       ]
      }
     ],
     "prompt_number": 22
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Exercise 1\n",
      "\n",
      "One thought that occurs is as follows. Do we need all the data to create this posterior estimate? Are any of the data points redundant? What happens to the model if you remove some data?\n",
      "\n",
      "*Hint:* \n",
      "```python\n",
      "X2 = np.delete(X,range(8),0)\n",
      "y2 = np.delete(y,range(8),0)\n",
      "```"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# Exercise 2 answer here"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 25
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Building the Low Rank Approximation\n",
      "\n",
      "Now we\u2019ll consider a GP that uses a low rank approximation to fit the data."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from IPython.display import display\n",
      "Z = np.random.rand(3,1)*12\n",
      "m = GPy.models.SparseGPRegression(X,y,Z=Z)\n",
      "display(m)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "creating /var/folders/22/6ls22g994bdfdpwx4f9gcmsw0000gn/T/scipy-neil-tKO4Yo/python27_intermediate/compiler_08f57918657bc3b6af9ef49d104a532a\n"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "clang: warning: argument unused during compilation: '-fopenmp'\n",
        "In file included from /Users/neil/.cache/scipy/python27_compiled/sc_1790bf65208b11355ffcfd4b65a5f1090.cpp:11:\n",
        "In file included from /Users/neil/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/scipy/weave/blitz/blitz/array.h:26:\n",
        "In file included from /Users/neil/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/scipy/weave/blitz/blitz/array-impl.h:37:\n",
        "/Users/neil/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/scipy/weave/blitz/blitz/range.h:120:34: warning: '&&' within '||' [-Wlogical-op-parentheses]\n",
        "        return ((first_ < last_) && (stride_ == 1) || (first_ == last_));\n",
        "                ~~~~~~~~~~~~~~~~~^~~~~~~~~~~~~~~~~ ~~\n",
        "/Users/neil/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/scipy/weave/blitz/blitz/range.h:120:34: note: place parentheses around the '&&' expression to silence this warning\n",
        "        return ((first_ < last_) && (stride_ == 1) || (first_ == last_));\n",
        "                                 ^\n",
        "                (                                 )\n",
        "In file included from /Users/neil/.cache/scipy/python27_compiled/sc_1790bf65208b11355ffcfd4b65a5f1090.cpp:23:\n",
        "In file included from /Users/neil/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/numpy/core/include/numpy/arrayobject.h:4:\n",
        "In file included from /Users/neil/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/numpy/core/include/numpy/ndarrayobject.h:17:\n",
        "In file included from /Users/neil/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/numpy/core/include/numpy/ndarraytypes.h:1761:\n",
        "/Users/neil/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/numpy/core/include/numpy/npy_1_7_deprecated_api.h:15:2: warning: \"Using deprecated NumPy API, disable it by \"          \"#defining NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION\" [-W#warnings]\n",
        "#warning \"Using deprecated NumPy API, disable it by \" \\\n",
        " ^\n",
        "/Users/neil/.cache/scipy/python27_compiled/sc_1790bf65208b11355ffcfd4b65a5f1090.cpp:24:10: fatal error: 'omp.h' file not found\n",
        "#include <omp.h>\n",
        "         ^\n",
        "2 warnings and 1 error generated.\n",
        "\n",
        " Weave compilation failed. Falling back to (slower) numpy implementation\n",
        "\n"
       ]
      },
      {
       "html": [
        "<style type=\"text/css\">\n",
        ".pd{\n",
        "    font-family:\"Courier New\", Courier, monospace !important;\n",
        "}\n",
        "</style>\n",
        "\n",
        "<p class=pd>\n",
        "<b>Model</b>: sparse gp mpi<br>\n",
        "<b>Log-likelihood</b>: -89.7838229455<br>\n",
        "<b>Number of Parameters</b>: 6<br>\n",
        "</p>\n",
        "<style type=\"text/css\">\n",
        ".tg  {border-collapse:collapse;border-spacing:0;border-color:#999;}\n",
        ".tg td{font-family:Arial, sans-serif;font-size:14px;padding:2px 3px;border-style:solid;border-width:1px;overflow:hidden;word-break:normal;border-color:#999;color:#444;background-color:#F7FDFA;}\n",
        ".tg th{font-family:Arial, sans-serif;font-size:14px;font-weight:normal;padding:2px 3px;border-style:solid;border-width:1px;overflow:hidden;word-break:normal;border-color:#999;color:#fff;background-color:#26ADE4;}\n",
        ".tg .tg-left{font-family:\"Courier New\", Courier, monospace !important;;text-align:left}\n",
        ".tg .tg-right{font-family:\"Courier New\", Courier, monospace !important;;text-align:right}\n",
        "</style>\n",
        "<table class=\"tg\">\n",
        "<tr>\n",
        "  <th><b>sparse_gp_mpi.</b></th>\n",
        "  <th><b>Value</b></th>\n",
        "  <th><b>Constraint</b></th>\n",
        "  <th><b>Prior</b></th>\n",
        "  <th><b>Tied to</b></th>\n",
        "</tr>\n",
        "<tr><td class=tg-left>inducing inputs        </td><td class=tg-right>(3, 1)</td><td class=tg-left>          </td><td class=tg-left>     </td><td class=tg-left>       </td></tr>\n",
        "<tr><td class=tg-left>rbf.variance           </td><td class=tg-right>   1.0</td><td class=tg-left>   +ve    </td><td class=tg-left>     </td><td class=tg-left>       </td></tr>\n",
        "<tr><td class=tg-left>rbf.lengthscale        </td><td class=tg-right>   1.0</td><td class=tg-left>   +ve    </td><td class=tg-left>     </td><td class=tg-left>       </td></tr>\n",
        "<tr><td class=tg-left>Gaussian_noise.variance</td><td class=tg-right>   1.0</td><td class=tg-left>   +ve    </td><td class=tg-left>     </td><td class=tg-left>       </td></tr>\n",
        "</table>"
       ],
       "metadata": {},
       "output_type": "display_data",
       "text": [
        "<GPy.models.sparse_gp_regression.SparseGPRegression at 0x111578b50>"
       ]
      }
     ],
     "prompt_number": 23
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "In GPy, the sparse inputs $\\mathbf{Z}$ are abbreviated 'iip' , for inducing input. Plot the posterior\n",
      "of $u$ in the same manner as for the full GP:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "mu, var = m._raw_predict(Z) \n",
      "plt.vlines(Z[:,0], mu[:,0]-2.*np.sqrt(var[:,0]), mu[:,0]+2.*np.sqrt(var[:,0]),color='r')"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 24,
       "text": [
        "<matplotlib.collections.LineCollection at 0x1115b99d0>"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAEACAYAAABVtcpZAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAEVlJREFUeJzt3X+sZGV9x/H3x12xBYmEaEF+WIpCBZsQqNluNIZpNGaB\nFqSpRdJEhARJK6lNjSJosrd/WWpIGjRaQtTwhxZbQ80SMbI0TCWxRbf8cEVWWOumC+LWFsFfNfLj\n2z/ugVwvc+/OnTM7c3ef9yu5uefHc87z3Wfmzuee58zcTVUhSWrTi+ZdgCRpfgwBSWqYISBJDTME\nJKlhhoAkNcwQkKSG9Q6BJJ9Osi/JzlXaXJ/k4ST3Jzmzb5+SpOmYxpXAZ4AtK+1Mci7wmqo6BXg3\n8Mkp9ClJmoLeIVBVdwE/WqXJ+cBNXdu7gaOSHNO3X0lSf7O4J3A8sHfJ+iPACTPoV5K0H7O6MZxl\n6/6tCklaBzbOoI9HgROXrJ/QbfsVSQwGSZpAVS3/RXtss7gS2Aa8EyDJZuCJqto3qmFVrauvrVu3\nzr0Gazq06rIma5r2V1+9rwSS/ANwNvDyJHuBrcCLAarqhqq6Lcm5SXYDPwMu7dunJGk6eodAVV08\nRpsr+/YjSZo+PzG8isFgMO8SXsCaxrce67Km8VjT7GQac0rTkKTWSy2SdLBIQq3zG8OSpHXKEJCk\nhhkCktQwQ0CSGmYISFLDDAFJapghIEkNMwQkqWGGgCQ1zBCQpIYZApLUMENAkhpmCEhSwwwBSWqY\nISBJDTMEJKlhhoAkNcwQkKSGGQKS1DBDQJIaZghIUsN6h0CSLUl2JXk4yVUj9g+SPJnk3u7rw337\n1CFsYWHeFUhNSVVNfnCyAfgO8BbgUeAbwMVV9eCSNgPgr6rq/P2cq/rUokNEAj4PpLEloaoy6fF9\nrwQ2Aburak9VPQXcDFwwot3EBUqSDpy+IXA8sHfJ+iPdtqUKeEOS+5PcluT0nn1KkqZkY8/jx7lu\nvwc4sap+nuQc4IvAqaMaLiyZDx4MBgwGg57lSdKhZTgcMhwOp3a+vvcENgMLVbWlW78aeLaqrl3l\nmO8Bv1tVjy/b7j0BeU9AWqN53xPYAZyS5KQkhwEXAduWFXhMknTLm1gMnsdfeCpJ0qz1mg6qqqeT\nXAl8BdgAfKqqHkxyRbf/BuCPgT9L8jTwc+AdPWuWJE1Jr+mgaXI6SIDTQdIazXs6SJJ0EDMEJKlh\nhoAkNcwQkKSGGQKS1DBDQJIaZghIUsMMAUlqmCEgSQ0zBCSpYYaAJDXMEJCkhhkCktQwQ0CSGmYI\nSFLDDAFJapghIEkNMwQkqWGGgCQ1zBCQpIYZApLUMENAkhpmCEhSw3qHQJItSXYleTjJVSu0ub7b\nf3+SM/v2KUmajl4hkGQD8HFgC3A6cHGS05a1ORd4TVWdArwb+GSfPiVJ09P3SmATsLuq9lTVU8DN\nwAXL2pwP3ARQVXcDRyU5pme/kqQp6BsCxwN7l6w/0m3bX5sTevYrSZqCjT2PrzHbZZzjFhYWnl8e\nDAYMBoOJipKkQ9VwOGQ4HE7tfKka93V8xMHJZmChqrZ061cDz1bVtUva/D0wrKqbu/VdwNlVtW/Z\nuapPLTpEJODzQBpbEqpq+S/aY+s7HbQDOCXJSUkOAy4Cti1rsw14JzwfGk8sDwBJ0nz0mg6qqqeT\nXAl8BdgAfKqqHkxyRbf/hqq6Lcm5SXYDPwMu7V21JGkqek0HTZPTQQKcDpLWaN7TQZKkg5ghIEkN\nW38hsORtopKkA2v93RNwTrhtPv7SmnhPQJI0MUNAkhpmCEhSwwwBSWqYISBJDTMEJKlhhoAkNcwQ\nkKSGGQKS1DBDQJIaZghIUsMMAUlqmCEgSQ0zBCSpYYaAJDXMEJCkhhkCktQwQ0CSGmYISFLDNk56\nYJKjgc8DvwnsAf6kqp4Y0W4P8GPgGeCpqto0aZ+SpOnqcyXwQWB7VZ0K/Eu3PkoBg6o60wCQpPWl\nTwicD9zULd8EvG2VtunRjyTpAOkTAsdU1b5ueR9wzArtCrgjyY4kl/foT5I0ZaveE0iyHTh2xK4P\nLV2pqkpSK5zmjVX1WJJXANuT7Kqqu0Y1XFhYeG6BwWDAYDBYvXpJasxwOGQ4HE7tfKla6bV7Pwcm\nu1ic6/9BklcCd1bVa/dzzFbgp1V13Yh9VVWQwIQ16RDg4y+tSRKqauIp9z7TQduAS7rlS4AvLm+Q\n5PAkR3bLRwBvBXb26FOSNEV9rgSOBv4ReBVL3iKa5Djgxqo6L8nJwC3dIRuBz1bVR1Y4n1cC8vGX\n1qjvlcDEITBthoAAH39pjeY5HSRJOsgZApLUMENAkhpmCEhSwwwBSWqYISBJDTMEJKlhhoAkNcwQ\nkKSGGQKS1DBDQJIaZghIUsMMAUlqmCEgSQ0zBCSpYYaAJDXMEJCkhhkCktQwQ0CSGmYISFLDDAFJ\napghIEkNMwQkqWETh0CStyd5IMkzSc5apd2WJLuSPJzkqkn7kyRNX58rgZ3AhcBXV2qQZAPwcWAL\ncDpwcZLTevQpSZqijZMeWFW7AJKs1mwTsLuq9nRtbwYuAB6ctF9J0vQc6HsCxwN7l6w/0m2TJK0D\nq14JJNkOHDti1zVVdesY56+1FLOwsPDcAoPBgMFgsJbDJemQNxwOGQ6HUztfqtb0Ov3CEyR3Au+r\nqntG7NsMLFTVlm79auDZqrp2RNuqKkigZ006iPn4S2uShKpadV5+NdOaDlqpgB3AKUlOSnIYcBGw\nbUp96lC0deu8K5CaMvGVQJILgeuBlwNPAvdW1TlJjgNurKrzunbnAH8HbAA+VVUfWeF8XglI0hr1\nvRLoPR00LYaApEPKwsLi1wFmCEjSejSj17L1ck9AknQQMgQkqWGGgCQ1zBCQpIYZApLUMENAkhpm\nCEhSwwwBSWqYISBJDTMEJKlhhoAkNcwQkKSGGQKS1DBDQJIaZghIUsMMAUlqmCEgSQ0zBCSpYYaA\nJDXMEJCkhhkCktQwQ0CSGjZxCCR5e5IHkjyT5KxV2u1J8s0k9yb5+qT9SZKmb2OPY3cCFwI37Kdd\nAYOqerxHX5KkA2DiEKiqXQBJxmk+ViNJ0mzN4p5AAXck2ZHk8hn0J0ka06pXAkm2A8eO2HVNVd06\nZh9vrKrHkrwC2J5kV1XdNarhwsLCcwsMBgMGg8GYXUhSG4bDIcPhcGrnS1X1O0FyJ/C+qrpnjLZb\ngZ9W1XUj9lVVQQI9a5KkuZvRa1kSqmriKfdpTQeNLCDJ4UmO7JaPAN7K4g3llW3dOqWSJEn7M/GV\nQJILgeuBlwNPAvdW1TlJjgNurKrzkpwM3NIdshH4bFV9ZIXzVd+rEklaNw6SK4He00HTYghIOqQc\nJCHgJ4YlqWGGgCQ1zBCQpIYZApLUMENAkhpmCEhSwwwBSWqYISBJDTMEJKlhhoAkNcwQkKSGGQKS\n1DBDQJIaZghIUsMMAUlqmCEgSQ0zBCSpYYaAJDXMEJCkhhkCktQwQ0CSGmYISFLDDAFJatjEIZDk\no0keTHJ/kluSvGyFdluS7ErycJKrJi9VkjRtfa4EbgdeV1VnAA8BVy9vkGQD8HFgC3A6cHGS03r0\nKUmaoolDoKq2V9Wz3erdwAkjmm0CdlfVnqp6CrgZuGDSPiVJ0zWtewKXAbeN2H48sHfJ+iPdNknS\nOrBxtZ1JtgPHjth1TVXd2rX5EPDLqvrciHa1lmIWFhaeXx4MBgwGg7UcLkmHvOFwyHA4nNr5UrWm\n1+lfPTh5F3A58Oaq+sWI/ZuBhara0q1fDTxbVdeOaFt9apGkdSWBGbymJaGqMunxfd4dtAV4P3DB\nqADo7ABOSXJSksOAi4Btk/YpSZquPvcEPga8FNie5N4knwBIclySLwFU1dPAlcBXgG8Dn6+qB3vW\nLEmakl7TQdPkdJCkQ8qhPh0kSTr4GQKS1DBDQJIaZghIUsMMAUlqmCEgSQ0zBCSpYYaAJDXMEJCk\nhhkCktQwQ0CSGmYISFLDDAFJapghIEkNMwQkqWGGgCQ1zBCQpIYZApLUMENAkhpmCEhSwwwBSWqY\nISBJDds46YFJPgr8AfBL4LvApVX15Ih2e4AfA88AT1XVpkn7lCRNV58rgduB11XVGcBDwNUrtCtg\nUFVnHmwBMBwO513CC1jT+NZjXdY0HmuanYlDoKq2V9Wz3erdwAmrNM+k/czTenzQrWl867EuaxqP\nNc3OtO4JXAbctsK+Au5IsiPJ5VPqT5I0BaveE0iyHTh2xK5rqurWrs2HgF9W1edWOM0bq+qxJK8A\ntifZVVV39apakta7rVvnXcFYUlWTH5y8C7gceHNV/WKM9luBn1bVdSP2TV6IJDWsqiaecu/z7qAt\nwPuBs1cKgCSHAxuq6idJjgDeCvz1qLZ9/hGSpMlMfCWQ5GHgMODxbtO/VdWfJzkOuLGqzktyMnBL\nt38j8Nmq+kjfoiVJ09FrOkiSdHCb+SeGkxyV5AtJHkzy7SSbl+0fJHkyyb3d14cPcD2/vaSve7u+\n/2JEu+uTPJzk/iRnzrumWY9T1+fVSR5IsjPJ55K8ZESbmY3TODXNaZze29XzrSTvXaHNTMdpnLpm\nMVZJPp1kX5KdS7YdnWR7koeS3J7kqBWO3ZJkVzduV62TmvYk+WY3Xl8/wDW9vXuuP5PkrFWOXds4\nVdVMv4CbgMu65Y3Ay5btHwDbZl1X1/eLgMeAE5dtPxe4rVv+PeDf10FNMx0n4CTgP4GXdOufBy6Z\n5ziNWdOsx+l3gJ3ArwEbgO3Aq+f9fBqzrgM+VsCbgDOBnUu2/S3wgW75KuBvRhy3AdjdPeYvBu4D\nTptnTd2+7wFHz2icXgucCtwJnLXCcWsep5leCSR5GfCmqvo0QFU9XSP+1ATz+3DZW4DvVtXeZdvP\nZzG8qKq7gaOSHDPnmmC24/Rj4Cng8CQbgcOBR5e1mfU4jVMTzHacXgvcXVW/qKpngH8F/mhZm3k8\nn8apCw7wWNXi28N/tGzz8+PRfX/biEM3Aburak9VPQXcDFww55qeM/UxG1VTVe2qqof2c+iax2nW\n00G/BfwwyWeS3JPkxu4dREsV8IbuMvm2JKfPsL53AKM+73A8sPRF+BFW/4T0NK1U00zHqaoeB64D\n/gv4PvBEVd2xrNlMx2nMmmb9fPoW8KZuOuFw4DxeOAbzeD6NU9e8fvaOqap93fI+YFQgjhqz4+dc\nE6y/D8OueZxmHQIbgbOAT1TVWcDPgA8ua3MPi1MfZwAfA744i8KSHAb8IfBPKzVZtn7A76jvp6aZ\njlOSVwN/yeJl5nHAS5P86aimy9YP2DiNWdNMx6mqdgHXsvi3tb4M3As8O6LpTJ9PY9Y1l5+9pWpx\nTmPUWMztHSyr1ASLH4Y9EzgHeE+SN82uspHWPE6zDoFHgEeq6hvd+hdYDIXnVdVPqurn3fKXgRcn\nOXoGtZ0D/EdV/XDEvkeBE5esn8DoaYeZ1TSHcXo98LWq+t+qeprFt/6+YVmbWY/Tfmuax/Opqj5d\nVa+vqrOBJ4DvLGsyl+fT/uqa48/eviTHAiR5JfDfI9osH7MTWXw9mWdNVNVj3fcfAv/M4nTMPK15\nnGYaAlX1A2BvklO7TW8BHljaJskxSdItb2LxbayPc+BdDPzDCvu2Ae/satrM4rTDvhXazqSmOYzT\nLmBzkl/v+n0L8O1lbWY9TvutaR7PpyS/0X1/FXAhL5zOm8vzaX91zfFnbxtwSbd8CaOvQHYApyQ5\nqbtCvqg7bm41JTk8yZHd8nMfht25vN0BstJ9iLWP07Tvao9x1/sM4BvA/Sz+5nYUcAVwRbf/PSzO\nX94HfA3YPIOajgD+Bzhyybbna+rWP87iXff7WeHO/CxrmtM4fYDF0N7J4s2yw9bBOK1a05zG6atd\nTfcBv78enk/j1DWLsWLxl5rvs/j/kOwFLgWOBu5g8U/S3w4c1bU9DvjSkmPPYfHqZTdw9bxrAk7u\nxuq+btwOZE2XsXhzei/wf8APgC9PY5z8sJgkNcz/XlKSGmYISFLDDAFJapghIEkNMwQkqWGGgCQ1\nzBCQpIYZApLUsP8HbTRqKm8biDkAAAAASUVORK5CYII=\n",
       "text": [
        "<matplotlib.figure.Figure at 0x1115c5ad0>"
       ]
      }
     ],
     "prompt_number": 24
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Exercise 2\n",
      "\n",
      "a) Optimise and plot the model. The inducing inputs are marked \u2013 how\n",
      "are they placed? You can move them around with e.g. `m['iip_2_0'] = 100` . What\n",
      "happens to the likelihood? What happens to the fit if you remove an input?"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# Exercise 3 a answer"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 28
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "b) How does the fit of the sparse compare with the full GP? Play around\n",
      "with the number of inducing inputs, the fit should improve as $M$ increases. How many\n",
      "inducing points are needed? What do you think happens in higher dimensions?"
     ]
    },
    {
     "cell_type": "raw",
     "metadata": {},
     "source": [
      "# Exercise 3 b answer"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Exercise 3\n",
      "\n",
      "Can you build a low rank Gaussian process with the intrinsic model of coregionalization? Do you have to treat the 2nd input (which specifies the event number) in a special way?"
     ]
    }
   ],
   "metadata": {}
  }
 ]
}