From 9d7b670160684d760136737b18237ae5405c5c97 Mon Sep 17 00:00:00 2001 From: Alan Saul Date: Thu, 19 Sep 2013 15:56:18 +0100 Subject: [PATCH] Tests setup but not fitting properly yet --- GPy/examples/laplace_approximations.py | 87 +++++++++++++++++++------- 1 file changed, 65 insertions(+), 22 deletions(-) diff --git a/GPy/examples/laplace_approximations.py b/GPy/examples/laplace_approximations.py index 3e24c89f..1ad4eb38 100644 --- a/GPy/examples/laplace_approximations.py +++ b/GPy/examples/laplace_approximations.py @@ -659,9 +659,10 @@ def boston_example(): Y = data['Y'].copy() Y = Y-Y.mean() Y = Y/Y.std() - num_folds = 2 + import ipdb; ipdb.set_trace() # XXX BREAKPOINT + num_folds = 10 kf = KFold(len(Y), n_folds=num_folds, indices=True) - score_folds = np.zeros((3, num_folds)) + score_folds = np.zeros((4, num_folds)) def rmse(Y, Ystar): return np.sqrt(np.mean((Y-Ystar)**2)) #for train, test in kf: @@ -673,56 +674,98 @@ def boston_example(): #Gaussian GP print "Gauss GP" - kernelgp = GPy.kern.rbf(X.shape[1]) #+ GPy.kern.white(X.shape[1]) + kernelgp = GPy.kern.rbf(X.shape[1]) + GPy.kern.white(X.shape[1], variance=0.01) mgp = GPy.models.GPRegression(X_train.copy(), Y_train.copy(), kernel=kernelgp) mgp.ensure_default_constraints() mgp['noise'] = noise + mgp.constrain_fixed('white', 0.01) + print mgp mgp.optimize(messages=1) Y_test_pred = mgp.predict(X_test) score_folds[0, n] = rmse(Y_test, Y_test_pred[0]) - plt.figure() - plt.scatter(X_test[:, 0], Y_test_pred[0]) - plt.scatter(X_test[:, 0], Y_test, c='r', marker='x') + print mgp print score_folds - plt.title('GP gauss') + #plt.figure() + #plt.scatter(X_test[:, 0], Y_test_pred[0]) + #plt.scatter(X_test[:, 0], Y_test, c='r', marker='x') + #plt.title('GP gauss') print "Gaussian Laplace GP" sigma2_start = 1 - kernelstu = GPy.kern.rbf(X.shape[1]) #+ GPy.kern.white(X.shape[1], variance=0.01) + kernelstu = GPy.kern.rbf(X.shape[1]) + GPy.kern.white(X.shape[1], variance=0.1) N, D = Y_train.shape g_distribution = GPy.likelihoods.functions.Gaussian(variance=noise, N=N, D=D) g_likelihood = GPy.likelihoods.Laplace(Y_train.copy(), g_distribution, opt='rasm') mg = GPy.models.GPRegression(X_train.copy(), Y_train.copy(), kernel=kernelstu, likelihood=g_likelihood) mg.ensure_default_constraints() mg.constrain_positive('noise_variance') - mg.optimize(messages=1) + mg.constrain_fixed('white', 0.01) + mg['noise'] = noise + print mg + try: + mg.optimize(messages=1) + except Exception: + print "Blew up" Y_test_pred = mg.predict(X_test) score_folds[1, n] = rmse(Y_test, Y_test_pred[0]) print score_folds - plt.figure() - plt.scatter(X_test[:, 0], Y_test_pred[0]) - plt.scatter(X_test[:, 0], Y_test, c='r', marker='x') - plt.title('Lap gauss') + print mg + #plt.figure() + #plt.scatter(X_test[:, 0], Y_test_pred[0]) + #plt.scatter(X_test[:, 0], Y_test, c='r', marker='x') + #plt.title('Lap gauss') #Student t likelihood - print "Student-T GP" deg_free = 5 - kernelstu = GPy.kern.rbf(X.shape[1]) #+ GPy.kern.white(X.shape[1], variance=0.01) + print "Student-T GP {}df".format(deg_free) + kernelstu = GPy.kern.rbf(X.shape[1]) + GPy.kern.white(X.shape[1], variance=0.1) t_distribution = GPy.likelihoods.functions.StudentT(deg_free, sigma2=noise) stu_t_likelihood = GPy.likelihoods.Laplace(Y_train.copy(), t_distribution, opt='rasm') mstu_t = GPy.models.GPRegression(X_train.copy(), Y_train.copy(), kernel=kernelstu, likelihood=stu_t_likelihood) mstu_t.ensure_default_constraints() + mstu_t.constrain_fixed('white', 0.01) #mstu_t.constrain_positive('t_noise') - mstu_t.constrain_bounded('t_noise', 0.01, 1000) - mstu_t.optimize(messages=1) + mstu_t.constrain_bounded('t_noise', 0.001, 1000) + mstu_t['t_noise'] = noise + print mstu_t + try: + mstu_t.optimize(messages=1) + except Exception: + print "Blew up" Y_test_pred = mstu_t.predict(X_test) score_folds[2, n] = rmse(Y_test, Y_test_pred[0]) print score_folds - plt.figure() - plt.scatter(X_test[:, 0], Y_test_pred[0]) - plt.scatter(X_test[:, 0], Y_test, c='r', marker='x') - plt.title('Stu t') - import ipdb; ipdb.set_trace() # XXX BREAKPOINT + print mstu_t + #plt.figure() + #plt.scatter(X_test[:, 0], Y_test_pred[0]) + #plt.scatter(X_test[:, 0], Y_test, c='r', marker='x') + #plt.title('Stu t {}df'.format(deg_free)) + + deg_free = 3 + print "Student-T GP {}df".format(deg_free) + kernelstu = GPy.kern.rbf(X.shape[1]) + GPy.kern.white(X.shape[1], variance=0.1) + t_distribution = GPy.likelihoods.functions.StudentT(deg_free, sigma2=noise) + stu_t_likelihood = GPy.likelihoods.Laplace(Y_train.copy(), t_distribution, opt='rasm') + mstu_t = GPy.models.GPRegression(X_train.copy(), Y_train.copy(), kernel=kernelstu, likelihood=stu_t_likelihood) + mstu_t.ensure_default_constraints() + mstu_t.constrain_fixed('white', 0.01) + #mstu_t.constrain_positive('t_noise') + mstu_t.constrain_bounded('t_noise', 0.001, 1000) + mstu_t['t_noise'] = noise + print mstu_t + try: + mstu_t.optimize(messages=1) + except Exception: + print "Blew up" + mstu_t.optimize(messages=1) + Y_test_pred = mstu_t.predict(X_test) + score_folds[3, n] = rmse(Y_test, Y_test_pred[0]) + print score_folds + print mstu_t + #plt.figure() + #plt.scatter(X_test[:, 0], Y_test_pred[0]) + #plt.scatter(X_test[:, 0], Y_test, c='r', marker='x') + #plt.title('Stu t {}df'.format(deg_free)) def plot_f_approx(model):