gpr_model_selection0.ipynb

{
  "metadata": {
    "name": "Simple model selection: grid search for GPR"
  }, 
  "nbformat": 3, 
  "nbformat_minor": 0, 
  "worksheets": [
    {
      "cells": [
        {
          "cell_type": "heading", 
          "level": 1, 
          "metadata": {}, 
          "source": [
            "Simple model selection: grid search for GPR"
          ]
        }, 
        {
          "cell_type": "markdown", 
          "metadata": {}, 
          "source": [
            "\n\n", 
            "Run simple model selection (grid search over hyperparameters' space) of\nGaussian Process Regression (GPR) on a simple 1D example."
          ]
        }, 
        {
          "cell_type": "code", 
          "collapsed": false, 
          "input": [
            "__docformat__ = 'restructuredtext'\n\nimport numpy as np\nfrom mvpa2.suite import *\nimport pylab as pl\n\n# Generate train and test dataset:\ntrain_size = 40\ntest_size = 100\nF = 1\ndataset = data_generators.sin_modulated(train_size, F)\ndataset_test = data_generators.sin_modulated(test_size, F, flat=True)\n\nprint \"Looking for better hyperparameters: grid search\"\n\n# definition of the search grid:\nsigma_noise_steps = np.linspace(0.1, 0.5, num=20)\nlength_scale_steps = np.linspace(0.05, 0.6, num=20)\n\n# Evaluation of log maringal likelohood spanning the hyperparameters' grid:\nlml = np.zeros((len(sigma_noise_steps), len(length_scale_steps)))\nlml_best = -np.inf\nlength_scale_best = 0.0\nsigma_noise_best = 0.0\ni = 0\nfor x in sigma_noise_steps:\n    j = 0\n    for y in length_scale_steps:\n        kse = SquaredExponentialKernel(length_scale=y)\n        g = GPR(kse, sigma_noise=x)\n        g.ca.enable(\"log_marginal_likelihood\")\n        g.train(dataset)\n        lml[i, j] = g.ca.log_marginal_likelihood\n        if lml[i, j] > lml_best:\n            lml_best = lml[i, j]\n            length_scale_best = y\n            sigma_noise_best = x\n            # print x,y,lml_best\n            pass\n        j += 1\n        pass\n    i += 1\n    pass\n\n# Log marginal likelihood contour plot:\npl.figure()\nX = np.repeat(sigma_noise_steps[:, np.newaxis], sigma_noise_steps.size,\n             axis=1)\nY = np.repeat(length_scale_steps[np.newaxis, :], length_scale_steps.size,\n             axis=0)\nstep = (lml.max()-lml.min())/30\npl.contour(X, Y, lml, np.arange(lml.min(), lml.max()+step, step),\n              colors='k')\npl.plot([sigma_noise_best], [length_scale_best], \"k+\",\n           markeredgewidth=2, markersize=8)\npl.xlabel(\"noise standard deviation\")\npl.ylabel(\"characteristic length_scale\")\npl.title(\"log marginal likelihood\")\npl.axis(\"tight\")\nprint \"lml_best\", lml_best\nprint \"sigma_noise_best\", sigma_noise_best\nprint \"length_scale_best\", length_scale_best\nprint \"number of expected upcrossing on the unitary intervale:\", \\\n      1.0/(2*np.pi*length_scale_best)\n\n\n# TODO: reincarnate by providing a function within gpr.py\n#\n# Plot predicted values using best hyperparameters:\n# pl.figure()\n# compute_prediction(1.0, length_scale_best, sigma_noise_best, True, dataset,\n#                    dataset_test.samples, dataset_test.targets, F, True)"
          ], 
          "language": "python", 
          "metadata": {}, 
          "outputs": []
        }
      ], 
      "metadata": {}
    }
  ]
}