The issue is known, and it is still not resolved.
This happens because if you have the same points, your matrix is not reversible (singular). (this means that you cannot calculate A ^ -1 - this is part of the solution for the GP).
To solve this problem, simply add some small Gaussian noise to your examples or use another GP library.
You can always try to implement it, in fact it is not so difficult. The most important thing in GP is your kernel function, for example a Gaussian kernel:
exponential_kernel = lambda x, y, params: params[0] * \ np.exp( -0.5 * params[1] * np.sum((x - y)**2) )
Now we need to build a covariance matrix, for example:
covariance = lambda kernel, x, y, params: \ np.array([[kernel(xi, yi, params) for xi in x] for yi in y])
So, if you want to predict a new point x calculate its covariance:
sigma1 = covariance(exponential_kernel, x, x, theta)
and apply the following:
def predict(x, data, kernel, params, sigma, t): k = [kernel(x, y, params) for y in data] Sinv = np.linalg.inv(sigma) y_pred = np.dot(k, Sinv).dot(t) sigma_new = kernel(x, x, params) - np.dot(k, Sinv).dot(k) return y_pred, sigma_new
This is a very naive implementation, and for data with large sizes, the execution time will be high. The hardest part to calculate here is Sinv = np.linalg.inv(sigma) , which takes O(N^3) .