{ "metadata": { "name": "", "signature": "sha256:1d21006562db74bc6d75184f507a84bc08b02d3de17d4d78aaba614aeeb5a878" }, "nbformat": 3, "nbformat_minor": 0, "worksheets": [ { "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Gaussian Process Regression in GPy\n", "\n", "# Gaussian Process Summer School, Melbourne, Australia\n", "### 25th-27th February 2015\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": 5 }, { "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": 6, "text": [ "[]" ] }, { "metadata": {}, "output_type": "display_data", "png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEACAYAAAC08h1NAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAEpVJREFUeJzt3X+MbGddx/H31976B0SsTUkrbbGgLT+MkkK4VNB0DMjM\nXJNWFIT6C9EgIamZTDQWkHjXP5TwBxkGQWwaIDUxFKMEC+3OWPQONAaqtbRU6NVepcltwUYov8sf\nLfv1j52du93uzzm7Z2bneb+SSebMPHvOt092P/fpc55zTmQmkqRy/MCsC5Ak1cvgl6TCGPySVBiD\nX5IKY/BLUmEMfkkqTOXgj4gPRsTDEXHvFt83IuKbEfG58evtVY8pSZrekX3Yx4eAvwD+eps2n8rM\nq/bhWJKkiiqP+DPzduDrOzSLqseRJO2POub4E3hpRNwTEbdGxPNrOKYkaQv7MdWzk7uAizPz0Yho\nAx8DLqvhuJKkTRx48Gfmt9e9X46Iv4yIczPzkfXtIsKbBknSFDJzT9PpBz7VExHnR0SM3x8FYmPo\nr8lMX5kcP3585jXMy8u+sC/si+1f06g84o+IDwNXAudFxGngOHD2OMivB14NvDkiHgceBV5X9ZiS\npOlVDv7MvGaH798HvK/qcSRJ+8Mrd+dQo9GYdQlzw744w744w76oJqadI9pvEZHzUou0KAaDAc1m\nk/FpNjKT4XBIq9WacWXaLxFBztvJXUmzMRgMaLfbdLvdyYnAbrdLu91mMBjMujzNUB3r+CXNQLPZ\npNPp0O/3J5/1+306nQ7NZnOGlWnWDH5pQUUEvV4PYBL+nU6HXq83mfpRmZzqkaTCOOKXFtTanP7a\n9A6cGfk76i+bwS8tqOFwOAn9tSkfWA3/Vqvlyp6CuZxTWmAu51x80yznNPgl6RBzHb8kaUcGvyQV\nxuCXpMIY/JJUGINfkgpj8EtSYQx+SSqMwS9JhTH4JakwBr8kFcbgl6TCGPySVBiDX5IKY/BLUmEM\nfkkqjMEvSYUx+CWpMAa/JBWmcvBHxAcj4uGIuHebNu+JiPsj4p6IuLzqMSVJ09uPEf+HgC2f3BwR\nx4CfyMxLgd8D3r8Px5QkTaly8Gfm7cDXt2lyFXDjuO0dwDkRcX7V40qSplPHHP+FwOl12w8CF9Vw\nXEnSJuo6uRsbtrOm40qSNjhSwzEeAi5et33R+LMnWVpamrxvNBo0Go2DrEuSDp3RaMRoNKq0j8is\nPviOiEuAj2fmT23y3THg2sw8FhFXAO/OzCs2aZf7UYsklSQiyMyNsyrbqjzij4gPA1cC50XEaeA4\ncDZAZl6fmbdGxLGIOAV8F3hD1WNKkqa3LyP+/eCIX5L2bpoRv1fuSlJhDH5JKozBL0mFMfglqTAG\nvyQVxuCXDshgMGD9SrXMZDAYzLAiaZXBLx2AwWBAu92m2+2SmWQm3W6Xdrtt+Gvm6rhlg1ScZrNJ\np9Oh3+9PPuv3+3Q6HZrN5gwrkwx+6UBEBL1eD2AS/p1Oh16vR8SerrWR9p1TPZIOnOc75ovBLx2A\ntTn9temdtWmftTn/kni+Y/441SMdgOFwOAn9tSkfWJ32abVatFpbPq104Xi+Y/54kzbpgAwGA5rN\n5mROPzMZDodFhf6a9f8HBJ7v2E8zuS2zpM1tDPiIKDL0NX8MfkkHauP5Djiz0slR/2wY/JIOlOc7\n5o9z/JIOnOc7Ds40c/wGvyQdYj6BS5K0I4Nfkgpj8EtSjebh9hUGvyTVZF5uX+FyTkmqybzcvsLg\nl6SazMvtup3qkaTCOOKXpJrMy+0rDH5Jqsm83L7CK3clqUb7ffsKb9kgSYWZyS0bIqIVEScj4v6I\nuG6T7xsR8c2I+Nz49faqx5QkTa/SHH9EnAW8F3gF8BDwbxFxc2bet6HppzLzqirHkiTtj6oj/qPA\nqcx8IDMfA24Crt6knU9akKQ5UTX4LwROr9t+cPzZegm8NCLuiYhbI+L5FY8pSaqg6nLO3ZyNvQu4\nODMfjYg28DHgss0aLi0tTd43Gg0ajUbF8iRpsYxGI0ajUaV9VFrVExFXAEuZ2RpvvxVYycx3bvMz\nXwJelJmPbPjcVT2StEezWNVzJ3BpRFwSET8IvBa4eUNR58d4wWpEHGX1H5tHnrwrSVIdKk31ZObj\nEXEtMATOAj6QmfdFxJvG318PvBp4c0Q8DjwKvK5izZKkCryAS5IOMZ+5K0nakcEvSYUx+CWpMAa/\nJBXG4Jekwhj8klQYg1+SCmPwS1JhDH5JKozBL0mFMfglqTAGvyQVxuCXpMIY/JJUGINfkgpj8EtS\nYQx+SSqMwS9JhTH4JakwBr8kFcbgl6TCGPySVBiDX1IRBoMBmTnZzkwGg8EMK5odg1/SwhsMBrTb\nbbrdLplJZtLtdmm320WG/5FZFyBJB63ZbNLpdOj3+5PP+v0+nU6HZrM5w8pmw+CXtPAigl6vBzAJ\n/06nQ6/XIyJmWdpMONUjSYVxxC9p4a3N6a9N78CZkX+Jo/7KI/6IaEXEyYi4PyKu26LNe8bf3xMR\nl1c95rxxtYA034bD4ST0e70evV5vMuc/HA5nXV791s5wT/MCzgJOAZcAZwN3A8/b0OYYcOv4/UuA\nz26xrzyMlpeXE8hOp5MrKyu5srKSnU4ngVxeXp51eZLGlpeXc2VlZbK9srKyEH+j4+zcU3ZXneo5\nCpzKzAcAIuIm4GrgvnVtrgJuHCf7HRFxTkScn5kPVzz2XHC1gHQ4tFqtJ2xHxJM+K0XV4L8QOL1u\n+0FWR/U7tbkIWIjgd7WApMOmavDnzk0A2JiAm/7c0tLS5H2j0aDRaExVlCQtqtFoxGg0qrSPyNxt\ndm/ywxFXAEuZ2RpvvxVYycx3rmvzV8AoM28ab58Ertw41RMRWaWWWcktVgs46pdUh4ggM/cUNFVH\n/HcCl0bEJcCXgdcC12xoczNwLXDT+B+KbyzK/D48ebXAmn6/T6vVKnYOUdL8qjTiB4iINvBuVlf4\nfCAz3xERbwLIzOvHbd4LtIDvAm/IzLs22c+hHPHD6nLOZrM5Gd1nJsPh0NCXdOCmGfFXDv79cpiD\nX5JmZZrg95YNklQYg1+SCmPwS1JhDH5JKozBL0mFMfglqTAGvyQVxuCXpMIY/JJUGINfkgpj8EtS\nYQx+SSqMwS9JhTH4JakwBr8kFcbgl6TCGPySVBiDX5IKY/BLUmEOdfAPBgPWP6c3MxkMBjOsSJLm\n36EN/sFgQLvdptvtkplkJt1ul3a7bfhL0jaOzLqAaTWbTTqdDv1+f/JZv9+n0+nQbDZnWJkkzbdD\nG/wRQa/XA5iEf6fTodfrERGzLE2S5tqhneqRJE3n0I741+b016Z34MzI31G/JG3t0Ab/cDichP7a\nlA+shn+r1aLVas2wOkmaX7F+OeQsRUTutZbBYECz2ZyM7jOT4XBo6EsqRkSQmXua4jjUwS9JpZsm\n+Kee6omIc4GPAD8GPAD8amZ+Y5N2DwDfAr4PPJaZR6c9piSpuiqret4C3JaZlwH/NN7eTAKNzLzc\n0Jek2asS/FcBN47f3wj80jZtXWIjSXOiSvCfn5kPj98/DJy/RbsEPhkRd0bEGyscT5K0D7ad44+I\n24ALNvnqj9dvZGZGxFZnZl+WmV+JiKcDt0XEycy8fbOGS0tLk/eNRoNGo7FdeZJUnNFoxGg0qrSP\nqVf1RMRJVufu/zcifhQ4kZnP3eFnjgPfycx3bfKdq3okaY+mWdVTZarnZuD14/evBz62SUFPiYgf\nGr9/KvBK4N4Kx5QkVVRlxH8u8LfAM1m3nDMingHckJm/GBHPBj46/pEjwN9k5ju22J8jfknaIy/g\nkqTC1D3VI0k6hAx+SSqMwS9JhTH4JakwBv+CGAwGrD85npk+dF7Spgz+BTAYDGi323S7XTJz8nSy\ndrtt+Et6kkP7BC6d0Ww26XQ6k0dPApOnkzWbzRlWJmkeGfwLICImj59cC/+1R1L67GFJGznVI0mF\nccS/ANbm9Nemd+DMyN9Rv6SNDP4FMBwOJ6G/NuUDq+HfarV8+LykJ/BePQtiMBjQbDYno/vMZDgc\nGvrSgvNePQVrtVpPmNKJiGJD32sapO0Z/FooXtMg7cw5fi0Ur2mQdmbwa6F4TYO0M6d6JKkwjvi1\nULymQdqZwa+F4jUN0s5cx6+F4zUNKokPW5ekwngBlyRpRwa/JBXG4Jekwhj8klQYg1+SCmPwS1Jh\nDH5JKszUwR8Rr4mIL0TE9yPihdu0a0XEyYi4PyKum/Z4kqT9UWXEfy/wKuDTWzWIiLOA9wIt4PnA\nNRHxvArHlCRVNPW9ejLzJLDTTa+OAqcy84Fx25uAq4H7pj2uJKmag57jvxA4vW77wfFnkqQZ2XbE\nHxG3ARds8tXbMvPju9i/N9+RpDmzbfBn5i9U3P9DwMXrti9mddS/qaWlpcn7RqNBo9GoeHhJWiyj\n0YjRaFRpH5XvzhkRJ4A/zMx/3+S7I8B/Ai8Hvgz8K3BNZj5pjt+7c0rS3tV6d86IeFVEnAauAG6J\niOXx58+IiFsAMvNx4FpgCHwR+MhmoS9Jqo/345ekQ8z78UuSdmTwS1JhDH5JKozBL0mFMfglqTAG\nvyQVxuCXpMIY/JJUGINfkgpj8EtSYQx+SSqMwS9JhTH4JakwBr8kFcbgl6TCGPySVBiDX5IKY/BL\nUmEMfkkqjMEvSYUx+CWpMAa/JBXG4Jekwhj82leDwYDMnGxnJoPBYIYVSdrI4Ne+GQwGtNttut0u\nmUlm0u12abfbhr80R47MugAtjmazSafTod/vTz7r9/t0Oh2azeYMK5O0nsGvfRMR9Ho9gEn4dzod\ner0eETHL0iSt41SPJBVm6hF/RLwGWAKeC7w4M+/aot0DwLeA7wOPZebRaY+p+bY2p782vQNnRv6O\n+qU5snYSbq8vVgP/MuAE8MJt2n0JOHcX+0utOnHixKxLmMry8nIC2el0cmVlJVdWVrLT6SSQy8vL\nU+3zsPbFQbAvzrAvzhhn557ye+qpnsw8mZn/tcvmDvX2YDQazbqEqbRaLZaXlyej+7U5/+XlZVqt\n1lT7PKx9cRDsizPsi2rqmONP4JMRcWdEvLGG42mGWq3WE6Z0ImLq0Jd0MLad44+I24ALNvnqbZn5\n8V0e42WZ+ZWIeDpwW0SczMzb91qoJGl/RK67ynKqHUScAP4gtzi5u6HtceA7mfmuTb6rVogkFSoz\n9zSdvl/r+Dc9aEQ8BTgrM78dEU8FXgn86WZt91q4JGk6U8/xR8SrIuI0cAVwS0Qsjz9/RkTcMm52\nAXB7RNwN3AF8IjP/sWrRkqTpVZ7qkSQdLrVeuRsRrYg4GRH3R8R1W7R5z/j7eyLi8jrrq9NOfRER\nvz7ug89HxL9ExE/Pos467Ob3YtzuxRHxeET8cp311WmXfyONiPhcRPxHRIxqLrE2u/gbOS8iBhFx\n97gvfnsGZR64iPhgRDwcEfdu02ZvubnXhf/TvoCzgFPAJcDZwN3A8za0OQbcOn7/EuCzddVX52uX\nffEzwA+P37dK7ot17f4Z+ATwK7Oue4a/F+cAXwAuGm+fN+u6Z9gXS8A71voB+BpwZNa1H0Bf/Bxw\nOXDvFt/vOTfrHPEfBU5l5gOZ+RhwE3D1hjZXATcCZOYdwDkRcX6NNdZlx77IzM9k5jfHm3cAF9Vc\nY11283sB8PvA3wH/V2dxNdtNX/wa8PeZ+SBAZn615hrrspu++ArwtPH7pwFfy8zHa6yxFrm6/P3r\n2zTZc27WGfwXAqfXbT84/mynNosYeLvpi/V+F7j1QCuanR37IiIuZPWP/v3jjxb1xNRufi8uBc6N\niBPjiyJ/s7bq6rWbvrgB+MmI+DJwD9CpqbZ5s+fcrPO2zLv9Y924rHMR/8h3/d8UET8P/A7wsoMr\nZ6Z20xfvBt6SmRmrlwUv6tLf3fTF2cALgZcDTwE+ExGfzcz7D7Sy+u2mL94G3J2ZjYj4cVYvEH1B\nZn77gGubR3vKzTqD/yHg4nXbF7P6L9N2bS4af7ZodtMXjE/o3gC0MnO7/9U7zHbTFy8CbhrfCuI8\noB0Rj2XmzfWUWJvd9MVp4KuZ+T3gexHxaeAFwKIF/2764qXAnwFk5n9HxJeA5wB31lLh/NhzbtY5\n1XMncGlEXBIRPwi8Ftj4h3sz8FsAEXEF8I3MfLjGGuuyY19ExDOBjwK/kZmnZlBjXXbsi8x8dmY+\nKzOfxeo8/5sXMPRhd38j/wD8bEScNb5A8iXAF2uusw676YuTwCsAxnPazwH+p9Yq58Oec7O2EX9m\nPh4R1wJDVs/YfyAz74uIN42/vz4zb42IYxFxCvgu8Ia66qvTbvoC+BPgR4D3j0e6C/ksg132RRF2\n+TdyMiIGwOeBFeCGzFy44N/l78WfAx+KiHtYHcT+UWY+MrOiD0hEfBi4EjhvfNHscVan/KbOTS/g\nkqTC+OhFSSqMwS9JhTH4JakwBr8kFcbgl6TCGPySVBiDX5IKY/BLUmH+H2ePFB3cB54kAAAAAElF\nTkSuQmCC\n", "text": [ "" ] } ], "prompt_number": 6 }, { "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": 7 }, { "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": 8 }, { "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.3414898005
\n", "Number of Parameters: 3
\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": 9 }, { "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+/AAAIABJREFUeJzt3XlwnPd93/H373n2BrBYHIsbIAASvEkdNCnRoiTakmPF\n8dXEaepacWO3nSptMp1Op3GcHlamh5q2k8l06tSNE2eSJrFrO5dsx1YkRZRNiRIlUuJ9gARA3Mfi\nWuy9+zy//rEgRFIECRAL7PV9zYBcYB8svg+Oz/729/wOpbVGCCFEcTDyXYAQQoiVk9AWQogiIqEt\nhBBFREJbCCGKiIS2EEIUEQltIYQoIo71/gJKKRlTKIQQ90BrrW792Ia0tLXWBfH2la98Je81yDnJ\nOZXKm5zT+r4tR7pHhBCiiEhoCyFEESmr0D58+HC+S8g5OafiIOdUHIrhnNSd+k5y8gWU0uv9NYQQ\notQopdD5uhAphBAiNyS0hRCiiEhoCyFEEZHQFkKIIrKm0FZKeZRSbyql3lVKnVdKPZerwoQQQrzf\nmqaxa60TSqkPaa1jSikHcFQpdUhrfTRH9QkhhLjBmrtHtNaxxZsuwARm1vqYQgghbm/Noa2UMpRS\n7wITwCta6/NrL0sIIcTtrHmVP621DdyvlKoGXlBKHdZaH7nxmGeffXbp9uHDh4ti1pEQQmykI0eO\ncOTIkbsel9MZkUqpfw/Etdb/44aPyYxIIYRYpXWZEamUqldKBRZve4GPAO+s5TGFEEIsb63dI83A\nHymlDLJPAP9Xa/3y2ssSQtyr7HrMoIHrzTTDeF+DTRQpWTBKiCKSTFlcG59jaDJMOmNj2ZqMrbE1\n2HY2rG0AbZNtS4HWNigDczG3DUNhKDANhWEonKbC53LQ3lRNU20lLqeZt/MT71mue0RCW4gClUpb\nnLk6wdRcnGTGJm1lQ7ku4KMuUIHDzM2EZq01iXSG0EyUSDSBUgqHqXCYBi6Hoqm2kq3tdbhdEuYb\nSUJbiAJna03v0AxXR+ZIpC0ytqatMYC/wo1S+enesLVmZj7GRGgB01A4HQYuU1Hj97CzM0iVz5WX\nusqBhLYQBci2Nad6JxgKRUikLOoCFTTWVeYtpFdCa000nmZkYh6tbVwOA5fTIFjtZVdXELdr3fcL\nLwsS2kIUkL6RWc4PhIglLZobqqmt9uW7pDXRWhOOJhmdnMdQ4HaaeFwmrfWV9LTXYhqyNt1q2LbG\nNI3bhrY8JQqxQdIZm6Onh5gOJ6is9NLZHizoFvVqKKWorvRQXelZ+phla8bnY5ztv4LToXA5DCrc\nDjoa/WxqDpTMud8L29Yk0hbDE/MMTS2QTltLF5Uzlkbd4UlOWtpCrLPZhThHTw0TT9t0t9Xh9Tjz\nXVLepDM2ofkY8/NRnKbC5TTwuRxsaa+lucC7hVYrlbaYnI0wMDZPNJFZDGSbjK2xbI2BoqY6e1HZ\n6TDe97kPbq6T7hEhNtLIZJi3Lo1jaUXPpnrpIrgNrTXpjM3E9AKRaBKnQ+FwGLhNRU2Vh572Oip9\nLowCDPN0xmJyJsrgxDzRRJqUpbGsbGvZsjQ2Gr/PQ0Nddhjlap6QJLSF2EADY3OcvDyJ0+mgq622\npFqPG8XWmmg8xdhUmEzGxmEqTENhmgqnoXA7Ter8HoK1ldRUenCYRk4mENl29klkZiHO1GyUuYUE\naVtnuy80ZCwbW/NeKFd4qa/x4XY6cjqBSUJbiA0wNhPhjbOjuN1ONrXU5ruckqV1tpthIZpkPpIg\nnkihAQOFMsBUBtefJw2VnRl6PU81oDRY16eL6uwThL4+OWnxa1T6XPgrPfi8LhxGbp4QVuNOob0h\nFyL/9ngfOzrraa2vkum0ouTMLsR59dQwSpn0dDZIy3qdKaVwOkxqq31FP+rmXmxIaDc11HBmYJY3\nzo/jcig8LpMtLTV0tQQkxEXRSmdsXnyrn0TapmdTPYb0WYsNsCGhbZoGm5oDS+/bWtM3scA7V6YW\nrx4b3L+1ifoyfNYUxUdrzZvnRhmYCLOtq1HW6hAbakP6tM8Ozt3xmEzGZmB0hnTGwuM0aAh4eWBr\nE06H/DGIwnJ1ZJa3L43THAxQF5BGhlgfeb8QebfQvlU4mmRobBaXQ1HlcbJ/Vwt+n3udKhTi7mKJ\nNH97vB+320VHS02+yxElLu8XIlfLX+Fm15YmIDsW8qUTQ2DbeN0mu7vraQv65WKP2BBaa948P8bA\nRJgd3Y05W1lPiHtVkC3t5WitGRyfJxZP4nWa9LTXsKWtpiAH3oviNzkb45V3BmmsqyJYW5nvckQZ\nKbqW9nKUUksXNLXWXJtc4NTVKXzXA7y9VgJcrJlla145eY35WJrdPc35LkeImxRVaN9IKUVLg5+W\nBj8Ag1NhTvWFsgHeUUtPW410oYhVuzI8w/GL43S319PcKNdRROEpqu6RldBaMx6KMDsfxec22dUV\npLO5esO+vihOqbTFD9/sw+Fw0NkqsxlFfhXd6JFcud4HHo0mqPA42L+jmaAM0xK3eLd3ggvXZtje\n3YDLWbQvPkUJKZk+7dW6sQ/csm1eOzuKZVn4vU4O7mmj0itbJZWzcDTJC8f7qKmuZO+2lnyXI8SK\nlHRo38g0DHo2BYHsMMIfHR/AUFDv9/DwrlaZ1VZGtNa8dnaE0VCUHd1NMv1cFJWS7h5ZiVgiRf/Q\nDC6noj1Yyf09TZgyFrdkjYYWOHpmhOZg8W/xJUpX2fZpr9ZcOM7I5Bxep8m29lp6OmQt5FJhWTZ/\nd3KAaNJmS0e9/FxFQVu3Pm2lVDvwx0AD2aVqf09r/T/X8pj5FPB7Cfi92THgoQin+y5T4XVw/5YG\nmuur8l2euEeXB6c5eXmS7o56WuQ6hihya2ppK6WagCat9btKqUrgBPBprfWFG44pmpb27dhac210\njlQqhd/r5KGdLVRVyPjdYhBNpHjxrQG8HjdtTYG7f4IQBWLdWtpa63FgfPF2RCl1AWgBLtzxE4uI\noRRdrdkFgtIZmxdPDGEoTX2VmwNyAbNgvXFuhMHJBbZ3Nco1ClFSctanrZTqBF4FdmmtIzd8vKhb\n2suJJdIMjEzjdhq01Vdyf0+jjEIoAKOhBV47M0JjfbUsnSqK1rqP017sGvku8C9vDOzrvvrbzy3d\n3n/wEAcOPpqLL5tXPo+TnZuzKxHOhRN899VefC6TbR019LTX5bm68pNMW7x4vA8Lgx2bm+RCoyg6\nx4/9hLeOHQWy698sZ80tbaWUE/g+8EOt9e/c5v6SbGnfjtaaiekIM/NRqjwO7u9poKlOLmCuJ601\nx84OMxSKsbUziEs2zhAlYN2G/Klsc+aPgGmt9b9a5piyCe0b2bZmcHyOZDJFlcfJQ7ta8MsFzJzq\nHZrhZO8kbY0BAn5vvssRImfWM7QPAT8GTsPS7vNf1lr/6IZjyjK0b5SxbPqGp1FaE6hw8vDuNjyu\nspmMmnMjoTBvnh+jwuuRUSGiJMnkmgIST6bpH57BaSqC1R4O7GzF6ZALmCsxE45z9PQwyjTpapWJ\nT6J0SWgXqEgsxbXRadwOk2C1h/0S4LcVjiZ59Z1BLBSbO+plowtR8iS0i0A0nsoOIXSY1Pnd7N/R\ngrvMu1DGpxd46+IEGRs2d9RhypBKUSYktItMPJGmb2QGtwmVPhf7tzeX1UXM/tFZTl+dwljckEBa\n1qLclO162sXK63Gya3MjkJ2F+dLJIZS28bkdbNtUy6bG6pLrz81YNm+cG2F8JkZlhYctnY35LkmI\ngiQt7SJia83IxDwL0QRep0lNpYv7tzZT4XXmu7R71jcyy/lr08QSGTa11lHpkwWdhJCWdokwlKL9\nhiFu8WSaF94aQCmNx2nSVFvB7u4gzgKfYDIwPsf5/mxQV1V66WyTpVKFWClpaZcIrTUL0STDE/M4\nTHA7Tao8TnZvDlJTld+JJ4lUhlO9E4Tm48RTFpU+D61N1dJXLcQypKVdBpRS+Cs97Kz0LH0skc7w\n6qkRMhkLp8PA7TDxuAw6m6tpDfrXpUVu2ZrhyXmuDs8ST1skUza2hvbmGjrbK3P+9YQoNxLaJczj\ndCzti3ldOmPROxbl7UtTKEPhMBQOE0xlYJrgNA3cTpMav5cKjxOn04HTNDAMRTJtkUymiSYzzIXj\nJNMWqYxNRmusjE3a1liWprrKS1N9AIcsiSpEzklolxmnw6QlWEVL8P0LWWmtsWxNOmMzPp8kMRnD\nsmws20bbGqfTxOk0cTlNKit8VDpMTEO6OITYSBLaYolSCoepcJgGXrf8aghRiOT1qxAF5OiRl7hx\ncIDWmqNHXspjRaLQSGgLUSCOHnmJZz7/GX7rN7+M1hqtNb/1m1/mmc9/RoJbLJHXwEIUiEcef4Kn\nv/gMf/KNry197E++8TWe/uIzPPL4E3msTBQSCW0hCoRSii99Jbs13/XgfvqLz/Clrzwnk4/EEuke\nEUKIIiItbSEKxPU+7OtdIvBei1ta2+I6CW0hCsRrr768FNjXu0kgG9yHDj/JocNP5rE6UShk7REh\nCsjRIy/xyONPLLWqtda89urLEthlRjZBEEKIInKn0JYLkUIIUUQktIUQoohIaAshRBFZc2grpb6h\nlJpQSp3JRUFCCCGWl4uW9h8CT+XgcYQQQtzFmkNba/0TYDYHtQghhLgL6dMWQogisiEzIr/62+/N\n7tp/8BAHDj66EV9WCCGKxvFjP+GtY0eB7F6ry8nJ5BqlVCfwPa31ntvcJ5NrhBBiFWRyjRBrJDvK\niEKRiyF/3wReB7YqpYaUUl9Ye1lCFA7ZUUYUkjX3aWutP5uLQoQoVKW6o4wsTlWcZGlWIe6iFHeU\nuf7q4cZlYK+v5f21P/6uBHcBk9AWogyV6quHciChLcRdlOKOMqX46qFcSGgLcReyo4woJLIJghAr\nUOgX7dIZi0g0STSeJJWxSKctMpaNaShcTgdOh4nTaeJymlT63Dgd5m1fPRR7a7vQf04rdadx2tLS\nFmIFbv2jV0ptaBDYtmZyeoGBkWkGx2YYDy0wNbPAxPQCs+EY8UR6dY8338fga1+n8/6PYDcfpqG2\nikMfneZPvvE1dj3wQT72iY9jGsU1jaNcLq5KS1uIAhSJJblwdTz71jfOlWtTxJPLB7NhKKoq3FR4\n3bicJk6HicM0sWyb9GLLO5W2SKUzLESTWLZNbOoS3vqtN7VK46HL+ILbcJgGLY3VdLXW09VWR1db\nHZ2tddRW+wq2Fb7ctYdifPUge0QKUeBsW3NlcIoT5wY5eW6IS/0T2Lf8bdb4fXS21tLRUktLsJqG\nuioa66qoDVTg87gwjJWFkm1rIrEkc+EYcwtxZsMxJqcXGJuaZ3QyzNjUPNNz0dt+bsDvZVtnI9u6\nG9nW2UBPZwM+j2vN558rNwY3FO/FVekeEaIA2bbm0sAER09c5fWTfUzNRpbuc5gG2zsb2bG5iZ2b\nm9jW1UjA78vJ1zUMhb/Sg7/SQ8cyxySSaYbGZ+kfnn7vbSTEXDjOm6cHePP0AABKQUdzLdu6Gtm1\npYldPS001lUVXUgWE2lpC7HBpmYivPzGJV56/SLjofDSx+sCFRzYu4l9OzvYu721oFqwkG3Fjk2F\nudQ/sfTWNzSNZds3HResqWRXTzO7e1rY3dNMa2NgQ0K8XLpHpKUtxAbQWnPy/BDP/90ZTp4f5Hpb\nqS5QwaEHN3PoA5vZ1tm44i6OfFBK0dJQTUtDNR96aCsAyVSGvqEQF/rGOX9ljLO9Y0zNRjhyvJcj\nx3uBbJfK7i0t7OppZs/WFjqaa9flPMtlaKa0tIVYR6l0hiPHe/mrl08xOJrd4MnhMDh4XxcfeWQH\n921vLbpRGndi25rBsRnOXB7lXO8YZ3tHmVuI33RMVYWbXVua2dXTwp6eFrra63L2PSiHIX8S2kKs\ng0g0yfOvnOYHr55lfiEBQG11BZ/40G4+emgn/kpPnivcGFprRibmOLsY4Gd7RwnN3nyR0+dxsWNz\nE3u2ZrtTtmwK4jDNPFVcGCS0hdgg0XiS5//uDH/10imi8RQAm9vr+fST93Fo32acjvIOI601E9ML\nnL08uhjiYzf16wO4XQ52dDexezHEt3Y24HKWV0+uhLYQ6yyWSPH9V87wFy+eIhJLAnDf9lZ+4af3\nsWdrS1FdBNtoUzMRzvaOcu7KGGcvjzI8cXNeOB0m27oa2d3TzM4t2RCv9LnzVO3GkNAWYp1kLIsf\n/eQCf/b9twhHst0gu3ua+dwnDrBna0ueqytOs+EY53rHsv3iV0YZGJl53zHtTTVs625ke1cj27oa\n6WipKalrAxLaQuSY1pq3zw7yB3/+OsPj2d/v7d2NPP3JA9y3rVVa1jkUjiQ4f2WMM72jXOyb4OrQ\nFJnMzcMMPW4HPZsa2NbVSHd7PZvb62kOVhf0aJw7kdAWIocGRqb5/e+8zrsXhwFoDvr5ws8e5OD9\nXRLWGyCdtugbDi2NFb/YP8FEaOF9x3ncDrpa67Mh3lFPd1s9m1pqcToL/7qChLYQORBLpPjT59/i\ne6+cwdaaCp+Lz37sA/zM4d1lf4Ex32bDMS73T3D52hR9QyH6hkK3nYpvGgYtDdW0Nwdob66lvSlA\ne1MNrU0BPC5nHiq/PQltIdbozVP9fO1bR5majWAoxc8c3s1nf+YDZTN0rxjNL8S5uhjg199GJue4\nXeQpBQ11VXQ019LaGKC53k9T0E9TvZ+GuqoNf1KW0BbiHoVmI/zet4/y+jv9AGzpCPKrTz/O5o5g\nnisT9yKRTDM8McfQ2CzD47MMjs0yND7L2GT4fdPxrzOUor6mkqZgFU311TQF/QRrKqmvqaS+poL6\nQGXOu1wktIVYJcu2+eGPz/NHf/UG8UQaj9vBL37yIT7+od0lNUphIxTDLMV0xmJ8Kszg2Ayjk/OM\nh8KMTYUZnwoTmo28b8XFW1VXeagP3BDki6Fe4/cR8HsJVHmpqvSs+HdHQluIVegfnuZ//ekRLvVP\nAvDQ3k6e+QePEqytzHNlxacUNiZIZywmpxcYD4WZCIUZD4WZmo0yPRthajbCzFxs2Vb6jQyl8Fd5\nCFT5qPF7Cfh9BKq8VFd5l96vrvTir/TgcTv54M4mWTBKiDtJpNJ86wcn+MsXT2HZNrXVFTzzDw7J\nqJA1KIVd350Ok9bGAK2Ngdveb9k28wtxpmYiTM9FCc1GFt+izC3EmAvHmQvHCUcTS7cHRu69njWH\ntlLqKeB3ABP4fa31b631MctFMbxsLBcnzw/xu3/2Y8ZDYZSCn3l8N5//9AEqvKU98269lcOu76Zh\nUFtdQW11xR2Py1gW8wuJpc0n5sIxZsPxpdtz4TgL0QThSIL5SHzZx1lTaCulTOB/AU8CI8BbSqnn\ntdYX1vK45aAUXjaWgrlwjK9/53VefSu7jGhnay2/8rnH2d7dlOfKRKlxmCZ1gQrqAncOd8guebvv\nr3/t9o+zxjoOAFe01gMASqlvAZ8CJLTvohReNhYzrTUvvn6Rb/z5MSKxJG6ng89+/AN8+sm9Zb/C\nXC4ttzEBUFKt7Vy70/dlraHdCgzd8P4w8NAaH7MslMPLxkI1ND7LV//0Vc72jgHwwM52/sVnH6Mp\n6M9zZaWnXDYm2EhrDe0VDT35zC98AYdp4HCYtPfspXv7A3jcTtwuBx6X46bbPq+bSl/2rcLnwut2\nSoiJnEinLb79o5N854WTZDI21VUe/unPP8Lj+3vkd2ydHDr8JF/74+/edO3mS195TgL7No4f+wlv\nHTsKgGUvH61rGvKnlHoYeFZr/dTi+18G7BsvRiqldOfHnlvuIe7KNAwqfK6lIK+q8FBTnR0qU+P3\nUVPty/6/OGSmwusqij/AUtrPrhicuTzKV//01aVlP3/qke184WcPUlUhMxpF4VnPPSLfBnqUUp3A\nKPALwGdvPejZX/kYiWSGZDpDMpleup244XYymSaeTBNLpIjEkkSiKaKxJMl0hnAksbTs5d24nCbB\n2koaaqtoqKuise76/34a66uo8fsKIhDlZePGCEcS/OFfHOPF1y8C0NYY4FeefpzdPbJsqihOa55c\no5T6ad4b8vcHWuvnbrl/TZNr0mmLSDxJJJYkGksyH0ksDY+ZDcfee5vPDpuJJ9N3fDyvx0lbYw1t\nTQHaGgPZ/5tqaAlWb/jqXzLkb/3YtublNy7xh39xjHAkgcNh8PefepCf/+iDRbHKmyhvZTUjMpZI\nMTUTYXJ6gcnpBSZu+H88NM9CNHnbzzMNg/bmGja319PdXkd3e5CutrqS3yGjFA2MTPO73/wx56+M\nA7B3Wwu//NnHaG+qyXNlQqxMWYX23cxH4gyPzzE8Ppv9fyJ7ezwUvu3qX411VfRsamB7dyPbuxvZ\n3B6UllqBiifSfPMHb/FXL5/GtjWBKi//+DMf5PABudAoiouE9grEE2kGRqfpHwrRNzzN1cEQ10an\nSaWtm45zOAw2twfZ3t3Iju4mdvU0U+P35alqAdlupWPv9vN73z5KaDaKUvCxx3bzi586IK+URFGS\n0L5HlmUzND7L5YFJLvaNc7FvgsGx2fcd195Uw95tLezZ2srunmYCEuIbpm8oxO9/93VOX8ou5rCl\nI8i/+Nxj9GxqyHNlQtw7Ce1VutMFwkgsuRTi56+Mc6FvnGQqc9Pnd7TUsHdrKw/ubGfP1la8nsLZ\nEaNUzMxH+ZPnj/Pi6xfRGqoq3Dz9iQM89dhOWTpVFD0J7VVY7Zog6YxF78AkZy6PcvryCBevTpBM\nvxfiDtNg55ZmHtzZzoO72ulqrZP+1TVIJNP89cun+c4LJ0kkM5iGwccXd5GprJCuEFEaJLRXYa2T\nXtJpi0sDE5y6OMzJ80NcHpi86QJnbbWPB3ZkA/yBHe2yXdUKJVMZfvjjc3z3hXeYW8iugPbQ3k6+\n+HMHl10yU4hiJaG9SjcGN6xtTZBwJMG7iwF+8twgM/OxpfuUgp5NDezb1c6+XR30dDbIS/tbpNMW\nL7x2gW//8MTS925rZwP/6NMPcd/2tjxXJ8T6WM8ZkeIu/JUeHvvAFh77wBa01lwbneHEuUFOnh/i\n3JUxLg9Mcnlgkm/+4ARVFW4e2NHOvt0dPLizvaxHpcTiKf72tQv89cunmZqNALC5vZ7PfWI/+/ds\nki4mUbakpX2LjVwTJJFMc/rSCCfODXHi3CDjofBN92/pCLJvdwf7drazrasR0yz9Vvjk9ALPv3KG\nF46eJ57Izm7d1FLL5z6xX3aQEWVDWtqrsJFrgnjcTg7s7eTA3k601oxOzvP22UFOnh/kzOVRrgxO\ncWVwiv/3Nyeo9Ll5YEcb+3Z1sG9XBzXVpdMKt23N6csj/OjH53n9nb6lTVT3bG3h00/ex/7dmzAM\nCWshQFrat1UIa4IkUmnOXh7lxLkh3j47yNjU/E33d7fXs29XO3u3trJjcxMed/ENKxyfCvPK8cu8\n9PpFJqYXgOxyAo/t38KnntjLlo5gnisUIj/kQmQJGJ2c58S5QU6cG+TMpdGbhhWahkFPZ5DdPS3s\n2drCjs1N+DyuPFa7vJGJOV5/p4+jJ69ydTC09PFgbSVPHtzORw/toL5Gdj0X5U1Cu8Sk0hnO9o7x\nzvkhzvaOcnUwtNSlAGAYis3t9fRsamDLpiBbOxtob6rJS594OJLg/JWx7OiZ80M39dt7F7uHPvLB\n7ezd1ipdIEIsktAucdF4kvNXxjnbO8rZ3lF6r01h37LzRWbmCjsf/CAdLbXZJWkbqxnvP80nPvnJ\nnIS51prZcIzB0VkGx2a4MjjFxb4JRidv7tap9LnZt6uDR/dt5sFd7biccllFiFtJaJeZeCJN77VJ\neq9NcuXaFMdfP8Kll3+Xqs4PUrvj4wDMXPg+CwOv07T/C7RsuZ+aah+11RXU+LO7Ank8zqWt4K6z\nbU06YxGNp4gtrnEemo0SmosQmokSS6TeV4vLadKzqYH7tmen9ctYdCHuTkK7zGmt+Y//7tf49v/9\nOoc++vMsRJOcOvo8jdsP4+n6qZwNo6vwudjUXEtHSy2bWmrZ3t1IV1ud7G4uxCrJkL8yp5Ti3/+n\n/4bLab5vlmfGspkLx5mZjzI7H2MmHGN+IU4ytbgdXCoDWqOUQilwOkx8XhcVixsw1wUqqKupoD5Q\nSVWFe0PHURfCKB8hNpqEdplzOrJ7agZri2vExmoX9hKiVEhol4HlZnkCRbvz+yOPP8HTX3xm6Tzg\nvZmrjzz+RB4rE2J9SWiXgVLc+V0ptXQuuVjYS4hiIaFdBg4dfpKv/fF3b+r//dJXnivawBainElo\nl4lbw1kpVdSBXYpdPkKshIS2KEql2OUjxErc8zhtpdTPA88C24H9WuuTyxwn47TFupAhf6JUrcvk\nGqXUdsAG/g/wryW0hRAiN9Zlco3W+iIgfYdCCLGBZBEIIYQoIndsaSulXgSabnPXb2itv7fSL/LV\n337vQtH+g4c4cPDRFRcohBDl4Pixn/DWsaMAWPby3dZrXjBKKfUK0qcthBA5c6c+7Vx1j0jHthBC\nbIB7Dm2l1N9TSg0BDwM/UEr9MHdlCSGEuJ21jB75S+Avc1iLEEKIu5DRI0IIUUQktIUQoohIaAsh\nRBGR0BZCiCIioS2EEEVEQlvk3NEjL3HjpC2tNUePvJTHioQoHRLaIqeub7j7W7/5ZbTWS5sVPPP5\nz0hwC5EDsgmCyCnZcFeI9SWhLXJKNtwVYn1J94gQQhQRaWmLnJINd4VYXxLaIqdkw10h1tea19O+\n6xeQ9bTLjmy4u3601iTTFhOhBSKxJIZSmKbCUGAohWFkb2vAXPz+W1qjAFuDbevs/1ov3tY4HQ5q\nq734Kz04HWZez09krcvGvisloS3EvdFaM7eQYGxyHtNUuBwGLqdJhctBd2uAYE0FDnNtl6VsWxOJ\npxiZCjM5GyWRsklbNhlLk7Y0bpdJa2MAj0telG8kCW0hikQkluTa2BxOU+FxGrTUVbKjsx6Xc+Nb\nwPbik8aF/ikW4mmSGZt0RuNwmHS0BHA7JcjXi4S2EAVsaibC5EwEr8skGPDy4NamvIT0SmitmQkn\nOHN1gmgyQzJlYwObmmvxeZ35Lq9kSGgLUWAi8SQDIzN4nSZdzdXs6g5iFOnImmg8zemrE8wsJEmk\nMqAUna3UF2GhAAANvUlEQVS1eFwS4vdKQluIAqC1ZmB0jlQySWOtj4d2tq65T7oQReNp3rk8xlwk\nRSJj4/O4aGusxizBc10vdwpt6ZQSYp2lMha916bwOk0e2NpAW9Cf75LWVYXXyaH7OoDsE9VoaIFz\n/SGiiQwZGza11FDhdeW5yuIloS3EOokn0vQNT1Ptc/CJg9143eXXXaCUojXop3XxiSqZynDy0hj9\nU/PEUxb1tVUEa3wy6WoVJLSFyLGFaJJro9M0Bnz87KM9OBzSLXCd2+Xg4J52IDs65dLANP3D08RT\nGdweJx3NNZiGfL/uRPq0hciReCLN1aEQrXU+Du5pk/BZpbHpCKevZEelOBxOOltqyrYfXPq0hVhH\nqbRF77VJ6v0ePvP41rINmrVqrqukua4SgNBcjJOXxogkLQzDpKuttiQv2t4LaWkLcY8s26b3Wgif\ny+DD+zpxF+jY6mI3E45z4uIYC4kMKIPutjqcJd7ltC4tbaXUfwc+DqSAq8AXtNbz916mEMVjaHye\neDzBhx/ooKbKk+9ySlqt38tHDnQDMBdJ8Nb5URbiaUynk64y7EK555a2UuojwMtaa1sp9V8BtNa/\nfpvjpKUtSsZsOM7Y5DwP9ATZ0lab73LK2tRiF0o0mcHn9dDaWF20E5RutS4tba31ize8+ybwc/f6\nWEIUOsuyudg/QXuwks8c3ipD1ApAMODjow9tBmBwYp6z/SFiiQw1gUoaaytK9meUqwuRXwS+maPH\nEqKgDE/Mk4wn+fjBzfg85TfWuhh0NFbT0ViN1prLgzNcGQoRTWVoaQhQ4/fmu7ycumNoK6VeBJpu\nc9dvaK2/t3jMvwVSWus/W+5xvvrb7y2Gv//gIQ4cfPTeqhViA8Xjaa4OTfFATwNbOzryXY5YAaUU\n2zbVsW1THbatOdU7wZWBSZJpm672Wnyewp2JefzYT3jr2FEALHv5bus1jR5RSv0S8E+BJ7TWiWWO\nkT5tUVS01vQPz+BxKj68b5OMty4BqbTFm+dGmAonsDX0dAQLetLTeo0eeQr4N8DjywW2EMVmNhxn\nfGqeQ3taaKqrync5IkdcTpNH78++WlqIJTl2doRwLI3b5aSrrbao+r/XMnqkF3ABM4sfOqa1/ue3\nOU5a2qLgWZZN77UpWmt9PLynLd/liA0yElrg1OVJoqkMgSovLQ3V+S4JWL/RIz1rK0uIwjAeCrOw\nEOOn9ndS5XPnuxyxgVrrq2itr8JevIB5aWCSRNqirSlAoKowL2DKjEhRtlIZi8v9E+zsqGPPloZ8\nlyMKRMayOXFxnJHpCBlLs7m9Ds8Gr9Aoa48IcYuB0Vl0JsPfe3RrwW7tJfLDYRo8tKsFgHgyw7Gz\nw/QvJFFK0bMpmPcZmBLaoqxEY0n6hkIc2NFMd2tNvssRBc7rdvDhfZ0ATM/HefP8aHYGpsdJZ2td\nXmqS7hFRFrTW9F4LUeUx+ZAM4xNrYGvN8OQC716ZJJGyqAtU0FSf25FG0j0iylpoLsbE1ByPPdBO\nU01lvssRRc5Qio5GPx2Nfixbc3EgxJUNvIApLW1RsjKWzYWr43Q3V3NgZ0u+yxElLp2xeOfyBKPT\nUVIZm+62Orz3uOyBtLRF2RkcnSWZTPGpQ1vKcm9GsfGcDnOpcZBMZXjz/CgjE0kSKZvO1hoqczSc\nVEJblJRINHuh8QM7mulpkwuNIj/cLgePLc7ATKUtTlwa4+q17GbGLQ3V1Fb77vmxJbRFSbBsm0v9\nU9RVufj7T2yXC42iYLicJgd3Z2fZWrbm7NVJ+ganSKQtqqt8tASrVjWNXvq0RdG7NjpDIpHip/Z3\nUekr3FXchLiR1pr+0TkuDc0QT1lYNnS11uL1OMuzT9u2NemMzXwkzlw4Tsay0YCpwFBgKANlAFpj\nGNnvi1Jw43OYZQMaNNmWnK2z32iP20nA78Vf4cZhGkW12EwpmZ6LMR6a5wNbm+huDeS7HCFWRSlF\nd2vN0nyBWDLDiQujjE4myVh3+LxibWlrrUmkMkyEIkTjSUxD4TAUpqkwDYXTNPC4TBpqKmiqq8Tn\ncWIaak0Bq7XGsjXz0STj0xFCczGSaYu0rclkbDJ29ski4PfSUFspu0evk2g8Sd/QNF3N1RzY0SxP\nmqLkaK0xDOO2Le2CD23b1oSjCcZDC2hb4zQVDoeB01RUeJx0twQI1lQUTECmMxbDk2EGxuaJpy3S\nGZt0RuPxOGhtCMiU6TVIpjNc7p+kMeDlsfs78j6dWIj1pJQq7NC2tSYcSTI2FQatcTqy4ew2FQ0B\nH5vb66go0q2etNaE5uKc658kmrJIp23SlqYl6Ke6yiMtxbtIZywu90/ir3Dy4Qc75YlPlIWCCW2t\nNZFYipHJedAah2HgcCjcDoNgwEdPRx1el1nyQZZKW1wYCDE2EyGZtkmmLVqC1QT83pI/95VKJDP0\nXpvE73PxoQc7ZLy1KCt5De3vHLmEwzBwOsBhmtT53WztqKfS6yyZLe/XKmPZnO+fYnQ6QiKd7VJp\na6ymqsJddiEejiQYHJulzu/msfs6pGUtylJeQ9uybQnnVcpYNmeuTjA2HSOZsbAsTVtzTcku0q+1\nZmh8nlgsQXOdj4d3tS2N6hGiHOU1tNf7a5SDZMriVO84k/Nxkuns8MNNLTVUeIt7XHI0nqR/eAav\ny2Rvd5DOFhm6JwRIaJecWCLNO5fHmV1Iksxkx6BvaqnB5yn8EI8n0/QPT+MwDYLVHh7e1Vowo3+E\nKBQS2iUumkhzqjcb4ikr2ydeF6igoa4y711TWmumZiJMzUXxOA0CFW7272jB6y7ZuV1CrJmEdpnJ\nWDb9o3P0j81lJwBZmowNPo+TlgY/buf6BWYynWFkIkw8mcJlGridBp1N1WzrqJN+aiFWSEJbYNma\nqdkol4emiSUzWLYmk4GMbWPp7OLuPo+TqkoPFV7X0gzS6781mvdmhS5EE8wvJIgn09llAQyF08iO\nra9wOdjeWUcwUCEhLcQ9ktAWd2TZmmQqw1wkQWg+xvxCcmm9Fk32H6XANMA0DGr9XoI1FVRXuHE4\njLx3wQhRanIe2kqp/wh8kuzf9DTwS1rrodscJ6EthBCrtB6hXaW1Xli8/avAfVrrf3Kb4yS0hRBi\nlZYL7XseZ3U9sBdVAqF7fSwhhBArs6YhBEqp/wz8IhADHs5JRUIIIZZ1x+4RpdSLQNNt7voNrfX3\nbjju14FtWusv3OYxpHtECCFWabnukTu2tLXWH1nh4/8Z8DfL3fnss88u3T58+DCHDx9e4cMKIUR5\nOHLkCEeOHLnrcWu5ENmjte5dvP2rwAGt9S/e5jhpaQshxCrdU0v7Lp5TSm0DLOAq8MtreCwhhBAr\nIJNrhBCiAOV8yJ8QQoiNJ6EthBBFREJbCCGKiIS2EEIUEQltIYQoImUV2isZuF5s5JyKg5xTcSiG\nc5LQLnJyTsVBzqk4FMM5lVVoCyFEsZPQFkKIIrIhMyLX9QsIIUSJyssekUIIIXJHukeEEKKISGgL\nIUQRKbnQVko9pZS6qJTqVUp9aZlj/ufi/aeUUg9sdI2rdbdzUkp9bvFcTiulXlNK7c1Hnauxkp/T\n4nH7lVIZpdTPbmR992qFv3+HlVLvKKXOKqWObHCJq7aC3796pdSPlFLvLp7TL+WhzBVTSn1DKTWh\nlDpzh2MKNyO01iXzBpjAFaATcALvAjtuOeZjwN8s3n4IeCPfdefgnA4C1Yu3nyqFc7rhuL8Dvg/8\nXL7rztHPKgCcA9oW36/Pd905OKdngeeunw8wDTjyXfsdzulR4AHgzDL3F3RGlFpL+wBwRWs9oLVO\nA98CPnXLMZ8E/ghAa/0mEFBKNW5smaty13PSWh/TWs8vvvsm0LbBNa7WSn5OAL8KfBeY2sji1mAl\n5/UPgT/XWg8DaK1DG1zjaq3knMYA/+JtPzCttc5sYI2rorX+CTB7h0MKOiNKLbRbgaEb3h9e/Njd\njinkkFvJOd3oH3OH/ToLxF3PSSnVSjYc/vfih4phmNNKflY9QK1S6hWl1NtKqfdt0VdgVnJOXwd2\nKaVGgVPAv9yg2tZLQWfEWrYbK0Qr/cO+dexjIQfCimtTSn0I+CLwyPqVkxMrOaffAX5da62VUor3\n/8wK0UrOywk8CDwB+IBjSqk39OJ+qwVoJef0G8C7WuvDSqnNwItKqfu01gvrXNt6KtiMKLXQHgHa\nb3i/neyz5J2OaVv8WKFayTmxePHx68BTWus7vfQrBCs5p33At7J5TT3w00qptNb6+Y0p8Z6s5LyG\ngJDWOg7ElVI/Bu4DCjW0V3JOHwT+M4DW+qpSqh/YBry9IRXmXkFnRKl1j7wN9CilOpVSLuAXgFv/\nyJ8HPg+glHoYmNNaT2xsmaty13NSSnUAfwE8rbW+kocaV+uu56S17tZad2mtu8j2a/9ygQc2rOz3\n76+BQ0opUynlI3uh6/wG17kaKzmni8CTAIt9v9uAvg2tMrcKOiNKqqWttc4opX4FeIHsVe8/0Fpf\nUEr9s8X7/4/W+m+UUh9TSl0BosAX8ljyXa3knID/ANQA/3uxZZrWWh/IV813s8JzKjor/P27qJT6\nEXAasIGva60LNrRX+LP6L8AfKqVOkW0I/prWeiZvRd+FUuqbwONAvVJqCPgK2W6rosgImcYuhBBF\npNS6R4QQoqRJaAshRBGR0BZCiCIioS2EEEVEQlsIIYqIhLYQQhQRCW0hhCgiEtpCCFFE/j/bG9/z\nlI7RoAAAAABJRU5ErkJggg==\n", "text": [ "" ] } ], "prompt_number": 11 }, { "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": {} } ] }