I tried to implement a Gaussian process in octave.
As a starting point I used the algorithm described on page 19 of Rasmussens GP book (http://www.gaussianprocess.org/gpml/).
As a covariance matrix I used the squared exponential function (as it is used in the book as well):
function re = k(x1, x2) re = exp(-(1/2.0) * abs(x1.-x2).^2); endfunction
And calculate the covariance matrix (of the training inputs with):
# Calculate covariance matrix s = size(X); K = ; for i = 1:s for j = 1:i re = k(X(i), X(j)); K(i,j) = re; K(j,i) = re; endfor endfor
But for some reason the resulting covariance matrix K sometimes is not positive definite (depending on inputs X).
So can anyone tell me what I'm doing wrong here, please? And is there a way to test whether a covariance function results in a positive definite covariance matrix? Since the squared exponential function seems to be a covariance function, I assumed it should create a positive definite matrix.