{ "metadata": { "name": "", "signature": "sha256:0b9da3d8d8d7aa47f4390cebfe573c7132f48d0cec1ea5dbe1e117f7d389cdb7" }, "nbformat": 3, "nbformat_minor": 0, "worksheets": [ { "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Gaussian Process Regression in GPy\n", "\n", "## Machine Learning Summer School, Sydney, Australia\n", "\n", "### 23rd February 2015\n", "\n", "### Neil D. Lawrence and Nicolas Durrande\n", "\n", "In the introduction to `GPy` you saw how it was possible to build covariance functions with the GPy software. The covariance function contains the assumptions about the data in it. In the Gaussian process, the covariance funciton *encodes* the model. However, to make predictions, you need to combine the model with data. \n", "\n", "If the data we are given, $\\mathbf{y}$ is *real* valued, then the problem is known as a regression problem. If a Gaussian noise model is used, then this is known as Gaussian process regression.\n" ] }, { "cell_type": "code", "collapsed": false, "input": [ "%matplotlib inline\n", "import numpy as np\n", "import pods\n", "import pylab as plt\n", "import GPy\n", "from IPython.display import display" ], "language": "python", "metadata": {}, "outputs": [ { "output_type": "stream", "stream": "stderr", "text": [ "/Users/neil/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages/pytz/__init__.py:29: UserWarning: Module pods was already imported from /Users/neil/sods/ods/pods/__init__.pyc, but /Users/neil/Library/Enthought/Canopy_64bit/User/lib/python2.7/site-packages is being added to sys.path\n", " from pkg_resources import resource_stream\n" ] } ], "prompt_number": 1 }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will now combine the Gaussian process prior with some data to form a GP regression model with GPy. We will generate data from the function \n", "$$\n", "f(x) = \u2212 \\cos(\\pi x ) + \\sin(4\\pi x )\n", "$$ \n", "over $[0, 1]$, adding some noise to give $y(x) = f(x) + \\epsilon$, with the noise being Gaussian distributed, $\\epsilon \\sim \\mathcal{N}(0, 0.01)$. " ] }, { "cell_type": "code", "collapsed": false, "input": [ "X = np.linspace(0.05,0.95,10)[:,None]\n", "Y = -np.cos(np.pi*X) + np.sin(4*np.pi*X) + np.random.normal(loc=0.0, scale=0.1, size=(10,1)) \n", "plt.figure()\n", "plt.plot(X,Y,'kx',mew=1.5)" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "pyout", "prompt_number": 2, "text": [ "[]" ] }, { "metadata": {}, "output_type": "display_data", "png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEACAYAAAC08h1NAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAEpxJREFUeJzt3X+sZGddx/HPx93yB0SsTclW2sWCtPwwSgphWUHTMSBz\nZk1aURDqL0SDhKTmZKKxgMS9/qGEP8h0sFibBkhJDMUowUI7Zyy6A42BQi1dKnSxq2yyLbgRyu/y\nR5f5+sc9d3b29v6YmXPnzNz7vF/JJHNmnp3nu8/e+cyzz3nmHkeEAADp+LFFFwAAqBfBDwCJIfgB\nIDEEPwAkhuAHgMQQ/ACQmMrBb/sDts/YfnCT5xu2v2P7C+XtnVX7BADMbv8OvMYHJf2NpA9t0eZT\nEXHNDvQFAKio8ow/Iu6R9K1tmrlqPwCAnVHHGn9Iernt47bvsv3CGvoEAGxiJ5Z6tnO/pIMR8bjt\nlqSPSbqyhn4BABuYe/BHxPfG7vds/63tiyLisfF2tvmlQQAwg4iYajl97ks9tg/Ydnn/kCSvD/01\nEcEtQkePHl14DctyYyyqjcVwOFSe56P3WJ7nGg6HC/+78HOxc7dZVJ7x2/6wpKslXWz7tKSjki4o\ng/wWSa+V9FbbZyU9LukNVfsEAMyucvBHxHXbPP8+Se+r2g+A6USE2u22ut3uaNbf7XYlSZ1OR+V/\nxJGgOk7uYkqNRmPRJSwNxuKcacei3++PQr/T6Ywe73a7yrJMWZbtcIX14eeiGs+6RrTTbMey1ALs\nFUVRqNlsjmb3EaF+v7+rQx/ns62Y8uQuwQ8Au9gswc8vaQOAxBD8AJAYgh8AEkPwA0BiCH4ASAzB\nDwCJIfgBIDEEPwAkhuAHgMQQ/ACQGIIfABJD8ANAYgh+AEgMwQ8AiSH4ASAxBD8AJIbgB4DEEPwA\nkBiCHwASQ/ADQGIIfgBIDMEPAIkh+AEgMQQ/ACSG4AeAxFQOftsfsH3G9oNbtHmv7YdtH7d9VdU+\nAQCz24kZ/wclZZs9afuIpOdGxBWS/kjSzTvQJwBgRpWDPyLukfStLZpcI+m2su29ki60faBqvwCA\n2dSxxn+ppNNjx49IuqyGfgEAG6jr5K7XHUdN/QIA1tlfQx+PSjo4dnxZ+diTrKysjO43Gg01Go15\n1gUAu85gMNBgMKj0Go6oPvm2fbmkj0fEz23w3BFJ10fEEduHJd0YEYc3aBc7UQsApMS2ImL9qsqW\nKs/4bX9Y0tWSLrZ9WtJRSRdIUkTcEhF32T5i+6SkH0h6U9U+AQCz25EZ/05gxg8A05tlxs83dwEg\nMQQ/ACSG4AeAxBD8AJAYgh8AEkPwA0BiCH4ASAzBDwCJIfgBIDEEPwAkhuAHgMQQ/ACQGIIfABJD\n8AOYu6IoNP7bdyNCRVEssKK0EfwA5qooCrVaLbXbbUWEIkLtdlutVovwX5A6Lr0IIGHNZlN5nqvb\n7Y4e63a7yvNczWZzgZWli+AHMFe21el0JGkU/nmeq9PpyJ7q+iHYISz1AEBimPEDmKu1Nf215R3p\n3MyfWf9iEPwA5qrf749Cf23JR1oN/yzLlGXZAqtLExdbB+akKAo1m83RjDYi1O/3kww6xmJ+uNg6\nsCTYwni+LMvOW9KxTegvEEs9wBywhRHLjOAH5oAtjFhmLPUAQGKY8QNzwBZGLDOCH5gDtjBimbGd\nE5gTtjCiDrNs5yT4AWAXW8g+ftuZ7RO2H7Z9wwbPN2x/x/YXyts7q/YJAJhdpTV+2/sk3STpVZIe\nlfR523dExEPrmn4qIq6p0hcAYGdUnfEfknQyIk5FxBOSbpd07Qbt2MIAAEuiavBfKun02PEj5WPj\nQtLLbR+3fZftF1bsEwBQQdXtnJOcjb1f0sGIeNx2S9LHJF25UcOVlZXR/UajoUajUbE8ANhbBoOB\nBoNBpdeotKvH9mFJKxGRlcdvlzSMiHdv8We+KuklEfHYusfZ1QMAU1rErp77JF1h+3LbT5H0ekl3\nrCvqgMuNzLYPafXD5rEnvxQAoA6Vlnoi4qzt6yX1Je2T9P6IeMj2W8rnb5H0WklvtX1W0uOS3lCx\nZgBABXyBCwB2MS7EAgDYFsEPAIkh+AEgMQQ/ACSG4AeAxBD8AJAYgh8AEkPwA0BiCH4ASAzBDwCJ\nIfgBIDEEPwAkhuAHgMQQ/ACQGIIfABJD8ANAYgh+AEgMwQ8AiSH4ASAxBD8A1KgoCo1fXzwiVBRF\nrTUQ/DtgGf4hASy/oijUarXUbrcVEYoItdtttVqtejNjrfNF31ZL2X16vV5IijzPYzgcxnA4jDzP\nQ1L0er1FlwdgiYznQ57n590fDoczvWaZnVPl7f76PmL2pmazqTzP1e12R491u13lea5ms7nAygAs\nG9vqdDqSNMqMPM/V6XRku7Y6CP6KluUfEgAmRfADQE2iXNNfWxWQzk0Y65wsEvwVLcs/JIDl1+/3\nR1mxtlIgrWZGlmXKsqyeQqY9KTCvmzi5CyABvV7vvBO5w+GwUlZohpO7jrFtiItkO5allmkVRaFm\nszma3UeE+v1+fZ/eALa1V9+nthURUy0tVN7HbzuzfcL2w7Zv2KTNe8vnj9u+qmqfyybLsvOWdGzv\n+h8mYC9Zmv3zS6LSGr/tfZJukvQqSY9K+rztOyLiobE2RyQ9NyKusP0ySTdLOlylXwCYBtuuz1f1\n5O4hSScj4pQk2b5d0rWSHhprc42k2yQpIu61faHtAxFxpmLfADARtl2fr+pSz6WSTo8dP1I+tl2b\nyyr2CwCYUdUZ/6RnY9d/pG7451ZWVkb3G42GGo3GTEUBwLi1Nf29sO16MBhoMBhUeo1Ku3psH5a0\nEhFZefx2ScOIePdYm7+TNIiI28vjE5KuXr/Us5t39QBYbmsnd8f3z699EPR6vV29GWOWXT1Vg3+/\npK9IeqWkr0n6nKTrNji5e31EHCk/KG6MiCed3CX4AcwT2znH/kzVsLXdknSjpH2S3h8R77L9FkmK\niFvKNjdJyiT9QNKbIuL+DV6H4AeAKS0k+HcKwQ8A01vIF7gAALsLwQ8AiSH4ASAxBD8AJIbgB4DE\nEPwAkBiCHwASQ/ADQGIIfgBIDMEPAIkh+AEgMQQ/ACSG4AeAxBD8AJAYgh8AEkPwA0BiCH4ASAzB\nDwCJIfgBIDG7OviLotD4dXojQkVRLLAiAFh+uzb4i6JQq9VSu91WRCgi1G631Wq1CH8A2ML+RRcw\nq2azqTzP1e12R491u13lea5ms7nAygBgue3a4LetTqcjSaPwz/NcnU5HthdZGgAstV271AMAmM2u\nnfGvremvLe9I52b+zPoBYHO7Nvj7/f4o9NeWfKTV8M+yTFmWLbA6AFheHt8OuUi2Y9paiqJQs9kc\nze4jQv1+n9AHkAzbioipljh2dfADQOpmCf6Zl3psXyTpI5J+WtIpSb8ZEd/eoN0pSd+V9CNJT0TE\noVn7BABUV2VXz9sk3R0RV0r61/J4IyGpERFXEfoAsHhVgv8aSbeV92+T9GtbtGWLDQAsiSrBfyAi\nzpT3z0g6sEm7kPRJ2/fZfnOF/gAAO2DLNX7bd0u6ZIOn/nz8ICLC9mZnZl8REV+3/QxJd9s+ERH3\nbNRwZWVldL/RaKjRaGxVHgAkZzAYaDAYVHqNmXf12D6h1bX7/7X9U5KORcTzt/kzRyV9PyLes8Fz\n7OoBgCnNsqunylLPHZLeWN5/o6SPbVDQU23/eHn/aZJeLenBCn0CACqqMuO/SNI/SHqWxrZz2n6m\npFsj4ldtP0fSR8s/sl/S30fEuzZ5PWb8ADAlvsAFAImpe6kHALALEfwAkBiCHwASQ/ADQGIIfgBI\nDMGPPacoCo3vEIsIFUWxwIqA5ULwY08pikKtVkvtdlsRMbpEZ6vVIvyB0q699CLOx9XIVjWbTeV5\nPrr+sqTRJTqbzeYCKwOWyNqsaNG31VIwi16vF5Iiz/MYDocxHA4jz/OQFL1eb9Hl1W787z8+LsBe\nVGbnVHnLjH8PYJYLYBoE/x5gW51OR5JG4Z/nuTqdzmjpJxVRrumvffBJ58YkxfEANkLwY0/p9/uj\n0F/7MJRWwz/LsuTOeQAb4Ze07QGbzXJTnfVzohsp4bdzJmptC+P4LHftg6DX6xF4wB5G8CeMWS6Q\nJoIfABLD7+MHAGyL4AeAxBD8AJAYgh8AEkPwA0BiCH4ASAzBDwCJIfgBIDEEPwAkhuAHgMQQ/ACQ\nGIIfABIzc/Dbfp3tL9n+ke0Xb9Eus33C9sO2b5i1PwDAzqgy439Q0mskfXqzBrb3SbpJUibphZKu\ns/2CCn0CACqa+dKLEXFC0nZXdzok6WREnCrb3i7pWkkPzdovAKCaea/xXyrp9NjxI+VjAIAF2XLG\nb/tuSZds8NQ7IuLjE7w+V1YBgCWzZfBHxK9UfP1HJR0cOz6o1Vn/hlZWVkb3G42GGo1Gxe4BYG8Z\nDAYaDAaVXqPypRdtH5P0pxHxHxs8t1/SVyS9UtLXJH1O0nUR8aQ1fi69CADTq/XSi7ZfY/u0pMOS\n7rTdKx9/pu07JSkizkq6XlJf0pclfWSj0AcA1IeLrQPALsbF1gEA2yL4ASAxBD8AJIbgB4DEEPwA\nkBiCHwASQ/ADQGIIfgBIDMEPAIkh+AEgMQQ/ACSG4AeAxBD8AJAYgh8AEkPwA0BiCH4ASAzBDwCJ\nIfgBIDEEPwAkhuAHgMQQ/ACQGIIfABJD8ANAYgh+7KiiKBQRo+OIUFEUC6wIwHoEP3ZMURRqtVpq\nt9uKCEWE2u22Wq0W4Q8skf2LLgB7R7PZVJ7n6na7o8e63a7yPFez2VxgZQDGEfzYMbbV6XQkaRT+\neZ6r0+nI9iJLAzCGpR4ASMzMM37br5O0Iun5kl4aEfdv0u6UpO9K+pGkJyLi0Kx9YrmtremvLe9I\n52b+zPqBJbJ2Em7am1YD/0pJxyS9eIt2X5V00QSvF1h17NixRZcwk16vF5Iiz/MYDocxHA4jz/OQ\nFL1eb6bX3K1jMQ+MxTmMxTlldk6V3zMv9UTEiYj4rwmbM9WbwmAwWHQJM8myTL1ebzS7X1vz7/V6\nyrJsptfcrWMxD4zFOYxFNXWs8YekT9q+z/aba+gPC5Rl2XlLOrZnDn0A87HlGr/tuyVdssFT74iI\nj0/Yxysi4uu2nyHpbtsnIuKeaQsFAOwMx9i3LGd6AfuYpD+JTU7urmt7VNL3I+I9GzxXrRAASFRE\nTLWcvlP7+Dfs1PZTJe2LiO/ZfpqkV0v6y43aTls4AGA2M6/x236N7dOSDku603avfPyZtu8sm10i\n6R7bD0i6V9InIuJfqhYNAJhd5aUeAMDuUus3d21ntk/Yftj2DZu0eW/5/HHbV9VZX522Gwvbv12O\nwRdt/7vtn19EnXWY5OeibPdS22dt/3qd9dVpwvdIw/YXbP+n7UHNJdZmgvfIxbYL2w+UY/H7Cyhz\n7mx/wPYZ2w9u0Wa63Jx24/+sN0n7JJ2UdLmkCyQ9IOkF69ockXRXef9lkj5bV3113iYci1+Q9BPl\n/SzlsRhr92+SPiHpNxZd9wJ/Li6U9CVJl5XHFy+67gWOxYqkd62Ng6RvStq/6NrnMBa/JOkqSQ9u\n8vzUuVnnjP+QpJMRcSoinpB0u6Rr17W5RtJtkhQR90q60PaBGmusy7ZjERGfiYjvlIf3Srqs5hrr\nMsnPhST9saR/lPR/dRZXs0nG4rck/VNEPCJJEfGNmmusyyRj8XVJTy/vP13SNyPibI011iJWt79/\na4smU+dmncF/qaTTY8ePlI9t12YvBt4kYzHuDyXdNdeKFmfbsbB9qVbf9DeXD+3VE1OT/FxcIeki\n28fKL0X+bm3V1WuSsbhV0s/a/pqk45LymmpbNlPnZp2/lnnSN+v6bZ178U0+8d/J9i9L+gNJr5hf\nOQs1yVjcKOltERFe/VrwXt36O8lYXCDpxZJeKempkj5j+7MR8fBcK6vfJGPxDkkPRETD9s9o9Qui\nL4qI7825tmU0VW7WGfyPSjo4dnxQq59MW7W5rHxsr5lkLFSe0L1VUhYRW/1XbzebZCxeIun28ldB\nXCypZfuJiLijnhJrM8lYnJb0jYj4oaQf2v60pBdJ2mvBP8lYvFzSX0lSRPy37a9Kep6k+2qpcHlM\nnZt1LvXcJ+kK25fbfoqk10ta/8a9Q9LvSZLtw5K+HRFnaqyxLtuOhe1nSfqopN+JiJMLqLEu245F\nRDwnIp4dEc/W6jr/W/dg6EuTvUf+WdIv2t5XfkHyZZK+XHOddZhkLE5IepUklWvaz5P0P7VWuRym\nzs3aZvwRcdb29ZL6Wj1j//6IeMj2W8rnb4mIu2wfsX1S0g8kvamu+uo0yVhI+gtJPynp5nKmuyev\nZTDhWCRhwvfICduFpC9KGkq6NSL2XPBP+HPx15I+aPu4ViexfxYRjy2s6Dmx/WFJV0u6uPzS7FGt\nLvnNnJt8gQsAEsOlFwEgMQQ/ACSG4AeAxBD8AJAYgh8AEkPwA0BiCH4ASAzBDwCJ+X9CrSz25Rrd\nMgAAAABJRU5ErkJggg==\n", "text": [ "" ] } ], "prompt_number": 2 }, { "cell_type": "markdown", "metadata": {}, "source": [ "A GP regression model based on an exponentiated quadratic covariance function can be defined by first defining a covariance function, " ] }, { "cell_type": "code", "collapsed": false, "input": [ "k = GPy.kern.RBF(input_dim=1, variance=1., lengthscale=0.2)" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 3 }, { "cell_type": "markdown", "metadata": {}, "source": [ "And then combining it with the data to form a Gaussian process model," ] }, { "cell_type": "code", "collapsed": false, "input": [ "model = GPy.models.GPRegression(X,Y,k)" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 6 }, { "cell_type": "markdown", "metadata": {}, "source": [ "Just as for the covariance function object, we can find out about the model using the command `display(m)`. " ] }, { "cell_type": "code", "collapsed": false, "input": [ "display(model)" ], "language": "python", "metadata": {}, "outputs": [ { "html": [ "\n", "\n", "

\n", "Model: GP regression
\n", "Log-likelihood: -13.2481164705
\n", "Number of Parameters: 3
\n", "Updates: True
\n", "

\n", "\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "
GP_regression.ValueConstraintPriorTied to
rbf.variance 1.0 +ve
rbf.lengthscale 0.2 +ve
Gaussian_noise.variance 1.0 +ve
" ], "metadata": {}, "output_type": "display_data", "text": [ "" ] } ], "prompt_number": 7 }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note that by default the model includes some observation noise\n", "with variance 1. We can see the posterior mean prediction and visualize the marginal posterior variances using \n", "```python\n", "model.plot()\n", "```" ] }, { "cell_type": "code", "collapsed": false, "input": [ "_ = model.plot()" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "display_data", "png": "iVBORw0KGgoAAAANSUhEUgAAAW0AAAD7CAYAAAChScXIAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl03Nd93/33/c0+Awz2fSe4L5K4aiElUhZdS3E2x3Ia\nx6ofR00aNU968jztaeMkPbVz+uRx0566PWl9oiSN0zhy7LR2XFuP9SiWFFMmKYqLKO4EV+z7MgBm\nX36/2z8GBEGKJABigFnwfZ2DgwHmx5nvjxh85uLe+7tXaa0RQgiRH4xsFyCEEGLhJLSFECKPSGgL\nIUQekdAWQog8IqEthBB5REJbCCHyiH25n0ApJXMKhRDiIWit1d3fW5GWttY6Jz6+9KUvZb0GOSc5\np0L5kHNa3o/7ke4RIYTIIxLaQgiRR1ZVaB84cCDbJWScnFN+kHPKD/lwTupBfScZeQKl9HI/hxBC\nFBqlFDpbA5FCCCEyQ0JbCCHyyJJCWynlVkodV0qdUUpdUkp9JVOFCSGE+KglXVyjtY4ppZ7VWkeU\nUnbgiFJqn9b6SIbqE0IIMceSu0e01pGZm07ABkws9TGFEELc25JDWyllKKXOAMPAj7XWl5ZelhBC\niHvJREvb0lo/BjQCzyilDiy5KiGEEPeUsQWjtNZTSqkfAruAQ3Pv+/KXvzx7+8CBA3kxgV0IIVbS\noUOHOHTo0LzHLeniGqVUJZDSWk8qpTzA3wG/r7V+Z84xcnGNEEIs0v0urllqS7sO+EullEG6q+Wv\n5ga2EEKIzJLL2IUQIgctV0tbCJEDLCvdMNKAUmCoj/yuiwIhoS1EjosnTG70TzA4HiaZsjC1RcqE\nlGWBBSYAGgM1c0sDCkOBYSgM0p9tNoXDpihyO9nQUkFJkRubIeGeb6R7RIgcYlmaa70TdA1NEUuZ\nJFMaraGyzEdlqQ9jiSGrtSYSTzIwPEUyZWEzwGk3cNoNWmpKWNNQisNuy9DZiKW4X/eIhLYQWTYV\ninGqY4hQLEk8aVFZ5qO6vAi1gl0cltaMToQYC4SxGwqXw6C82MX29XW4XfIHeTZIaAuRQ4KRBMcu\n9BGKJUHZWNNYgcOeW4tuhqIJegcDGIDbaaOhqogtbVXYbblVZ6GS0BYiyyxLc/rKIL1jYSytaW+q\nyrmgvh+tNYHpKIMjU7icNoo8DvZsqqPY68p2aQVLQluILAlHk7x7podQLEV9dQllfk+2S1qyRNLk\nZt8EaBOP087m1kpaav0r2qVT6CS0hVhhQxNhTlwcIJ6yWN9ajT1PWtWLpbWmb3iaYDiK22GjvsLH\nY+tqsEk3ypJIaAuxQroGJ/nw2giGzU57U/mqa31OBWP0jUzidhiUeJ3s2lSHz+3Mdll5J6uhPRoI\nU1HiWXUvXrG69I9M8/7lQTxuJy315dkuJyfEEyk6+8YxFPhcNra111BXWZTtsnJaPGFyczDA5taq\n7IX2949cJ55M4XQYuO021jeW01JXsuQ5p0LkgtHJMEfO9WOz22lvqsh2OTnLsjS9g5PE4gncThuN\n1UVsbate1TkQjiW5dHOE8ek4SdMikbLQQFW5nwOP1GUvtC/0TM5+bWnN4Og0k8EoLruB22Gwua2K\nhqpiufRW5JVoPMnbJ7tIaVjXUiV/SS6C1pqJqQhD40HcDgOPw8bGlgoaqwt3MHM6HOdS5yiBUJyE\naZFKaVCKhpoSijzOO847kTTZ0V6RG6F9N9Oy6BuaIhSJ43YYFHkd7FxfS0mRe1nrEuJhWZbmJ2d7\nGZ6MsLG1pmAHGFeSZWn6R9I54LQrfC47m1orqatY2YuMMsGyNCOBMFd6xgnHUiRNi2TKwrDbaKop\nxeOyz3tOOR3ad0ukZqYSWSYeh521TaWsbShf1X9CidxxsXOUC51jtNSXU+yThsVySZkWfUOTRGKJ\n2cvsy4pcbG6rotjrzIkgtyzNVDjO1e4xpiIJkimLhKlJmRYet5P6aj8ux8NdTZpXoT2X1pqhsSCB\nqQguh0FFiZvdG+txOmRtBLGyAsEo73zQTUmxj/pqf7bLWXW01kRiSfqHp0iZJk6bgd2msNsMvC47\njdXF1FUU43LaMhbolqWJJ01GJkL0jQaJxJMkTU1qJphTlsblsFNfVYLHPX/reTHyNrTvFo7G6eoP\n4LApij0O9myux++TK7LE8kl3hfQwMhVnU1sVhiFdIblEa008ZTI5FSUwHcHSYFNgtxkYgDLSU+cM\nBVqDoQAFamZFRFOn10REpxfmMq30uJtlaUytsRkGfp+b8hIvTqdtxcbdCia050qmLG72jgMmPped\nXRvrqCjxZvx5xOrVNTTFsQsDNNeXUVqc/1cxivzxoNDO2+W7HHaDDW1VAJimxeFzA5iWiddlZ8f6\nWmrKfVmuUOSrZMrkzeM3QRk8sqE+2+UIcYe8De25bDaDda0zAW5ZHO8YJpFI4nPb2b2xjspSaYGL\nhTl7fYSLXWNsbK3BJUuSihxUcK9Km2HMXuBgWhZHLgxgmibFHgdPbGmQPnBxT9PhOD862Ym/2Muj\nGxqyXY4Q91VwoT2XzTBY15JugadSFm990INCU+Jz8uSWBrxuR5YrFNmmtebYhQF6R4Nsbq+RgUaR\n8/J2IHIp4skUN3rGsBuKCr+bJ7Y0yDTCVWhwPMRPzvZRX+WnvFTGQETuKMjZI5kSjibo6pvA5TCo\nr/SxY32tLClZ4EzT4u9PdxOKm6yf+UtMiFwiob1Ak8Eo/UOTuF021taXsam1IieuvBKZc70vwKmO\nIdY0V+DzyPiGyE3LNuVPKdUEfAOoBjTwp1rrP1rKY2ZTabGH0mIPWmuGJ8JcPHQVn8vG5rZKWutK\ns12eWIJYIsVbJzqxOx1sk2l8Io8tqaWtlKoFarXWZ5RSRcAHwM9rrS/POSZvWtr3kt6VY4pwJIbP\nZWfnhlqqyqT/M5+c6hjk+sAUG9uqcdhl7ELkvmVraWuth4ChmdshpdRloB64/MB/mEeUUjTVplvZ\npmlx7PIwqVSKEq+DxzfXUyQbm+askUCIn5ztp6K8mG3r6rJdjhAZkbEpf0qpVmA7cDxTj5lrbLbb\nc8ATSZO3TvdiaE2F38WezTIDJVeYlsWPP+giGLfYuKZGxiVEQclIaM90jXwH+C2tdSgTj5nrnA7b\n7MyDaDzJ60dv4HQo6sq97NhQJ/N9s+Ry1xhnr4+ytqWSOtmXUBSgJYe2UsoBfBd4TWv9v+51zNe+\n+pXZ27uf3MeeJ59e6tPmFI/Lwcb2GrTWBMMJvvvuNbwuG621frasqc52eavC+FSEd8/04vW5ZaBR\n5KUTxw5z8tgRAEzr/mONSx2IVMBfAuNa6//7Psfk9UDkw9JaMz4VYXQ8SJHbzsbmclrry7JdVsFJ\npizeOdVJJKlZ31IpXSGiICzbPG2l1D7gJ8A50lP+AH5Ha/3mnGNWZWjPZVmaofEQwVCEIreDR9ur\nqKmQHamXQmvN8Yv9dI+EWNtSidspSxKIwiEX1+QQ09L0D00SjSfwexzs2lgrazUv0tnrw1zpCVBf\nU0qZX/7vROEpyPW085XNUDTPdJMkUybvnhtAWxZ+j50ntjbKIlYPcLl7nItdY1SUFLFFpvCJVUpC\nO4scdhtrmyuB9CJWbxzvwmGDUp+Tx7c04HbKjwfSM0IudU/gL/KwaU1ttssRIqukeyQHReNJuvon\ncBhQ5HGwe9Pq2wtTa80HHYN0DQcpK/FRW1mc7ZKEWDHSp53HkimLrv5xLMvC57Kxrb2a+gIOsETS\n5Mi5XsanY9RVlVAm+36KVUhCu0BYlqZ3aJJINI7HZaeu3Msja2uwF8BSsp0DAS52jhFJmKxpqsTj\nkr59sXrJQGSBMAxFy5y53lPBKN87fA2n3cDjtLGppYLGan/ezFUen45y6vIAwWiKIp+b1qaqvKld\niGyR0M5jJcUeSmamC1qW5nx3gPcvDeFyGnicdtY1ltJcW4qRI0GotaZnaIqOngnC8RR2u53W+nIa\nC+AvBSFWioR2gTAMRUvd7Va4ZWmuDgQ5dWUUu03htBsUu+1saK2kstS7IkFuac3A6DQd3RNEEini\nSYtin4eGuvKceSMRIt9IaBcow1A0VPtpqPbPfi+eNDlxZYRILIHDMHDaFXabgcthUFHqpbW2FK/b\ngc1YXKAmUybj0zG6BgJMRxIkTYtkUpO0ND6vk4bq0oLodxciF8hA5Cqntca0NNOhGKMTIVKWxgCU\nAYZSGEqhVPo4QyksDWiNRmFqjbY0Fhqf20VVRREel11a0Utw5NDb7N3/3Gzfvtaao+++w74DB7Nc\nmVhJMhAp7ksphd2mKC/xUi7T67LqyKG3eeXzL/LSy6/w219Kr4z5h7//O7z29Vd59RvfkeAWgIS2\nEDlj7/7neOnlV3jt66/Ofu+1r7/KSy+/wt79z2WxMpFLJLSFyBFKqdkW9q3gvtXqlqmQ4hYZHRJC\niDwiLW0hcoTWerYP+6WXXwFut7iltS1ukdAWIkccffed2cC+1U0C6eDed+CgDEQKQKb8CZFTZMqf\nAFkwSggh8sqDQlsGIoUQIo9IaAshRB6R0BZCiDwioS2EEHlEQlsIIfKIhLYQQuSRJYe2UurrSqlh\npdT5TBQkhBDi/jLR0v4L4PkMPI4QQoh5LDm0tdaHgUAGahFCCDEP6dMWQog8siILRn3tq7cXv9n9\n5D72PPn0SjytEELkjRPHDnPy2BEATOv+y4tkZO0RpVQr8LrWets97pO1R4QQYhFk7REhhCgQmZjy\n9y3gPWC9UqpXKfUrSy9LCCHEvSy5T1tr/dlMFCKEEGJ+0j0ihBB5REJbiFXqyKG3mTsRQWvNkUNv\nZ7EisRAS2kKsQkcOvc0rn3+RP/z930FrPbup8Cuff1GCO8fJxr5CrEJ79z/HSy+/MrvbOzC7qfDe\n/c9lsTIxHwltIRag0DbcVUrN7vh+K7hv7QJ/6xxFbpLuESHmIV0JIpdIS1uIeeRSV4JlaaZDUcYn\nw0wGo0RjSaLxJLF4EsvSGEa6leyw2/B5nHg9Toq8LipKfZT5vdhs6XbarTeeW+dx65wAaW3nOAlt\nIeaRja4ErTWDo9Nc7Rqme2CC3sEAPYMBRiaCpFLWQz2moRSlfi/11X70dCc//ItX+cSnXuLX/q9/\nTZnfC6TPb9+Bg3nb7bMaSGgLkQO01vQNT3L6Yg/nrgzQ0TnEVDB2z2Pntpy9HgcelxO3y4FhqNnu\nm0TSJBxNEIkmCEZiTExGmAxGmJgKMzEVBtxU7/4CHbH1/B9f/CtK/R42tO7nc7/VRmnDZmKJJG6n\nY2X/E8SCSGgLMY/l6kowLYuL1wY58sENTl3oYWQieMf9JcVuNrbV0tZYQVNdGS315dRW+nG7Hi5M\nkymTickwfcOT9A0F6BncRO9ggO7BCSanoxw/1wXA0f/4fWyGwdqWSh7b1MT2TY1saKvBYbc91POK\nzJLQFmIeR999Zzawb3WTwMN1JWitudY9yjvHrnD0wxtMTkdn7yspdrNjczOPbWpkc3sttZX+jHa/\nOOw2air91FT62bml+Y6aBkam6OgcpuPmEB03h+nun+BK5whXOkf4mzc+wO2ys219Azs2N/H4o61U\nlxdnrC6xOBlZmvWBTyBLs4oCsNQpf6FInB8fv8qPjl6ms2989vt1VX6e3rmWJ7e30d5UNTuQmG2R\nWIILVwc409HHmct99AzeuTnVmqZKHn+klScebWVNU6UMXGbYg5ZmldAWYhkNj03z/b8/x4+OXiYW\nTwHg97l59on1PPv4etrzJPDGAiHOXO7jxPluTl/qmT0XgKqyIvY80sreHWvYsq4OmyEziZdKQluI\nFXa9e5TvvnWGox/cwJr5HXtkQwMvPL2ZJx5tw+HI3/7hRDLFuSv9vH+2ixPnupiYiszeV+r3sHd7\nO0/vamdze13O/OWQbyS0hVghHTeHeO0HJznT0QeAzTB4ZvdaPvXxR1nTWJnl6jLPsjTXukc4dqaT\no6dvMDg6PXtfeYmXvTvaeXpnOxvX1EqAL4KEthDL7GbfGH/1/ROcPN8NgMft4IWnN/Mzzz5CVXlR\nlqtbGVprbvSOceSDGxw+dZ3h8duzYSpKfezbkW6Bb2iryYsuoWyS0BZimfQPT/La6yc5fOo6AG6X\nnZ999hF+4eOPUeRzZbm67Lk1S+bIB9c5/MENRidCs/dVlRexb0c7+3a2s761WgL8HiS0hciwUDjO\nX//wJD88dBHTsrDbDX7qmS185vkds1cXijStNVc6hzn8wQ2OfHCD8cnw7H1VZUU8tX0N+3amW+DS\nhZImoS1EhpimxZuHL/HN108yHY6hFHz8qY189pO7V003yFJYVjrAf3LqOu99ePOOAC8v8bF3xxr2\n7ljDpvbaVT0LRUJbiAw429HHn/6Po3QPTACwdV09/+QX97KmqfAGGFeCZWmudA1z9PRNjp6+swul\nzO/lqe1t7N3Rzpa1dbMLXa0WEtpCLMFUMMp/+857/Pj4VQBqKov5x59+iicfa5P+2Ay51Qd+9PQN\njpy+wfDY7UHMIq+LnVua2fNICzs3Nz9wrKBQ1j2X0M5RhfICK1Raa945doU//+57BMNxnA4bv/RT\nO/n5g4/idMgKEMvl1iyUo6dvcOzDTvqGb+eHYSg2t9exZ1sLu7e10FhbOvv7c2vd87nLDdxaM+bV\nb3wnr36vJLRzUCG9wApR//AkX/vrdzl3ZQCA7Zsa+Y1ffoa6qpIsV7b69A9PcvJ8NyfOd3Px2iCm\ndXtp2vISH49ubOCxjQ08sqGBP/+jP/jIwl75uCOPhHYOut/Kcfn4AiskyZTJd390hr954wOSKRN/\nkZtf+8xeDuxZJz+THBCKxPnwUi8nznfz4aVeJoPRO+6vr/YTuvYm546+DuTvFmoPCu0l/42nlHoe\n+M+ADfhvWus/XOpjrgayR1/uuXR9kP/6zXdnF0c6+NRGXv6FJ/EXubNcWX7LZDdgkdfF07vW8vSu\ntWit6R6Y4GxHH2c6+rlwdYD+4Skm+m8vyPWjo5ex/vTv2Limlk1ramlvrsz7rq0lVa+UsgH/FTgI\n9AMnlVI/0FpfzkRxQqyEUCTOf//e+7x5+BIA9dUl/Obn9vPIhoYsV5b/lrMbUClFa0MFrQ0V/Nxz\nj5JMpfi9f/nPeaPrPdp3fILAdISRK+/y+jeTHN300yilsNsMmuvLaakvp62hgpaGctoaKyjze/Om\nsbTUt5w9wHWtdReAUurbwM8BEtrzkD36sk9rzdHTN/mTvzlCYDqC3Wbw4ie284sv7Mj71liuWMn9\nNY8fOcQb3/3G7BuEZWn+zRf/Bd//m79g34GDROwN9AxOcLN3jJu9Y/x4zr/1F7lpqS+nvrqEuqoS\naiv91Fb5qavy4/Pk1pWtS31lNgC9c77uAx5f4mOuCplcWF8s3sh4kD/+9uHZtUI2t9fym5/bT3N9\neZYrKywr2Q2478BBXv3Gd2a7Ymw2xf/z77/KC5/85OzvUySWoHtggq6+cbr6J+jqH6drYJzpUIzz\nVwc4f3XgI4/r97mpqSymrMRLeUl6m7fykvRHWYmPkiI3Pq8Lr9u5Ild0LjW0FzSK+eu//psYhoHN\npli3ZSebHtmN02HD6bDf+dlpx+ty4vE48LqdBb290d0vMEi3sCWwl5dpWvzgx+f55usniMVT+DxO\nvvCpJ/jEvs1yCXUBuPt3Ryl1x/e8biebZvq3b9FaMxYI0z0wztBYkMHRKYZGpxkcm2J4NMh0OMZ0\n+N77dc5lKIXX48TndVLkdeHzuHA57bgcdhwOGy6n/SO5Z7cZGIbCMAw6O05zs+MMigcH61JDux9o\nmvN1E+nW9p0HqS3pKlLQdTbMW2cPLejBnQ4bXrcTj9uJz+PE63Hg87goKfakP4rcc257KCl24y9y\nY7flR9jP9wITmXW9Z5T/8tohbvSMAfD0znZ+7Rf3Ul7iy3JlhSsfugGVUlSVF91zGQKtNRNTEUYn\ngkxMRZiYihCYCt++PR0hGIoRisaJxpKEInFCkTjDBO/xTAuxbt4jlhrap4B1SqlWYAD4h8Bn7z7o\ncz+zm2TSJJ5MkUyaJJIpEnd9jifSn2PxJJFYgnA0MXNf9CPTeuZT5vdSWVZEVblv5nMxVWXpH0pD\ndemqXn1tNYrEEnzz9ZO8/vfnsbSmqryI3/jsM+ze1pLt0gpevncDKqWoKPVRUTr/G7tpWoSjCUKR\nOOFonHAkQTyRvJ1xSZPkzOdEIp19KdPEsnT6Q9/+nEqZ/Pc37lPTUudpK6Ve4PaUvz/XWn/lrvsf\nap621pp4MkU0miQcjROJpcM8FI4zFYwyFYqmPwdjt2+HogRD8dmdQu7HX+SmvrqEhupSGmpKqK8u\npbE2/ZEvrXQxP601x8508qf/4whjgTCGUvzMx7bx0s/sweN+uB3NxeLJlb+Lt6ourjFNi4mpCGOB\nECMTQcYCIcYCYUYnQgyPTzMwMkU8kbrnv7XbDZrrymltKKetsZK2hgpaG8oplaU2887dA43rWqr4\nPz+3n7XNVVmuTIj5rarQno/WmvHJMAMjUzMfk/QNT9I7GLhjq6S5qsqK2NBWw/q2aja21dDeXIXL\nKVPCclHKNPlfb5/jWz88RTyRwut28vmff5wXntm8qpf6FPlFQnuBIrEE3f0TdPaP09k3RldfekpQ\nNJ684zibYdDWWMGGthq2rqtj6/p6Wfg+B5y/OsCr3z48u3TqM7vW8qufeUoGGkXekdBeAtOy6Bua\n5Ern8OxHz0DgI/3mjbWlbFtfz7b1DWxbV09ZiYT4Shkem+brf3uMo6dvAlBX5eeVX3qanVuas1yZ\nEA9HQjvDIrEE17tHuXxjiAvXBrh0Y+gj/eSNtaXs2NzEjs1NbF1fj9spA1+ZFo0l+c7fneZv3zpL\nMmXicth58fnt/MLHH5PuK5HXJLQXabGj3cmUyfWe0dkrqi7fGCQWvx3iDruNLWvr0iG+pYmW+vKc\nmJ+ar0zT4p33r/DaD04yMZXermr/7nV84VNPyJZfoiBIaC9CJha4SZkmVzpHOH2xh9OXerneM8rc\n/+byEh87tjSxa0szj21qpMgr88YXwrI0R0/f4LXXT9A/PAWkZ4X8k1/cx6b22nn+tRD5Q0J7EZZj\nneupYJQzHX18cLGHDy/1EZiOzN5nGIpNa2rZubWZXVuaaWuskFb4XbTWnLrQw1/94AQ3e9NXM9ZV\n+fnln97N/t3r5PJzUXAktBdpbnBDZhe40VrT2T/O6Yu9nLrQw+UbQx/ZiWPnliZ2bU23wnNthbGV\nZFoWxz7s5Ls/OsO17hEAKkp9fPaTuzj41Aa5EEoUrGXdBEEsjlKKNY2VrGms5MVPbCccjXPmcroV\nfupCDxNTYd56r4O33uvAZhhsaq9l19Zmdm1tXjV94fFEineOdfC9t8/Ozp33F7l58RPb+eT+rTLI\nKFY1aWnfJZvbgN1qhZ+60MMHF3q4fHMIy7r986ks87FzSzM7Z/rCvW7nstWSDcNj0/zovQ7ePHyR\nqWB6VbXaSj+fOvgozz21QWbgiFVDukcWIZc23A1F0q3wUxfTIT63L9xuM9i8tpZdW1rYubWZ5rqy\nvGyFJ5IpTpzr5q33LnP6Uu/sgO3a5io+/Q8e46nta7DZ5EpGsbpIaC9SLi5wY1mazr6xdCv8Yg8d\nN4fvuMCnqryIxzY2snVdHVvW1VNTUZyzIZ5Mmpy72s+7J69x7Ewn0Vj6ilOH3cbeHWt4/unNbFlb\nl7P1C7HcJLQLUDAc48PLfZy60M3pix/dlbqyzMeWtfVsaq9lXUsVbY0VWdtCS2vN6ESIDy/3cvJ8\nD2c7+u5YGqC9uZJn96znY09skE10hUBCu+BZluZm7xjnrw1w8doAF68PEgzH7zjGZhi0NJSztjkd\n4A01pTTWlFJZVpTxKXOhSJyu/nE6+8bpuDnExeuDjAXCdxzT2lDO3h3tPLNrLQ01pRl9fiHynYT2\nKmNZmp7BCS5eG+Rq1wjXekboHQxwrx+1y2GnrrqE2spiSv1eyvweSou9lJV4KfK6ZrdFcthtmJaF\naVqkTGt2l45gOJbe9Xo8yPB4kKHRaUYDoY88j8/rZMvaOnZtbWHX1maqy4uXfJ652I0lRCZIaAve\n+dGb1LQ9wvXuUXoGA/QNBeg4+z5WUWvGn8vpsNFSX05rQwXrWqrYvLaO5rryjLboc2nAWIhMk3na\nq9yRQ2/zW7/6Sx8JuJs/+RP+0599i+b1OxgLhAhMR5mcjjAZjBKYisxs+TazLVLKxLApbIaBzTDw\nehwUed0UeV2UFnuoriimuqKYmopiaiv9yz7jY+/+53jp5VdmL4CC21Mz9+5/blmfW4hsktBeBR4U\ncAf/wfMopdjQVpPFChdPKTX7BrQcV64KkasktFcBCTghCoeEtshL97tyFZA3I1HQJLRXgUIMuKPv\nvnPH8gK3vPb1V9l34KAMRIqCJaG9ChRiwO07cJBXv/GdO6b8/faXvpK35yPEQsmUv1VC5jQLkT9k\nnrYQQuSRB4X2Q0+mVUp9Ril1USllKqV2LK1EIYQQC7GUKyDOA58CfpKhWoQQQszjoQcitdYdQF7O\nPBBCiHwlq8sLIUQeeWBLWyn1FlB7j7t+V2v9+kKf5GtfvT3NbPeT+9jz5NMLLlAIIVaDE8cOc/LY\nEQBM6/4TRJY8e0Qp9WPgX2itT9/nfpk9IoQQi7Ass0fuIh3bQgixApYy5e9TSqle4Angh0qp/z9z\nZQkhhLiXpcwe+R7wvQzWIoQQYh4ye0QIIfKIhLYQQuQRCW0hhMgjEtpCCJFHJLSFECKPSGiLjDty\n6G3mXrSltebIobezWJEQhUNCW2TUkUNv88rnX+QPf/930FrPbnX2yudflOAWIgNkuzGRUXv3P8dL\nL78yuwclMLvV2d79z2WxMiEKg4S2yCil1Ow+lLeC+9belLKMrxBLJ90jQgiRR6SlLTLqVh/2rS4R\nuN3ilta2EEsnoS0y6ui778wG9q1uEkgH974DB2X3dyGWSHZjFxl35NDb7N3/3GyrWmvN0XffkcAW\nYoEetJ62hLYQeUZrjaUhFk8RjiWIxhLEEym0pVFKoRTY7Ta8HidFXhcuhw2boaRrKo88KLSle0SI\nHGVpzVjVcRK6AAAO5UlEQVQgzOhECJuhsNsUdruBXSmcdgOv20F9iYviuiJ8XheGUsCtQE8yHYoT\nCMUIBJMkTYukqTFNTcrSpEyNy2mntsqP12WXQM8jEtpC5AjTtOgZmiQWT+KyG7idBk3VJezbUovD\nblvUY5X4XNSUF933fq01gWCMqz3j9E8kSKQs4ikLu91GU20pbqdEQ66Sn4wQWRSNJbnZP4HLrvA6\nbexcV01dxf3DNlOUUpT7PTyxtfGO7weCUS7cGGUomiCesMAwaG0ow+WQqMgV0qctxApLpExudI9h\ntykq/S52b6rHlaMt2+lwnNNXBglGk0QTFhWlPmoqiqQ7ZZlJn7YQOaB/ZJrpYITSIicvPN6Gz+PI\ndknz8vtcHNjRCqT72K/3THC9b5xoIoXH46KprnSmL12sFAltIZaRaVrc7BsHrdnaWsHanU3ZLumh\nGUqxvqWC9S0VAAyMBTl3Y4RwLIXX7aJRAnxFSGgLsQwSKZNr3aP4XDaeeaSBcr8n2yVlXH1lMfWV\nxQD0j05z7sYokVgKv99HXaV0oSwXCW0hMiiRNLnaNUKpz8HPPtW+amZhNFT5aajyA3Ctd4IrvWNE\n4ynqqkspK8A3rGySgUghMiCZsrjWNYLfa+fZHa04HYuboleILK05c3WI3pEQKQvamytx2GWNuoWQ\nKyKFWCamZXG1cxSvy8azO5rxuHJ/cDEbwrEExy70MxVOUOT10FDjl+6TB1iW2SNKqf8A/DSQAG4A\nv6K1nnr4MoXIH1prugYCpBJJPrG7hWKfK9sl5TSf28nBXW0AdA4GuHhzlGjCpLmhnCKP/N8txkO3\ntJVSHwfe0VpbSql/B6C1/uI9jpOWtigoIxMhxgMh9myqpbmmJNvl5K2UaXH8Uj/DgSgul5OWulJp\nfc9Ylpa21vqtOV8eBz79sI8lRD6IxBLc6BljQ2MZzx7YkO1y8p7dZrB3W3oK5MBYkA+vDhOOp2it\nr8DndWa5utyVqaHtl4FvZeixhMgpWmtu9k7gdsBnDmzAZpPBtEy7NX3QNC3ev9TPleEAbqeT5npp\nfd/tgaGtlHoLqL3HXb+rtX595pjfAxJa67++3+N87au3F8Pf/eQ+9jz59MNVK8QKG5+KMDI2xb5t\nDdRWFGe7nIJnm9P67huZ5uz19MU7rU0V+NyF3fo+cewwJ48dAcC07t9tvaTZI0qpLwC/BjyntY7d\n5xjp0xZ5J2VaXOsepbHSxxNbGrJdzqqWTFkcu9DHyGQUr8dF8yro+16u2SPPA/8S2H+/wBYiHw2M\nTBOORHlhTys+T2G37vKBw27wzGPNAPQMT3H22giRpEl7YyUe9+qbYrmU2SPXACcwMfOtY1rr37jH\ncdLSFnkhnkhyvXuUbWuq2NRame1yxAMkUxbvne9lbCqGpwBb33JxjRAPoLWmeyCADYuDu9uwy0Bj\nXukemuLs9REiCZM1jRUF8deRLM0qxH1MBWP0Dwd4fEsdzdUy5zoftdSW0FJbQiJp8t6FfnoGJgq6\n71tCW6xKltZc6Rqh2u/mxQMbCvKXe7VxOmwc2D7T9z00xZkbI0TjJm2N5fgK6KpLCW2x6gyPhwhM\nhji4s4XSYne2yxHLoLm2hObaEpKpdOu7dyCAy+2ktb4s79+gJbTFqpFMmXR0jrChoZSPPSZXNK4G\nDruN/TMzT/pGpzlzdZhIwqS5vpxib362viW0xarQPRDATCb5+X1rV80a1+JOjVV+Gqv8pEyLE5cG\nuDI8CcpgbVNFXl3lKq9eUdCC4Tjd/ePs3FDL2saybJcjcoDdZvDUtvQu9IFglOOXBghGUpT6vdRX\n+7Nc3fxkyp8oSFprrnSOUlrk4NkdLbJ3oXggrTXX+gJ0dI8TTZq01JVndbldmfInVpWhsSCBqTDP\nbm+iosSb7XJEHlBKsb6pnPVN5aRMi1OXB7jeNUU8ZdHeXJlTXWq5U4kQSxRLpLjWOcLGlnKe2y4D\njeLh2G0GT2xNd59E4ymOX+ynOxTH0oq1zRU47NndSm5Vh7bWGkunP2uAmdsAylAo0u/ASiF/Xucw\ny7K40jlKkdvOpw9skH0IRcZ4XHYO7GgBIBhJ8P7FPoKRJBg22hsrsGfhtVZwoa21JpmyGAuEmQpF\n0VpjGAY2BYahsBkKm0qHsc1QKENhN8BAzQY1gKXTF2BYWmOZmpQ1c1uDaVrpz5bG1OB12amu9ONx\n2jEMCfeV1DMYIBqN87EdTZT7pStELJ9ir5OP714DwGQoxvGL/YTjJqalWdO0cl0oeTsQmUiaDE+E\nCIZiGArsNoXdbuAwFB6HjYZqP3WVxbictmVtJZuWZmIqSudggOlwnHjKIpnSmKZGGQb1NX6KPM68\nn9Cfa8YCYYbHp3lsbTXrm8qzXY5YxSLxFCcupTctjicsGmpLl3zRVl4vGGVamvFAmNHJMDYDnDYD\nu03hcdpY01BGXUVR1vuY7icUTXKle5Sx6RiJpEUiZWF32GmsLcHtKLg/clbEVChK72CAtfWl7NhQ\nK2+GIqekTIuz14YZnAgTS1o4bDbaGssXPQ88L0Jba004mqB/eBrLsnDYFQ67gctho6mqmDUNZQWx\n+trEdJSLN0cJRpPEUyaplKaxrgy/7Ob9QOFonM6+CRqqinhqa4OMMYicp7VmbCrKmWtDROIm8ZRF\nVamPqvKieRsbORfaKdNicGyaYCiOMdN6dtoNKvxuNrRUUuRxrJoWVCJpcu76MCOTUeIpE9OClvpy\nvKtwcfd7mQpF6RuaorrUw9OPNObVlWtCzGVZmqs943QPTxNLmiRSFkVeFw01JdiMO1/XWQ/t/3no\nKmiNw57ud3bbbbQ3ltFYVVwQredMisSSnLk2TCAYI5a0sNlstDSU4czRLqDlMjQWZGIyRGNVEY9v\nbpABXlFwLEvTOzzN1b5xogmTREqjlKKpthS73cbObIZ2MBLH5149redMmgzFOH11iFA0SSxhzawT\n/NF35kJgWRY3+wKYqSTrGsvYuqZKXjNiVQlGElzsHCUaS/KxXW3ZC+3lfo7VwtKawbEgF26Ozrwz\nW5SX+KitLM7rcJsMRukbmsTnsrNnUx3V5b5slyRE1imlJLQLjWlpbvRPcL1vklgyRTIF9dV+yvye\nnA/xUDRBd98ELodBfaWPHetrpb9aiDkktFeBlGlxsXOUgbEQ8aRFwrTweVw01pRkfexAz/yVEJiM\n4HYaVJa42bWxHqdjdfXVC7FQEtqrkGVphsaDXO4eJzYz0JE0TZwOB811pTgdtmVrkceTKfqHp4jF\nkzhsRnpefX0p7Y1lMl1PiAWQ0BZAOsgngjE6ukYJx1IkTQtz9jL99DEet4PSYg8+rxObYXD3q0aT\n3gUmGI4TiiSIxhIz67OAzaZw2gy8TjubWiupLPNKSAvxECS0xbwsS5NImUyGYoxMhJkOx0nNrLOi\nNaBIL6JFehPVsmI3VaVeSorc2O2GhLMQGSShLYQQeeR+of3Qo1NKqX+rlDqrlDqjlHpHKdW0tBKF\nEELM56Fb2kqpYq11cOb2PwMe1Vr/6j2Ok5a2EEIsUsZb2rcCe0YRMPawjyWEEGJhlrQ+qFLqD4B/\nBESAJzJSkRBCiPt6YGgrpd4Cau9x1+9qrV/XWv8e8HtKqS8C/wn4lXs9zpe//OXZ2wcOHODAgQMP\nW68QQhSkQ4cOcejQoXmPy8jsEaVUM/CG1nrrPe6TPm0hhFik5Zg9sm7Olz8HfPiwjyWEEGJhltKn\n/RWl1AbABG4A/zQzJQkhhLgfubhGCCFyUMa7R4QQQqw8CW0hhMgjEtpCCJFHJLSFECKPrKrQXsjE\n9Xwj55Qf5JzyQz6ck4R2npNzyg9yTvkhH85pVYW2EELkOwltIYTIIytycc2yPoEQQhSorGw3JoQQ\nInOke0QIIfKIhLYQQuSRggttpdTzSqkOpdQ1pdRv3+eYP5q5/6xSavtK17hY852TUupzM+dyTil1\nVCn1SDbqXIyF/JxmjtutlEoppX5hJet7WAt8/R1QSn2olLqglDq0wiUu2gJef5VKqTdnNvm+oJT6\nQhbKXDCl1NeVUsNKqfMPOCZ3M0JrXTAfgA24DrQCDuAMsOmuY36K9IYNAI8D72e77gyc05NAyczt\n5wvhnOYc9/fA/wd8Ott1Z+hnVQpcBBpnvq7Mdt0ZOKcvA1+5dT7AOGDPdu0POKenge3A+fvcn9MZ\nUWgt7T3Ada11l9Y6CXyb9AYNc/0s8JcAWuvjQKlSqmZly1yUec9Ja31Maz018+VxoHGFa1yshfyc\nAP4Z8B1gdCWLW4KFnNcvA9/VWvcBaK1zfUPshZzTIOCfue0HxrXWqRWscVG01oeBwAMOyemMKLTQ\nbgB653zdN/O9+Y7J5ZBbyDnN9Y+BN5a1oqWb95yUUg2kw+GPZ76VD9OcFvKzWgeUK6V+rJQ6pZT6\nRytW3cNZyDn9GbBFKTUAnAV+a4VqWy45nRFL2o09By30F/vuuY+5HAgLrk0p9SzwMrB3+crJiIWc\n038Gvqi11kopxUd/ZrloIeflAHYAzwFe4JhS6n2t9bVlrezhLeScfhc4o7U+oJRqB95SSj2qtQ4u\nc23LKWczotBCux9omvN1E+l3yQcd0zjzvVy1kHNiZvDxz4DntdYP+tMvFyzknHYC307nNZXAC0qp\npNb6BytT4kNZyHn1AmNa6ygQVUr9BHgUyNXQXsg5PQX8AYDW+oZSqhPYAJxakQozL6czotC6R04B\n65RSrUopJ/APgbt/yX8AfB5AKfUEMKm1Hl7ZMhdl3nNSSjUDfwu8pLW+noUaF2vec9Jar9Fat2mt\n20j3a//THA9sWNjr7/vAPqWUTSnlJT3QdWmF61yMhZxTB3AQYKbvdwNwc0WrzKyczoiCamlrrVNK\nqd8E/o70qPefa60vK6V+feb+P9Fav6GU+iml1HUgDPxKFkue10LOCfg3QBnwxzMt06TWek+2ap7P\nAs8p7yzw9dehlHoTOAdYwJ9prXM2tBf4s/p/gb9QSp0l3RD8V1rriawVPQ+l1LeA/UClUqoX+BLp\nbqu8yAi5jF0IIfJIoXWPCCFEQZPQFkKIPCKhLYQQeURCWwgh8oiEthBC5BEJbSGEyCMS2kIIkUck\ntIUQIo/8b3VEnN2LXhrTAAAAAElFTkSuQmCC\n", "text": [ "" ] } ], "prompt_number": 8 }, { "cell_type": "markdown", "metadata": {}, "source": [ "The actual predictions of the model for a set of points `Xstar`\n", "(an $m \\times p$ array) can be computed using \n", "\n", "```python\n", "Ystar, Vstar, up95, lo95 = model.predict(Xstar)`\n", "```\n", "\n", "### Exercise 1\n", "\n", "What do you think about this first fit? Does the prior given by the GP seem to be\n", "appropriate?" ] }, { "cell_type": "raw", "metadata": {}, "source": [ "# Exercise 1 answer" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Exercise 2\n", "\n", "The parameters of the models can be modified using the parameter name, for example \n", "```python\n", "model.Gaussian_noise.variance = 0.001\n", "```\n", "Change the values of the parameters to obtain a better fit." ] }, { "cell_type": "code", "collapsed": false, "input": [ "# Exercise 2 answer here" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 12 }, { "cell_type": "markdown", "metadata": {}, "source": [ "As we saw when we introduced GPy and covariance functions, random sample paths from the conditional GP can be obtained using\n", "```python\n", "np.random.multivariate_normal(mu[:,0],C)\n", "``` \n", "Now you can sample paths from the *posterior* process by first obtaining the mean and covariance of the posterior, `mu` and `C`. These can be obtained from the `predict` method, \n", "```python\n", "mu, C, up95, lo95 = model.predict(Xp,full_cov=True)\n", "```\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Exercise 3\n", "\n", "Obtain 10 samples from the posterior sample and plot them alongside the data below." ] }, { "cell_type": "code", "collapsed": false, "input": [ "# Exercise 3" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 13 } ], "metadata": {} } ] }