# gpR: Gaussian processes in R
#
# Copyright (C) 2015 - 2016 Simon Dirmeier
#
# This file is part of gpR.
#
# gpR is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# gpR is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with gpR. If not, see <http://www.gnu.org/licenses/>.
#' Sample function values from a Gaussian process
#'
#' @noRd
#' @importFrom stats rnorm
#' @param x a vector of real values
#' @param m a vector of mean values, e.g. a vector of zeroes
#' @param K a symmetric positive semi-definite matrix of size \code{length(x)} x \code{length(x)}, e.g. calculated using \code{convariance.function}
#' @return a vector of normal-distributed values
#' @examples
#' sample.from.gp()
sample.from.gp <- function(x=seq(-10, 10, .25), m=NULL, K=NULL)
{
xlen <- length(x)
if (is.null(m)) m <- rep(0)
if (is.null(K))
K <- covariance.function(x, x,
list(var=1, inv.scale=2, gamma=2,
noise=0, kernel="gamma.exp"))
# small pseudo-count for numerical stability
C <- 1e-5 * diag(xlen)
# standard transformation of Gaussian random variables: X <- \mu + \sigma Z
# the cholesky decomposition of a matrix K is its square root
# the square root of the covariance yields the standard deviation
# K has to be symmetric and semi-definite.
# Since this is a requirement for a kernel-function,
# we are safe to use chol().
# thus the 'sampling' from a Gaussian process is just generating
# xlen random variables and transforming them accordingly
y <- m + t(chol(K + C)) %*% stats::rnorm(xlen)
as.vector(y)
}
#' Sample a Gaussian random variable
#'
#' @noRd
#' @param m the mean of the Gausian
#' @param var the variance of the Gaussian
#' @return a sample from a Gaussian with mean \code{m} and variance \code{var}
.sample.gaussian <- function(m, var) m + var * rnorm(1)
#' Sample class label probabilities from a logistic mapping of a Gaussian process
#'
#' @noRd
#' @param x a vector of real values
#' @param m a vector of mean values, e.g. a vector of zeroes
#' @param K a symmetric positive semi-definite matrix of size \code{length(x)} x \code{length(x)}, e.g. calculated using \code{convariance.function}
#' @return a vector of class mapping probabilities
#' @examples
#' x=seq(-10, 10, .25)
#' m=rep(0, length(x))
#' K=covariance.function(x, x, list(var=1, inv.scale=2, gamma=2, noise=0, kernel="gamma.exp"))
#' .sample.from.sigmoid.gp(x, m, K)
.sample.from.sigmoid.gp <-
function(x=seq(-10, 10, .25),m=NULL, K=NULL)
{
xlen <- length(x)
if (is.null(m)) m <- rep(0, xlen)
if (is.null(K))
K <- covariance.function(
x, x, list(var=1, inv.scale=2, gamma=2, noise=0, kernel="gamma.exp"))
y.latent <- sample.from.gp(x, m, K)
y <- .sigmoid(y.latent)
as.vector(y)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.