{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib inline"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n====================================================================\nProbabilistic predictions with Gaussian process classification (GPC)\n====================================================================\n\nThis example illustrates the predicted probability of GPC for an RBF kernel\nwith different choices of the hyperparameters. The first figure shows the\npredicted probability of GPC with arbitrarily chosen hyperparameters and with\nthe hyperparameters corresponding to the maximum log-marginal-likelihood (LML).\n\nWhile the hyperparameters chosen by optimizing LML have a considerable larger\nLML, they perform slightly worse according to the log-loss on test data. The\nfigure shows that this is because they exhibit a steep change of the class\nprobabilities at the class boundaries (which is good) but have predicted\nprobabilities close to 0.5 far away from the class boundaries (which is bad)\nThis undesirable effect is caused by the Laplace approximation used\ninternally by GPC.\n\nThe second figure shows the log-marginal-likelihood for different choices of\nthe kernel's hyperparameters, highlighting the two choices of the\nhyperparameters used in the first figure by black dots.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "print(__doc__)\n\n# Authors: Jan Hendrik Metzen <jhm@informatik.uni-bremen.de>\n#\n# License: BSD 3 clause\n\nimport numpy as np\n\nfrom matplotlib import pyplot as plt\n\nfrom sklearn.metrics.classification import accuracy_score, log_loss\nfrom sklearn.gaussian_process import GaussianProcessClassifier\nfrom sklearn.gaussian_process.kernels import RBF\n\n\n# Generate data\ntrain_size = 50\nrng = np.random.RandomState(0)\nX = rng.uniform(0, 5, 100)[:, np.newaxis]\ny = np.array(X[:, 0] > 2.5, dtype=int)\n\n# Specify Gaussian Processes with fixed and optimized hyperparameters\ngp_fix = GaussianProcessClassifier(kernel=1.0 * RBF(length_scale=1.0),\n                                   optimizer=None)\ngp_fix.fit(X[:train_size], y[:train_size])\n\ngp_opt = GaussianProcessClassifier(kernel=1.0 * RBF(length_scale=1.0))\ngp_opt.fit(X[:train_size], y[:train_size])\n\nprint(\"Log Marginal Likelihood (initial): %.3f\"\n      % gp_fix.log_marginal_likelihood(gp_fix.kernel_.theta))\nprint(\"Log Marginal Likelihood (optimized): %.3f\"\n      % gp_opt.log_marginal_likelihood(gp_opt.kernel_.theta))\n\nprint(\"Accuracy: %.3f (initial) %.3f (optimized)\"\n      % (accuracy_score(y[:train_size], gp_fix.predict(X[:train_size])),\n         accuracy_score(y[:train_size], gp_opt.predict(X[:train_size]))))\nprint(\"Log-loss: %.3f (initial) %.3f (optimized)\"\n      % (log_loss(y[:train_size], gp_fix.predict_proba(X[:train_size])[:, 1]),\n         log_loss(y[:train_size], gp_opt.predict_proba(X[:train_size])[:, 1])))\n\n\n# Plot posteriors\nplt.figure(0)\nplt.scatter(X[:train_size, 0], y[:train_size], c='k', label=\"Train data\",\n            edgecolors=(0, 0, 0))\nplt.scatter(X[train_size:, 0], y[train_size:], c='g', label=\"Test data\",\n            edgecolors=(0, 0, 0))\nX_ = np.linspace(0, 5, 100)\nplt.plot(X_, gp_fix.predict_proba(X_[:, np.newaxis])[:, 1], 'r',\n         label=\"Initial kernel: %s\" % gp_fix.kernel_)\nplt.plot(X_, gp_opt.predict_proba(X_[:, np.newaxis])[:, 1], 'b',\n         label=\"Optimized kernel: %s\" % gp_opt.kernel_)\nplt.xlabel(\"Feature\")\nplt.ylabel(\"Class 1 probability\")\nplt.xlim(0, 5)\nplt.ylim(-0.25, 1.5)\nplt.legend(loc=\"best\")\n\n# Plot LML landscape\nplt.figure(1)\ntheta0 = np.logspace(0, 8, 30)\ntheta1 = np.logspace(-1, 1, 29)\nTheta0, Theta1 = np.meshgrid(theta0, theta1)\nLML = [[gp_opt.log_marginal_likelihood(np.log([Theta0[i, j], Theta1[i, j]]))\n        for i in range(Theta0.shape[0])] for j in range(Theta0.shape[1])]\nLML = np.array(LML).T\nplt.plot(np.exp(gp_fix.kernel_.theta)[0], np.exp(gp_fix.kernel_.theta)[1],\n         'ko', zorder=10)\nplt.plot(np.exp(gp_opt.kernel_.theta)[0], np.exp(gp_opt.kernel_.theta)[1],\n         'ko', zorder=10)\nplt.pcolor(Theta0, Theta1, LML)\nplt.xscale(\"log\")\nplt.yscale(\"log\")\nplt.colorbar()\nplt.xlabel(\"Magnitude\")\nplt.ylabel(\"Length-scale\")\nplt.title(\"Log-marginal-likelihood\")\n\nplt.show()"
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.6.7"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}