Transcript Gaussian Process Regression for Dummies

Gaussian Process Regression for
Dummies
Greg Cox
Rich Shiffrin
Overview
• We often wish to describe data by a function.
• One-dimensional example:
– Forgetting as function of delay from study to test
• Multi-dimensional example:
– A hand movement trajectory varies in time and in
several spatial dimensions
• Traditionally one characterizes such a relation
with a parameterized function class:
– E.g. Forgetting might be characterized by a one
parameter exponential function or a two parameter
power function.
• Gaussian Process Regression (GPR) provides a
different way of characterizing functions that
does not require committing to a particular
function class, but instead to the relation that
different points on the function have to each
other.
• GPR can be used to characterize parameterized
functions as a special case, but offers much more
flexibility.
• GPR offers powerful tools, and there are
programs to apply it and an excellent book
(both freely downloadable) by Rasmussen and
Williams (Gaussian Processes for Machine
Learning, the MIT Press, 2006, ISBN
026218253X).
• However, most psychologists (and others) find
it hard to understand conceptually what the
approach does: The various articles, texts, and
lectures inhibit understanding by ‘explaining’
the method with matrix algebra.
• The ‘Gaussian’ part of GPR refers to the fact that
all uncertainty, about any variables, or
combinations of variables, is characterized by
Gaussian distributions. For a one-dimensional
distribution the Gaussian is characterized by a
mean and variance. A multi-dimensional
distribution is characterized by a multidimensional mean, and a multi-dimensional
covariance. We can visualize the multidimensional Gaussian as an ellipse. E.g. in three
dimensions this would be a ‘cigar’ shape.
• It is easiest to explain and illustrate GPR by
using GPR to characterize a simple
parameterized function.
• We will do so for the simplest function
available: a linear function Y = aX + b + ε
a is the slope, b the intercept and ε Gaussian error
• In turn we illustrate the cases in which only
the intercept b varies, only the slope a varies,
or both vary.
• When only the intercept varies, we assume
that a = 1.0 and ε = 1.0.
Linear regression
y=x
y = x +e
• The Bayesian approach is traditional: the
parameter b has a prior distribution, and as
data is collected (pairs-- Xi,Yi) we use Bayes
Theorem to form a posterior distribution for b.
• This approach is likely familiar and illustrated in
two ways: On the left of the next figure, at the
top, is the parameterized function. Just below we
show the prior probability distribution of b
(Gaussian) with a mean 0 and standard deviation
ε = 1.0.
• Below the prior we have the posterior
distribution of b after observing the single data
point Y = 2 given X = 3, derived from Bayes
Theorem: The posterior is proportional to the
prior times likelihood:
pd(b|2,3) ≈ p0(b)p(2,3|b)
The posterior mean has dropped to -1/2,
because we observed Y less than X, and the
variance has dropped, because the data
produces an improved estimate for b.
All the figures in the left column give the
parameter space representation.
Bayesian linear regression
• The second column shows the same idea in a
function space representation: In this case imagine
the probability density rising out of the figure, with a
darker region depicting higher probability.
• We see a Gaussian spread of linear functions with
slope 1.0. At X = 0, these are centered at Y = 0, as
specified by the prior.
• Below is the posterior: It has a lower variance
Gaussian spread of linear functions. At X = 0 these are
now centered at Y = -1/2.
• We have marked with vertical lines the distributions
at the (arbitrarily chosen) points X=2, X=4, X= 7, and
X=9, because it helps us move to data space.
• Note that the first column shows the distributions of
the parameter defining the function, and the second
column depicts the distributions of the functions
themselves. These are duals of each other.
• GPR uses yet another representation, one giving the
joint distribution of selected Y points on the function.
We have arbitrarily chosen to show this joint
distribution at two pairs of Y values, one for the Y
values at X=2 and X=9, and the other for Y values at
X=4 and X=7.
• (The technical treatises on GPR prove that not only this
distribution, but all joint distributions of Y values (for N
X values) are Gaussian, with some mean and
covariance structure.)
• The third column shows these two joint distributions,
one for Y at X=2 and Y at X=9, and the other for Y at
X=4 and Y at X = 7. We have placed the larger Y value
on the vertical axis. Again imagine that probability
density rises from the figure, and darker shading
depicts higher probability.
• These are ‘predictive’ distributions. For each pair of Y
points, the top panel gives what two Y values are
expected based on the prior, below is the same for the
posterior.
• This representation will not be familiar to viewers of
this lecture, so some detailed explanation will aid
understanding.
• The vertical lines in the second column show the Y distributions at
each of the four X values. The mean values at these points, 2, 4, 7,
and 9 for the prior, and 1.5, 3.5, 6.5, and 8.5 for the posterior, must
exactly match the center of the corresponding ellipse in the third
column. E.g. 2,9 has a posterior ellipse center at 1.5,8.5.
• But note that the distributions at various X values in the second
column are not independent: The Y values are constrained to lie on
a line of slope 1.0. Thus a high value for one Y must go with a high
value for other Y’s, and low values must also go together. This is
seen in the third column as distributions lying along a diagonal with
slope 1. Without observation noise this would be a uni-dimensional
distribution along the diagonal. However the two Y values have
independent Gaussian noise added to them; this noise reduces the
correlation and is seen as a spread of probability orthogonal to the
diagonal, thereby producing an ellipse. Thus variation along the
major diagonal reflects the uncertainty about b. Variation
orthogonally represents observation noise.
• Each of the ellipses in the third column depict the joint Gaussian
distribution of Y at two chosen X values. These distributions are not
equivalent to the parameter space and function space
representations, because they only specify the relation between
two X,Y points, for two different sets of points.
• However when we add the fact that the same shaped distribution
applies for any two Y values, with a center at X1-1/2 and X2-1/2,
then all three representations become equivalent.
• Choosing a random point from one of the GPR distribution(s) in the
right column is equivalent to choosing a random point from the b
distribution in column 1, using that intercept to predict mean Y
values at the specified X values, and then adding independent
samples of Gaussian noise to those two points. It is also equivalent
to choosing random points from the distributions in column 2, at
the specified X values (a pair of vertical lines).
• The third column depicts the actual joint
distribution. There is no room on the figure to
give more columns so the next figure gives what
would have been the fourth, fifth and sixth
columns. These give the mean and covariance
matrices corresponding to ellipses in the third
column in the previous figure. The first for the
pair of points---, the second for the pair of points
---, and the third for all four points (not
represented in the prior figure because the
ellipse would be four dimensional).
• These matrices are the usual way that GPR is
represented.
• In this example the shape of the GPR
distribution happens to be identical, save for
the center, at every pair of points. This is not a
general result, and in most cases the center
and shape of the distribution will change for
each different pair of points selected. (See the
next example in which the slope rather than
the intercept varies.) Nonetheless, knowing
the distributions for each pair of points does
provide a compete and equivalent description
to the Bayesian regression.
Example 2
• The next figure illustrates the same situation
for the case with slope, a, varying (with
intercept, b, fixed at 0). For the same
observation as in the first example, Y=2 at
X=3, the first column shows the prior and
posterior distributions for the parameter a,
and the second column shows the
distributions of the functions, which go
through the origin and spread out as distance
from the origin increases.
• The third column gives the joint GPR distributions
for the two pairs of Y points, for X=(2,9) and
X=(4,7). The distributions have different shapes
(as well as different centers), in large part
because variance increases with distance from
the origin.
• The final example is linear regression with both a
and b as variables. In this case, one observation is
not very informative, so we add another row to
the figure, showing the posterior distributions
when there is a second observation at (1,1).
• Finish discussion of a,b
• Our example was very simple: The parameterized
functions were linear. One can think of linear
regression as a sum of ‘weights’ (the two
parameters) times two fixed ‘basis’ functions (a
constant ‘1’ and y=x), plus added noise.
• Nothing critical changes if we say Y equals a sum of
N weights times N fully specified but arbitrarily
shaped basis functions, plus noise. In GPR data
space we still have a dual representation fully
defined by a mean and covariance function defined
for pairs of Y points.
• Such a regression approach is highly flexible
and can be used in many applications, and
one might ask what is gained by developing a
dual representation in data space.
• There are a number of good answers, some of
which will become clear later, but one is the
fact that GPR does not require specifying
particular functions but rather the relation of
points on the desired function to each other.
• Instead of specifying parameterized functions, let us start in data
space directly, and assume some mean and covariance structure.
• In the usual situation we have N observed data points presumed to
lie on some unknown function. Let us assume we do not have prior
knowledge that this function is in some simply stated
parameterized class.
• It is possible we might have prior knowledge which would allow us
to make an educated guess about the mean of the GPR Gaussian(s).
E.g. if we know that Y should be monotonically increasing we could
incorporate that knowledge in the prior for the means.
• However, in many situations we do not have such prior knowledge,
or that knowledge is so weak we are willing to allow the data points
to ‘speak for themselves’. In such cases we assume an
uninformative prior for the means; for any pair of Y points the prior
mean for that Gaussian (i.e. ‘ellipse’) is usually set to zero.
• To use GPR effectively it is however necessary
to make assumptions about the covariance or
correlation structure. For example, we might
want to impose a constraint that when X
values are near each other, then the
corresponding Y values are near each other.
This tends to produce smoothly varying
functions. We might go further and assume
that this ‘nearness’ constraint falls off
exponentially with distance of the X values
from each other.
• These constraints can be specified by setting
the covariance of points Y(X1) and Y(X2) to be:
• This is called a ‘kernel’. Note that this kernel
function is parameterized. For example the
parameter – controls how fast the Y values can
vary as distance increases. When the observed
Y values for nearby points vary considerably
then this parameter will --------
• In almost any application we are dealing with
a finite and limited set of data points, and
wish to predict only a finite set of points
(perhaps one at a time), making the
specification of the GPR model quite feasible.
For N data points we would have N(N-1)
covariance terms to represent the GPR
Gaussian. Typically these would have some
parameterized structure, reducing the number
of to be estimated quantities.
• Before going further, it would be wise to
illustrate the use of GPR with a simple
example. ???
• The example we used to explain the connection
of Bayesian regression to GPR was as simple as
one could make it. Conceptual understanding is
increased by considering a few generalizations.
• ??Consider first non-linear functions with one
parameter. The logic is the same as before: For a
given value of the parameter, knowing one point
would tell you the values at all others (except for
the added noise). However, the center and shape
of the column 3 distribution would be different,
and would be different for different pairs of
points. ???
• Next consider functions defined by more than
one parameter. One simple example would be
linear regression with both intercept and
slope as parameters.
• Note that GPR is used (usually) when we do
not want to commit to a particular
parameterized class of functions. Instead we
deal with the joint distribution of some finite
set of data points whose form is Gaussian, but
whose mean and covariance comes not from a
parameterized set of assumed functions, but
instead from an assumed mean and
covariance.
• Given we assume no particular functional form
we generally do not know how the mean changes
from one x value to the next, so only specify the
covariance. For N points this could be given as an
N by N matrix of pairwise covariances, but that
would be silly. Instead some simple
parameterized covariance structure is assumed.
• An example is an exponential distance covariance
by which points covary highly when close and less
so as a function of distance along the X axis.
There might be two or three parameters of this
assumed covariance structure and such
parameters are usually termed hyperparameters.
• The assumed covariance structure is termed
the ‘kernel’.
• Note that we do not know the implicitly
defined class of functions that could give rise
to the class of data Gaussians defined by the
(parameterized) kernel. In fact there will
generally be an infinite set of such functions
consistent with some kernel. This is not a
problem, but rather seen as a merit of the
GPR approach: We use GPR when we do not
know how to characterize the functions
explicitly.
• Once we define a covariance kernel, we can
apply Bayes Theorem directly in the data
space, producing both a posterior joint
distribution for the data points, and joint
distribution for the hyperparameters. Doing
this requires a bit of matrix algebra.
• In our simple regression example, we can write the
covariance and mean as:
• In this case we of course see b in the expression since b
is the only unknown parameter in the regression.
• Suppose we did not start with regression, but simply
operated in data space and for some (unknown) reason
assumed both the prior for b and this mean and
covariance as a function of b. If so, then nothing would
change, and the hyperparameter (as we now call it) has
the same prior and posterior as we had in the simple
regression.
• This shows the main point: We can operate
solely in data space and assume a covariance
(and mean) kernel. We can also assume a
prior for either the data point Gaussian or for
the hyperparameters (these are duals of each
other; given the kernel one implies the other).
Given some data we then apply Bayesian
analysis as usual and produce either a
posterior data Gaussian or a posterior for the
hyperparameters (these are duals of each
other, again).
• Suppose we operate solely in data space, the
usual GPR approach. We now see that the
analysis depends completely on the assumed
kernel. So far we have mentioned just two,
one corresponding to our simple linear
regression and one based on exponential
distance (we will say more about this in a
moment).
• Many kernels have been identified (with
different properties that are useful for
different settings), such as:
• A real world example: A few years ago, Dan
Little and I decided to investigate how people
(scientists) formed guesses about the form of
causal relation that underlies 2D noisy data.
Experiment 1:
Participants responded at
various locations on the graph
Responses
• We assumed that the guesses about
functional form could be of two types. One
was parameterized functions. We tried low
order polynomials (e.g.
Y = a + bX + cX2 + dX3 + eX4 + z
• But often people would choose to simply draw
smooth curves through the data as in the
example below:
Similarity-based Functions
• To model both these cases we used GPR. The
next slides show:
• the exponential distance kernel we used for
modeling smooth tracking, and….
• the polynomial kernel we used to model
parameterized functions.
Function hypotheses
Polynomial Kernel (degrees 1 to 4):


 
k x, x*    x   p x*   n2

where   x   x 0 , x1 ,..., x d
Similarity Kernel:



k x, x   exp c x  x
*
2
f

* 2

2
n
Function hypotheses
Polynomial Kernel (degrees 1 to 4):


 
k x, x*    x   p x*   n2

where   x   x 0 , x1 ,..., x d
Similarity Kernel:



k x, x   exp c x  x
*
2
f

* 2

2
n
• We applied this ‘mixture’ model both to the
displayed noisy data pattern on a trial
(showing the Bayesian posterior implied by
the model for that data), and to the data
produced by the observer (showing the
posterior for the observer’s guessed function).
• There were many interesting findings. For
example, the analysis showed that some
people seemed to prefer low order
polynomials (linear mostly) and others
smooth tracking, but all used some of each,
depending on the particular set of data.
[There was much more to the study and
analysis, but this talk is about understanding
GPR].