From 5c7caef88564a480d880d815d017147243a7b175 Mon Sep 17 00:00:00 2001 From: James Hensman Date: Mon, 10 Mar 2014 08:34:03 +0000 Subject: [PATCH] tidying in likelihoods --- GPy/examples/non_gaussian.py | 7 ++- .../latent_function_inference/laplace.py | 12 ++--- GPy/likelihoods/likelihood.py | 50 ++++--------------- GPy/likelihoods/poisson.py | 7 ++- 4 files changed, 23 insertions(+), 53 deletions(-) diff --git a/GPy/examples/non_gaussian.py b/GPy/examples/non_gaussian.py index 1c3cab76..7d3e1171 100644 --- a/GPy/examples/non_gaussian.py +++ b/GPy/examples/non_gaussian.py @@ -60,10 +60,9 @@ def student_t_approx(optimize=True, plot=True): m3['.*white'].constrain_fixed(1e-5) m3.randomize() debug = True - print m3 - if debug: - m3.optimize(messages=1) - return m3 + + #TODO: remove + return m3 #Student t GP model on corrupt data t_distribution = GPy.likelihoods.StudentT(deg_free=deg_free, sigma2=edited_real_sd) diff --git a/GPy/inference/latent_function_inference/laplace.py b/GPy/inference/latent_function_inference/laplace.py index 96a47512..94603240 100644 --- a/GPy/inference/latent_function_inference/laplace.py +++ b/GPy/inference/latent_function_inference/laplace.py @@ -51,16 +51,16 @@ class Laplace(object): Ki_f_init = self._previous_Ki_fhat f_hat, Ki_fhat = self.rasm_mode(K, Y, likelihood, Ki_f_init, Y_metadata=Y_metadata) - self.f_hat = f_hat + #Compute hessian and other variables at mode - log_marginal, woodbury_vector, woodbury_inv, dL_dK, dL_dthetaL = self.mode_computations(f_hat, Ki_fhat, K, Y, likelihood, kern, Y_metadata) + log_marginal, woodbury_inv, dL_dK, dL_dthetaL = self.mode_computations(f_hat, Ki_fhat, K, Y, likelihood, kern, Y_metadata) kern.update_gradients_full(dL_dK, X) likelihood.update_gradients(dL_dthetaL) self._previous_Ki_fhat = Ki_fhat.copy() - return Posterior(woodbury_vector=woodbury_vector, woodbury_inv=woodbury_inv, K=K), log_marginal, {'dL_dK':dL_dK} + return Posterior(woodbury_vector=Ki_fhat, woodbury_inv=woodbury_inv, K=K), log_marginal, {'dL_dK':dL_dK} def rasm_mode(self, K, Y, likelihood, Ki_f_init, Y_metadata=None): """ @@ -137,7 +137,6 @@ class Laplace(object): At the mode, compute the hessian and effective covariance matrix. returns: logZ : approximation to the marginal likelihood - woodbury_vector : variable required for calculating the approximation to the covariance matrix woodbury_inv : variable required for calculating the approximation to the covariance matrix dL_dthetaL : array of derivatives (1 x num_kernel_params) dL_dthetaL : array of derivatives (1 x num_likelihood_params) @@ -156,7 +155,6 @@ class Laplace(object): #Compute vival matrices for derivatives dW_df = -likelihood.d3logpdf_df3(f_hat, Y, Y_metadata=Y_metadata) # -d3lik_d3fhat - woodbury_vector = likelihood.dlogpdf_df(f_hat, Y, Y_metadata=Y_metadata) dL_dfhat = -0.5*(np.diag(Ki_W_i)[:, None]*dW_df) #why isn't this -0.5? s2 in R&W p126 line 9. #BiK, _ = dpotrs(L, K, lower=1) #dL_dfhat = 0.5*np.diag(BiK)[:, None]*dW_df @@ -170,7 +168,7 @@ class Laplace(object): explicit_part = 0.5*(np.dot(Ki_f, Ki_f.T) - K_Wi_i) #Implicit - implicit_part = np.dot(woodbury_vector, dL_dfhat.T).dot(I_KW_i) + implicit_part = np.dot(Ki_f, dL_dfhat.T).dot(I_KW_i) dL_dK = explicit_part + implicit_part else: @@ -201,7 +199,7 @@ class Laplace(object): else: dL_dthetaL = np.zeros(likelihood.size) - return log_marginal, woodbury_vector, K_Wi_i, dL_dK, dL_dthetaL + return log_marginal, K_Wi_i, dL_dK, dL_dthetaL def _compute_B_statistics(self, K, W, log_concave): """ diff --git a/GPy/likelihoods/likelihood.py b/GPy/likelihoods/likelihood.py index 070875f5..35010184 100644 --- a/GPy/likelihoods/likelihood.py +++ b/GPy/likelihoods/likelihood.py @@ -120,7 +120,7 @@ class Likelihood(Parameterized): return z, mean, variance - def _predictive_mean(self,mu,variance): + def _predictive_mean(self, mu, variance): """ Quadrature calculation of the predictive mean: E(Y_star|Y) = E( E(Y_star|f_star, Y) ) @@ -135,6 +135,10 @@ class Likelihood(Parameterized): return mean + def _conditional_mean(self, f): + """Quadrature calculation of the conditional mean: E(Y_star|f)""" + raise NotImplementedError, "implement this function to make predictions" + def _predictive_variance(self,mu,variance,predictive_mean=None): """ Numerical approximation to the predictive variance: V(Y_star) @@ -358,50 +362,16 @@ class Likelihood(Parameterized): return dlogpdf_dtheta, dlogpdf_df_dtheta, d2logpdf_df2_dtheta - def predictive_values(self, mu, var, full_cov=False, sampling=False, num_samples=10000): + def predictive_values(self, mu, var): """ - Compute mean, variance and conficence interval (percentiles 5 and 95) of the prediction. + Compute mean, variance of the predictive distibution. :param mu: mean of the latent variable, f, of posterior :param var: variance of the latent variable, f, of posterior - :param full_cov: whether to use the full covariance or just the diagonal - :type full_cov: Boolean - :param num_samples: number of samples to use in computing quantiles and - possibly mean variance - :type num_samples: integer - :param sampling: Whether to use samples for mean and variances anyway - :type sampling: Boolean - """ - - if sampling: - #Get gp_samples f* using posterior mean and variance - if not full_cov: - gp_samples = np.random.multivariate_normal(mu.flatten(), np.diag(var.flatten()), - size=num_samples).T - else: - gp_samples = np.random.multivariate_normal(mu.flatten(), var, - size=num_samples).T - #Push gp samples (f*) through likelihood to give p(y*|f*) - samples = self.samples(gp_samples) - axis=-1 - - #Calculate mean, variance and precentiles from samples - print "WARNING: Using sampling to calculate mean, variance and predictive quantiles." - pred_mean = np.mean(samples, axis=axis)[:,None] - pred_var = np.var(samples, axis=axis)[:,None] - q1 = np.percentile(samples, 2.5, axis=axis)[:,None] - q3 = np.percentile(samples, 97.5, axis=axis)[:,None] - - else: - - pred_mean = self.predictive_mean(mu, var) - pred_var = self.predictive_variance(mu, var, pred_mean) - print "WARNING: Predictive quantiles are only computed when sampling." - q1 = np.repeat(np.nan,pred_mean.size)[:,None] - q3 = q1.copy() - - return pred_mean, pred_var, q1, q3 + pred_mean = self.predictive_mean(mu, var) + pred_var = self.predictive_variance(mu, var, pred_mean) + return pred_mean, pred_var def samples(self, gp): """ diff --git a/GPy/likelihoods/poisson.py b/GPy/likelihoods/poisson.py index 616447d9..a306cbf1 100644 --- a/GPy/likelihoods/poisson.py +++ b/GPy/likelihoods/poisson.py @@ -25,8 +25,11 @@ class Poisson(Likelihood): super(Poisson, self).__init__(gp_link, name='Poisson') - def _preprocess_values(self,Y): - return Y + def _conditional_mean(self, f): + """ + the expected value of y given a value of f + """ + return self.gp_link.transf(gp) def pdf_link(self, link_f, y, Y_metadata=None): """