Gaussian Process Prediction Uncertainty

Question

I am using Gaussian Process Regression to interpolate my input points. I would like to measure the total uncertainty of my prediction thus I sum up the GPR prediction variances at all the testing points $x_i, i = \overline{1, N}$. $$ V = \sum_{i = 0}^N \sigma^2(x_i) $$

I was expecting that this total variance $V$ would reduce when I add more samples, and this is the case when I don't adapt the hyperparameters. However, when I adapt the hyperparameters using Maximum Likelihood, $V$ can increase when new points are added.

I found this quite counter-intuitive. Isn't hyperparameter adapting supposed to reduce the variance? Could some one please help me clarify this point. Thanks a lot!

I attach a simple Matlab script, based on GPML for your testing (http://www.gaussianprocess.org/gpml/code/matlab/doc/index.html)

covfunc = @covSEiso; hyp.cov = [log(0.1); log(1.0)]; 
likfunc = @likGauss; hyp.lik = log(0.01);

z = linspace(-1.0, 1.0, 101); z = z';

figure; grid on;
axis([-1.1, 1.1, -2, 2]);
x = [];
y = [];

ADAPT_PARAM = false;

while 1
    [xx, yy, button] = ginput(1);

    if button == 3
        break;
    end

    x = [x; xx];
    y = [y; yy];

    if ADAPT_PARAM
    % Adapt the hyperparams
        hyp.cov = [log(0.1); log(1.0)]; 
        hyp = minimize(hyp, @gp, -100, @infExact, [], covfunc, likfunc, x, y);
    end

    [m, s2] = gp(hyp, @infExact, [], covfunc, likfunc, x, y, z);
    f = [m+2*sqrt(s2); flipdim(m-2*sqrt(s2),1)];

    fill([z; flipdim(z,1)], f, [7 7 7]/8); 
    axis([-1.1, 1.1, -2, 2]); grid on; hold on;
    plot(z, m, 'LineWidth', 2);
    plot(x, y, 'r+', 'MarkerSize', 12)
    hold off;

    disp(['Total variance ', num2str(sum(s2))]);
end

Answer 1

Isn't hyperparameter adapting supposed to reduce the variance?

Nope, not necessarily! The most common way to tune parameters is to maximise the log marginal likelihood

$$\log P(y|x,\theta) = -0.5y^TK^{-1}y - 0.5\log|K| - c$$

Notice that there are two terms at play here and a constant which I'll ignore. The first is the data fit term which is maximized when the data fits the model very well (lower $V$). The second is a penalty on the complexity of the model, ie smoother the better (higher $V$). When you optimise you try and find the balance between the two and this changes with the data you observe.

Check out Chapter 5 of GPML.