Merge branch 'params' of github.com:SheffieldML/GPy into params

Conflicts:
	GPy/likelihoods/mixed_noise.py
This commit is contained in:
Ricardo 2014-03-19 09:39:35 +00:00
commit fcb6742b60
21 changed files with 197 additions and 228 deletions

View file

@ -23,13 +23,10 @@ K = Bias.prod(Coreg,name='X')
#K.coregion.W = 0 #K.coregion.W = 0
#print K.coregion.W #print K.coregion.W
#print Bias.K(_X,_X) #print Bias.K(_X,_X)
#print K.K(X,X) #print K.K(X,X)
#pb.matshow(K.K(X,X)) #pb.matshow(K.K(X,X))
Mlist = [GPy.kern.Matern32(1,lengthscale=20.,name="Mat")] Mlist = [GPy.kern.Matern32(1,lengthscale=20.,name="Mat")]
kern = GPy.util.multioutput.LCM(input_dim=1,num_outputs=2,kernels_list=Mlist,name='H') kern = GPy.util.multioutput.LCM(input_dim=1,num_outputs=2,kernels_list=Mlist,name='H')
kern.B.W = 0 kern.B.W = 0
@ -37,16 +34,22 @@ kern.B.kappa = 1.
#kern.B.W.fix() #kern.B.W.fix()
#kern.B.kappa.fix() #kern.B.kappa.fix()
#m = GPy.models.GPCoregionalizedRegression(X_list=[X1,X2], Y_list=[Y1,Y2], kernel=kern) #m = GPy.models.GPCoregionalizedRegression(X_list=[X1,X2], Y_list=[Y1,Y2], kernel=kern)
m = GPy.models.SparseGPCoregionalizedRegression(X_list=[X1], Y_list=[Y1], kernel=kern)
Z1 = np.array([1.5,2.5])[:,None]
m = GPy.models.SparseGPCoregionalizedRegression(X_list=[X1], Y_list=[Y1], Z_list = [Z1], kernel=kern)
#m.optimize() #m.optimize()
m.checkgrad(verbose=1) m.checkgrad(verbose=1)
"""
fig = pb.figure() fig = pb.figure()
ax0 = fig.add_subplot(211) ax0 = fig.add_subplot(211)
ax1 = fig.add_subplot(212) ax1 = fig.add_subplot(212)
slices = GPy.util.multioutput.get_slices([Y1,Y2]) slices = GPy.util.multioutput.get_slices([Y1,Y2])
m.plot(fixed_inputs=[(1,0)],which_data_rows=slices[0],ax=ax0) m.plot(fixed_inputs=[(1,0)],which_data_rows=slices[0],ax=ax0)
#m.plot(fixed_inputs=[(1,1)],which_data_rows=slices[1],ax=ax1) #m.plot(fixed_inputs=[(1,1)],which_data_rows=slices[1],ax=ax1)
"""

View file

@ -158,7 +158,7 @@ def boston_example(optimize=True, plot=True):
#Gaussian GP #Gaussian GP
print "Gauss GP" print "Gauss GP"
mgp = GPy.models.GPRegression(X_train.copy(), Y_train.copy(), kernel=kernelgp.copy()) mgp = GPy.models.GPRegression(X_train.copy(), Y_train.copy(), kernel=kernelgp.copy())
mgp.constrain_fixed('white', 1e-5) mgp.constrain_fixed('.*white', 1e-5)
mgp['rbf_len'] = rbf_len mgp['rbf_len'] = rbf_len
mgp['noise'] = noise mgp['noise'] = noise
print mgp print mgp
@ -176,7 +176,7 @@ def boston_example(optimize=True, plot=True):
g_likelihood = GPy.likelihoods.Laplace(Y_train.copy(), g_distribution) g_likelihood = GPy.likelihoods.Laplace(Y_train.copy(), g_distribution)
mg = GPy.models.GPRegression(X_train.copy(), Y_train.copy(), kernel=kernelstu.copy(), likelihood=g_likelihood) mg = GPy.models.GPRegression(X_train.copy(), Y_train.copy(), kernel=kernelstu.copy(), likelihood=g_likelihood)
mg.constrain_positive('noise_variance') mg.constrain_positive('noise_variance')
mg.constrain_fixed('white', 1e-5) mg.constrain_fixed('.*white', 1e-5)
mg['rbf_len'] = rbf_len mg['rbf_len'] = rbf_len
mg['noise'] = noise mg['noise'] = noise
print mg print mg
@ -194,10 +194,10 @@ def boston_example(optimize=True, plot=True):
t_distribution = GPy.likelihoods.noise_model_constructors.student_t(deg_free=df, sigma2=noise) t_distribution = GPy.likelihoods.noise_model_constructors.student_t(deg_free=df, sigma2=noise)
stu_t_likelihood = GPy.likelihoods.Laplace(Y_train.copy(), t_distribution) stu_t_likelihood = GPy.likelihoods.Laplace(Y_train.copy(), t_distribution)
mstu_t = GPy.models.GPRegression(X_train.copy(), Y_train.copy(), kernel=kernelstu.copy(), likelihood=stu_t_likelihood) mstu_t = GPy.models.GPRegression(X_train.copy(), Y_train.copy(), kernel=kernelstu.copy(), likelihood=stu_t_likelihood)
mstu_t.constrain_fixed('white', 1e-5) mstu_t.constrain_fixed('.*white', 1e-5)
mstu_t.constrain_bounded('t_noise', 0.0001, 1000) mstu_t.constrain_bounded('.*t_noise', 0.0001, 1000)
mstu_t['rbf_len'] = rbf_len mstu_t['rbf_len'] = rbf_len
mstu_t['t_noise'] = noise mstu_t['.*t_noise'] = noise
print mstu_t print mstu_t
if optimize: if optimize:
mstu_t.optimize(optimizer=optimizer, messages=messages) mstu_t.optimize(optimizer=optimizer, messages=messages)

View file

@ -43,7 +43,7 @@ class ExactGaussianInference(object):
K = kern.K(X) K = kern.K(X)
Ky = K.copy() Ky = K.copy()
diag.add(Ky, likelihood.gaussian_variance(Y, Y_metadata)) diag.add(Ky, likelihood.gaussian_variance(Y_metadata))
Wi, LW, LWi, W_logdet = pdinv(Ky) Wi, LW, LWi, W_logdet = pdinv(Ky)
alpha, _ = dpotrs(LW, YYT_factor, lower=1) alpha, _ = dpotrs(LW, YYT_factor, lower=1)

View file

@ -51,12 +51,11 @@ class Laplace(object):
Ki_f_init = self._previous_Ki_fhat Ki_f_init = self._previous_Ki_fhat
f_hat, Ki_fhat = self.rasm_mode(K, Y, likelihood, Ki_f_init, Y_metadata=Y_metadata) f_hat, Ki_fhat = self.rasm_mode(K, Y, likelihood, Ki_f_init, Y_metadata=Y_metadata)
self.f_hat = f_hat self.f_hat = f_hat
self.Ki_fhat = Ki_fhat self.Ki_fhat = Ki_fhat
self.K = K.copy() self.K = K.copy()
#Compute hessian and other variables at mode #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)
self._previous_Ki_fhat = Ki_fhat.copy() self._previous_Ki_fhat = Ki_fhat.copy()
return Posterior(woodbury_vector=Ki_fhat, woodbury_inv=woodbury_inv, K=K), log_marginal, {'dL_dK':dL_dK, 'dL_dthetaL':dL_dthetaL} return Posterior(woodbury_vector=Ki_fhat, woodbury_inv=woodbury_inv, K=K), log_marginal, {'dL_dK':dL_dK, 'dL_dthetaL':dL_dthetaL}
@ -86,13 +85,13 @@ class Laplace(object):
#define the objective function (to be maximised) #define the objective function (to be maximised)
def obj(Ki_f, f): def obj(Ki_f, f):
return -0.5*np.dot(Ki_f.flatten(), f.flatten()) + likelihood.logpdf(f, Y, extra_data=Y_metadata) return -0.5*np.dot(Ki_f.flatten(), f.flatten()) + likelihood.logpdf(f, Y, Y_metadata=Y_metadata)
difference = np.inf difference = np.inf
iteration = 0 iteration = 0
while difference > self._mode_finding_tolerance and iteration < self._mode_finding_max_iter: while difference > self._mode_finding_tolerance and iteration < self._mode_finding_max_iter:
W = -likelihood.d2logpdf_df2(f, Y, extra_data=Y_metadata) W = -likelihood.d2logpdf_df2(f, Y, Y_metadata=Y_metadata)
grad = likelihood.dlogpdf_df(f, Y, extra_data=Y_metadata) grad = likelihood.dlogpdf_df(f, Y, Y_metadata=Y_metadata)
W_f = W*f W_f = W*f
@ -136,13 +135,12 @@ class Laplace(object):
At the mode, compute the hessian and effective covariance matrix. At the mode, compute the hessian and effective covariance matrix.
returns: logZ : approximation to the marginal likelihood 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 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_kernel_params)
dL_dthetaL : array of derivatives (1 x num_likelihood_params) dL_dthetaL : array of derivatives (1 x num_likelihood_params)
""" """
#At this point get the hessian matrix (or vector as W is diagonal) #At this point get the hessian matrix (or vector as W is diagonal)
W = -likelihood.d2logpdf_df2(f_hat, Y, extra_data=Y_metadata) W = -likelihood.d2logpdf_df2(f_hat, Y, Y_metadata=Y_metadata)
K_Wi_i, L, LiW12 = self._compute_B_statistics(K, W, likelihood.log_concave) K_Wi_i, L, LiW12 = self._compute_B_statistics(K, W, likelihood.log_concave)
@ -151,11 +149,10 @@ class Laplace(object):
Ki_W_i = K - C.T.dot(C) #Could this be wrong? Ki_W_i = K - C.T.dot(C) #Could this be wrong?
#compute the log marginal #compute the log marginal
log_marginal = -0.5*np.dot(Ki_f.flatten(), f_hat.flatten()) + likelihood.logpdf(f_hat, Y, extra_data=Y_metadata) - np.sum(np.log(np.diag(L))) log_marginal = -0.5*np.dot(Ki_f.flatten(), f_hat.flatten()) + likelihood.logpdf(f_hat, Y, Y_metadata=Y_metadata) - np.sum(np.log(np.diag(L)))
#Compute vival matrices for derivatives #Compute vival matrices for derivatives
dW_df = -likelihood.d3logpdf_df3(f_hat, Y, extra_data=Y_metadata) # -d3lik_d3fhat dW_df = -likelihood.d3logpdf_df3(f_hat, Y, Y_metadata=Y_metadata) # -d3lik_d3fhat
woodbury_vector = likelihood.dlogpdf_df(f_hat, Y, extra_data=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. 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) #BiK, _ = dpotrs(L, K, lower=1)
#dL_dfhat = 0.5*np.diag(BiK)[:, None]*dW_df #dL_dfhat = 0.5*np.diag(BiK)[:, None]*dW_df
@ -169,7 +166,7 @@ class Laplace(object):
explicit_part = 0.5*(np.dot(Ki_f, Ki_f.T) - K_Wi_i) explicit_part = 0.5*(np.dot(Ki_f, Ki_f.T) - K_Wi_i)
#Implicit #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 dL_dK = explicit_part + implicit_part
else: else:
@ -179,7 +176,7 @@ class Laplace(object):
#compute dL_dthetaL# #compute dL_dthetaL#
#################### ####################
if likelihood.size > 0 and not likelihood.is_fixed: if likelihood.size > 0 and not likelihood.is_fixed:
dlik_dthetaL, dlik_grad_dthetaL, dlik_hess_dthetaL = likelihood._laplace_gradients(f_hat, Y, extra_data=Y_metadata) dlik_dthetaL, dlik_grad_dthetaL, dlik_hess_dthetaL = likelihood._laplace_gradients(f_hat, Y, Y_metadata=Y_metadata)
num_params = likelihood.size num_params = likelihood.size
# make space for one derivative for each likelihood parameter # make space for one derivative for each likelihood parameter
@ -200,7 +197,7 @@ class Laplace(object):
else: else:
dL_dthetaL = np.zeros(likelihood.size) 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): def _compute_B_statistics(self, K, W, log_concave):
""" """

View file

@ -65,7 +65,7 @@ class VarDTC(object):
_, output_dim = Y.shape _, output_dim = Y.shape
#see whether we've got a different noise variance for each datum #see whether we've got a different noise variance for each datum
beta = 1./np.fmax(likelihood.gaussian_variance(Y, Y_metadata), 1e-6) beta = 1./np.fmax(likelihood.gaussian_variance(Y_metadata), 1e-6)
# VVT_factor is a matrix such that tdot(VVT_factor) = VVT...this is for efficiency! # VVT_factor is a matrix such that tdot(VVT_factor) = VVT...this is for efficiency!
#self.YYTfactor = self.get_YYTfactor(Y) #self.YYTfactor = self.get_YYTfactor(Y)
#VVT_factor = self.get_VVTfactor(self.YYTfactor, beta) #VVT_factor = self.get_VVTfactor(self.YYTfactor, beta)
@ -176,7 +176,6 @@ class VarDTC(object):
#construct a posterior object #construct a posterior object
post = Posterior(woodbury_inv=woodbury_inv, woodbury_vector=woodbury_vector, K=Kmm, mean=None, cov=None, K_chol=Lm) post = Posterior(woodbury_inv=woodbury_inv, woodbury_vector=woodbury_vector, K=Kmm, mean=None, cov=None, K_chol=Lm)
return post, log_marginal, grad_dict return post, log_marginal, grad_dict
class VarDTCMissingData(object): class VarDTCMissingData(object):
@ -365,7 +364,7 @@ class VarDTCMissingData(object):
return post, log_marginal, grad_dict return post, log_marginal, grad_dict
def _compute_dL_dpsi(num_inducing, num_data, output_dim, beta, Lm, VVT_factor, Cpsi1Vf, DBi_plus_BiPBi, psi1, het_noise, uncertain_inputs): def _compute_dL_dpsi(num_inducing, num_data, output_dim, beta, Lm, VVT_factor, Cpsi1Vf, DBi_plus_BiPBi, psi1, het_noise, uncertain_inputs):
dL_dpsi0 = -0.5 * output_dim * (beta * np.ones([num_data, 1])).flatten() dL_dpsi0 = -0.5 * output_dim * (beta[:,None] * np.ones([num_data, 1])).flatten()
dL_dpsi1 = np.dot(VVT_factor, Cpsi1Vf.T) dL_dpsi1 = np.dot(VVT_factor, Cpsi1Vf.T)
dL_dpsi2_beta = 0.5 * backsub_both_sides(Lm, output_dim * np.eye(num_inducing) - DBi_plus_BiPBi) dL_dpsi2_beta = 0.5 * backsub_both_sides(Lm, output_dim * np.eye(num_inducing) - DBi_plus_BiPBi)
if het_noise: if het_noise:

View file

@ -95,7 +95,7 @@ class Bernoulli(Likelihood):
else: else:
return np.nan return np.nan
def pdf_link(self, link_f, y, extra_data=None): def pdf_link(self, link_f, y, Y_metadata=None):
""" """
Likelihood function given link(f) Likelihood function given link(f)
@ -106,7 +106,7 @@ class Bernoulli(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in bernoulli :param Y_metadata: Y_metadata not used in bernoulli
:returns: likelihood evaluated for this point :returns: likelihood evaluated for this point
:rtype: float :rtype: float
@ -118,7 +118,7 @@ class Bernoulli(Likelihood):
objective = np.where(y, link_f, 1.-link_f) objective = np.where(y, link_f, 1.-link_f)
return np.exp(np.sum(np.log(objective))) return np.exp(np.sum(np.log(objective)))
def logpdf_link(self, link_f, y, extra_data=None): def logpdf_link(self, link_f, y, Y_metadata=None):
""" """
Log Likelihood function given link(f) Log Likelihood function given link(f)
@ -129,7 +129,7 @@ class Bernoulli(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in bernoulli :param Y_metadata: Y_metadata not used in bernoulli
:returns: log likelihood evaluated at points link(f) :returns: log likelihood evaluated at points link(f)
:rtype: float :rtype: float
""" """
@ -140,7 +140,7 @@ class Bernoulli(Likelihood):
np.seterr(**state) np.seterr(**state)
return np.sum(objective) return np.sum(objective)
def dlogpdf_dlink(self, link_f, y, extra_data=None): def dlogpdf_dlink(self, link_f, y, Y_metadata=None):
""" """
Gradient of the pdf at y, given link(f) w.r.t link(f) Gradient of the pdf at y, given link(f) w.r.t link(f)
@ -151,7 +151,7 @@ class Bernoulli(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in bernoulli :param Y_metadata: Y_metadata not used in bernoulli
:returns: gradient of log likelihood evaluated at points link(f) :returns: gradient of log likelihood evaluated at points link(f)
:rtype: Nx1 array :rtype: Nx1 array
""" """
@ -162,7 +162,7 @@ class Bernoulli(Likelihood):
np.seterr(**state) np.seterr(**state)
return grad return grad
def d2logpdf_dlink2(self, link_f, y, extra_data=None): def d2logpdf_dlink2(self, link_f, y, Y_metadata=None):
""" """
Hessian at y, given link_f, w.r.t link_f the hessian will be 0 unless i == j Hessian at y, given link_f, w.r.t link_f the hessian will be 0 unless i == j
i.e. second derivative logpdf at y given link(f_i) link(f_j) w.r.t link(f_i) and link(f_j) i.e. second derivative logpdf at y given link(f_i) link(f_j) w.r.t link(f_i) and link(f_j)
@ -175,7 +175,7 @@ class Bernoulli(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in bernoulli :param Y_metadata: Y_metadata not used in bernoulli
:returns: Diagonal of log hessian matrix (second derivative of log likelihood evaluated at points link(f)) :returns: Diagonal of log hessian matrix (second derivative of log likelihood evaluated at points link(f))
:rtype: Nx1 array :rtype: Nx1 array
@ -190,7 +190,7 @@ class Bernoulli(Likelihood):
np.seterr(**state) np.seterr(**state)
return d2logpdf_dlink2 return d2logpdf_dlink2
def d3logpdf_dlink3(self, link_f, y, extra_data=None): def d3logpdf_dlink3(self, link_f, y, Y_metadata=None):
""" """
Third order derivative log-likelihood function at y given link(f) w.r.t link(f) Third order derivative log-likelihood function at y given link(f) w.r.t link(f)
@ -201,7 +201,7 @@ class Bernoulli(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in bernoulli :param Y_metadata: Y_metadata not used in bernoulli
:returns: third derivative of log likelihood evaluated at points link(f) :returns: third derivative of log likelihood evaluated at points link(f)
:rtype: Nx1 array :rtype: Nx1 array
""" """

View file

@ -18,13 +18,12 @@ class Exponential(Likelihood):
L(x) = \exp(\lambda) * \lambda**Y_i / Y_i! L(x) = \exp(\lambda) * \lambda**Y_i / Y_i!
$$ $$
""" """
def __init__(self,gp_link=None,analytical_mean=False,analytical_variance=False): def __init__(self,gp_link=None):
super(Exponential, self).__init__(gp_link,analytical_mean,analytical_variance) if gp_link is None:
gp_link = link_functions.Log()
super(Exponential, self).__init__(gp_link, 'ExpLikelihood')
def _preprocess_values(self,Y): def pdf_link(self, link_f, y, Y_metadata=None):
return Y
def pdf_link(self, link_f, y, extra_data=None):
""" """
Likelihood function given link(f) Likelihood function given link(f)
@ -35,16 +34,15 @@ class Exponential(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in exponential distribution :param Y_metadata: Y_metadata which is not used in exponential distribution
:returns: likelihood evaluated for this point :returns: likelihood evaluated for this point
:rtype: float :rtype: float
""" """
assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape
log_objective = link_f*np.exp(-y*link_f) log_objective = link_f*np.exp(-y*link_f)
return np.exp(np.sum(np.log(log_objective))) return np.exp(np.sum(np.log(log_objective)))
#return np.exp(np.sum(-y/link_f - np.log(link_f) ))
def logpdf_link(self, link_f, y, extra_data=None): def logpdf_link(self, link_f, y, Y_metadata=None):
""" """
Log Likelihood Function given link(f) Log Likelihood Function given link(f)
@ -55,17 +53,16 @@ class Exponential(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in exponential distribution :param Y_metadata: Y_metadata which is not used in exponential distribution
:returns: likelihood evaluated for this point :returns: likelihood evaluated for this point
:rtype: float :rtype: float
""" """
assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape
log_objective = np.log(link_f) - y*link_f log_objective = np.log(link_f) - y*link_f
#logpdf_link = np.sum(-np.log(link_f) - y/link_f)
return np.sum(log_objective) return np.sum(log_objective)
def dlogpdf_dlink(self, link_f, y, extra_data=None): def dlogpdf_dlink(self, link_f, y, Y_metadata=None):
""" """
Gradient of the log likelihood function at y, given link(f) w.r.t link(f) Gradient of the log likelihood function at y, given link(f) w.r.t link(f)
@ -76,7 +73,7 @@ class Exponential(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in exponential distribution :param Y_metadata: Y_metadata which is not used in exponential distribution
:returns: gradient of likelihood evaluated at points :returns: gradient of likelihood evaluated at points
:rtype: Nx1 array :rtype: Nx1 array
@ -86,7 +83,7 @@ class Exponential(Likelihood):
#grad = y/(link_f**2) - 1./link_f #grad = y/(link_f**2) - 1./link_f
return grad return grad
def d2logpdf_dlink2(self, link_f, y, extra_data=None): def d2logpdf_dlink2(self, link_f, y, Y_metadata=None):
""" """
Hessian at y, given link(f), w.r.t link(f) Hessian at y, given link(f), w.r.t link(f)
i.e. second derivative logpdf at y given link(f_i) and link(f_j) w.r.t link(f_i) and link(f_j) i.e. second derivative logpdf at y given link(f_i) and link(f_j) w.r.t link(f_i) and link(f_j)
@ -99,7 +96,7 @@ class Exponential(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in exponential distribution :param Y_metadata: Y_metadata which is not used in exponential distribution
:returns: Diagonal of hessian matrix (second derivative of likelihood evaluated at points f) :returns: Diagonal of hessian matrix (second derivative of likelihood evaluated at points f)
:rtype: Nx1 array :rtype: Nx1 array
@ -112,7 +109,7 @@ class Exponential(Likelihood):
#hess = -2*y/(link_f**3) + 1/(link_f**2) #hess = -2*y/(link_f**3) + 1/(link_f**2)
return hess return hess
def d3logpdf_dlink3(self, link_f, y, extra_data=None): def d3logpdf_dlink3(self, link_f, y, Y_metadata=None):
""" """
Third order derivative log-likelihood function at y given link(f) w.r.t link(f) Third order derivative log-likelihood function at y given link(f) w.r.t link(f)
@ -123,7 +120,7 @@ class Exponential(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in exponential distribution :param Y_metadata: Y_metadata which is not used in exponential distribution
:returns: third derivative of likelihood evaluated at points f :returns: third derivative of likelihood evaluated at points f
:rtype: Nx1 array :rtype: Nx1 array
""" """
@ -132,18 +129,6 @@ class Exponential(Likelihood):
#d3lik_dlink3 = 6*y/(link_f**4) - 2./(link_f**3) #d3lik_dlink3 = 6*y/(link_f**4) - 2./(link_f**3)
return d3lik_dlink3 return d3lik_dlink3
def _mean(self,gp):
"""
Mass (or density) function
"""
return self.gp_link.transf(gp)
def _variance(self,gp):
"""
Mass (or density) function
"""
return self.gp_link.transf(gp)**2
def samples(self, gp): def samples(self, gp):
""" """
Returns a set of samples of observations based on a given value of the latent variable. Returns a set of samples of observations based on a given value of the latent variable.

View file

@ -1,11 +1,12 @@
# Copyright (c) 2012, 2013 Ricardo Andrade # Copyright (c) 2012 - 2014, GPy authors (see AUTHORS.txt).
# Licensed under the BSD 3-clause license (see LICENSE.txt) # Licensed under the BSD 3-clause license (see LICENSE.txt)
import numpy as np import numpy as np
from scipy import stats,special from scipy import stats,special
import scipy as sp import scipy as sp
from GPy.util.univariate_Gaussian import std_norm_pdf,std_norm_cdf from ..util.univariate_Gaussian import std_norm_pdf,std_norm_cdf
from ..core.parameterization import Param
import link_functions import link_functions
from likelihood import Likelihood from likelihood import Likelihood
@ -18,14 +19,16 @@ class Gamma(Likelihood):
\\alpha_{i} = \\beta y_{i} \\alpha_{i} = \\beta y_{i}
""" """
def __init__(self,gp_link=None,analytical_mean=False,analytical_variance=False,beta=1.): def __init__(self,gp_link=None,beta=1.):
self.beta = beta if gp_link is None:
super(Gamma, self).__init__(gp_link,analytical_mean,analytical_variance) gp_link = link_functions.Log()
super(Gamma, self).__init__(gp_link, 'Gamma')
def _preprocess_values(self,Y): self.beta = Param('beta', beta)
return Y self.add_parameter(self.beta)
self.beta.fix()#TODO: gradients!
def pdf_link(self, link_f, y, extra_data=None): def pdf_link(self, link_f, y, Y_metadata=None):
""" """
Likelihood function given link(f) Likelihood function given link(f)
@ -37,7 +40,7 @@ class Gamma(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in poisson distribution :param Y_metadata: Y_metadata which is not used in poisson distribution
:returns: likelihood evaluated for this point :returns: likelihood evaluated for this point
:rtype: float :rtype: float
""" """
@ -47,7 +50,7 @@ class Gamma(Likelihood):
objective = (y**(alpha - 1.) * np.exp(-self.beta*y) * self.beta**alpha)/ special.gamma(alpha) objective = (y**(alpha - 1.) * np.exp(-self.beta*y) * self.beta**alpha)/ special.gamma(alpha)
return np.exp(np.sum(np.log(objective))) return np.exp(np.sum(np.log(objective)))
def logpdf_link(self, link_f, y, extra_data=None): def logpdf_link(self, link_f, y, Y_metadata=None):
""" """
Log Likelihood Function given link(f) Log Likelihood Function given link(f)
@ -59,7 +62,7 @@ class Gamma(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in poisson distribution :param Y_metadata: Y_metadata which is not used in poisson distribution
:returns: likelihood evaluated for this point :returns: likelihood evaluated for this point
:rtype: float :rtype: float
@ -71,7 +74,7 @@ class Gamma(Likelihood):
log_objective = alpha*np.log(self.beta) - np.log(special.gamma(alpha)) + (alpha - 1)*np.log(y) - self.beta*y log_objective = alpha*np.log(self.beta) - np.log(special.gamma(alpha)) + (alpha - 1)*np.log(y) - self.beta*y
return np.sum(log_objective) return np.sum(log_objective)
def dlogpdf_dlink(self, link_f, y, extra_data=None): def dlogpdf_dlink(self, link_f, y, Y_metadata=None):
""" """
Gradient of the log likelihood function at y, given link(f) w.r.t link(f) Gradient of the log likelihood function at y, given link(f) w.r.t link(f)
@ -83,7 +86,7 @@ class Gamma(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in gamma distribution :param Y_metadata: Y_metadata which is not used in gamma distribution
:returns: gradient of likelihood evaluated at points :returns: gradient of likelihood evaluated at points
:rtype: Nx1 array :rtype: Nx1 array
@ -94,7 +97,7 @@ class Gamma(Likelihood):
#return -self.gp_link.dtransf_df(gp)*self.beta*np.log(obs) + special.psi(self.gp_link.transf(gp)*self.beta) * self.gp_link.dtransf_df(gp)*self.beta #return -self.gp_link.dtransf_df(gp)*self.beta*np.log(obs) + special.psi(self.gp_link.transf(gp)*self.beta) * self.gp_link.dtransf_df(gp)*self.beta
return grad return grad
def d2logpdf_dlink2(self, link_f, y, extra_data=None): def d2logpdf_dlink2(self, link_f, y, Y_metadata=None):
""" """
Hessian at y, given link(f), w.r.t link(f) Hessian at y, given link(f), w.r.t link(f)
i.e. second derivative logpdf at y given link(f_i) and link(f_j) w.r.t link(f_i) and link(f_j) i.e. second derivative logpdf at y given link(f_i) and link(f_j) w.r.t link(f_i) and link(f_j)
@ -108,7 +111,7 @@ class Gamma(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in gamma distribution :param Y_metadata: Y_metadata which is not used in gamma distribution
:returns: Diagonal of hessian matrix (second derivative of likelihood evaluated at points f) :returns: Diagonal of hessian matrix (second derivative of likelihood evaluated at points f)
:rtype: Nx1 array :rtype: Nx1 array
@ -122,7 +125,7 @@ class Gamma(Likelihood):
#return -self.gp_link.d2transf_df2(gp)*self.beta*np.log(obs) + special.polygamma(1,self.gp_link.transf(gp)*self.beta)*(self.gp_link.dtransf_df(gp)*self.beta)**2 + special.psi(self.gp_link.transf(gp)*self.beta)*self.gp_link.d2transf_df2(gp)*self.beta #return -self.gp_link.d2transf_df2(gp)*self.beta*np.log(obs) + special.polygamma(1,self.gp_link.transf(gp)*self.beta)*(self.gp_link.dtransf_df(gp)*self.beta)**2 + special.psi(self.gp_link.transf(gp)*self.beta)*self.gp_link.d2transf_df2(gp)*self.beta
return hess return hess
def d3logpdf_dlink3(self, link_f, y, extra_data=None): def d3logpdf_dlink3(self, link_f, y, Y_metadata=None):
""" """
Third order derivative log-likelihood function at y given link(f) w.r.t link(f) Third order derivative log-likelihood function at y given link(f) w.r.t link(f)
@ -134,22 +137,10 @@ class Gamma(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in gamma distribution :param Y_metadata: Y_metadata which is not used in gamma distribution
:returns: third derivative of likelihood evaluated at points f :returns: third derivative of likelihood evaluated at points f
:rtype: Nx1 array :rtype: Nx1 array
""" """
assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape
d3lik_dlink3 = -special.polygamma(2, self.beta*link_f)*(self.beta**3) d3lik_dlink3 = -special.polygamma(2, self.beta*link_f)*(self.beta**3)
return d3lik_dlink3 return d3lik_dlink3
def _mean(self,gp):
"""
Mass (or density) function
"""
return self.gp_link.transf(gp)
def _variance(self,gp):
"""
Mass (or density) function
"""
return self.gp_link.transf(gp)/self.beta

View file

@ -35,12 +35,7 @@ class Gaussian(Likelihood):
if gp_link is None: if gp_link is None:
gp_link = link_functions.Identity() gp_link = link_functions.Identity()
if isinstance(gp_link, link_functions.Identity): assert isinstance(gp_link, link_functions.Identity), "the likelihood only implemented for the identity link"
analytical_variance = True
analytical_mean = True
else:
analytical_variance = False
analytical_mean = False
super(Gaussian, self).__init__(gp_link, name=name) super(Gaussian, self).__init__(gp_link, name=name)
@ -100,7 +95,7 @@ class Gaussian(Likelihood):
def predictive_quantiles(self, mu, var, quantiles, Y_metadata=None): def predictive_quantiles(self, mu, var, quantiles, Y_metadata=None):
return [stats.norm.ppf(q/100.)*np.sqrt(var) + mu for q in quantiles] return [stats.norm.ppf(q/100.)*np.sqrt(var) + mu for q in quantiles]
def pdf_link(self, link_f, y, extra_data=None): def pdf_link(self, link_f, y, Y_metadata=None):
""" """
Likelihood function given link(f) Likelihood function given link(f)
@ -111,14 +106,14 @@ class Gaussian(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in gaussian :param Y_metadata: Y_metadata not used in gaussian
:returns: likelihood evaluated for this point :returns: likelihood evaluated for this point
:rtype: float :rtype: float
""" """
#Assumes no covariance, exp, sum, log for numerical stability #Assumes no covariance, exp, sum, log for numerical stability
return np.exp(np.sum(np.log(stats.norm.pdf(y, link_f, np.sqrt(self.variance))))) return np.exp(np.sum(np.log(stats.norm.pdf(y, link_f, np.sqrt(self.variance)))))
def logpdf_link(self, link_f, y, extra_data=None): def logpdf_link(self, link_f, y, Y_metadata=None):
""" """
Log likelihood function given link(f) Log likelihood function given link(f)
@ -129,7 +124,7 @@ class Gaussian(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in gaussian :param Y_metadata: Y_metadata not used in gaussian
:returns: log likelihood evaluated for this point :returns: log likelihood evaluated for this point
:rtype: float :rtype: float
""" """
@ -139,7 +134,7 @@ class Gaussian(Likelihood):
return -0.5*(np.sum((y-link_f)**2/self.variance) + ln_det_cov + N*np.log(2.*np.pi)) return -0.5*(np.sum((y-link_f)**2/self.variance) + ln_det_cov + N*np.log(2.*np.pi))
def dlogpdf_dlink(self, link_f, y, extra_data=None): def dlogpdf_dlink(self, link_f, y, Y_metadata=None):
""" """
Gradient of the pdf at y, given link(f) w.r.t link(f) Gradient of the pdf at y, given link(f) w.r.t link(f)
@ -150,7 +145,7 @@ class Gaussian(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in gaussian :param Y_metadata: Y_metadata not used in gaussian
:returns: gradient of log likelihood evaluated at points link(f) :returns: gradient of log likelihood evaluated at points link(f)
:rtype: Nx1 array :rtype: Nx1 array
""" """
@ -159,7 +154,7 @@ class Gaussian(Likelihood):
grad = s2_i*y - s2_i*link_f grad = s2_i*y - s2_i*link_f
return grad return grad
def d2logpdf_dlink2(self, link_f, y, extra_data=None): def d2logpdf_dlink2(self, link_f, y, Y_metadata=None):
""" """
Hessian at y, given link_f, w.r.t link_f. Hessian at y, given link_f, w.r.t link_f.
i.e. second derivative logpdf at y given link(f_i) link(f_j) w.r.t link(f_i) and link(f_j) i.e. second derivative logpdf at y given link(f_i) link(f_j) w.r.t link(f_i) and link(f_j)
@ -173,7 +168,7 @@ class Gaussian(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in gaussian :param Y_metadata: Y_metadata not used in gaussian
:returns: Diagonal of log hessian matrix (second derivative of log likelihood evaluated at points link(f)) :returns: Diagonal of log hessian matrix (second derivative of log likelihood evaluated at points link(f))
:rtype: Nx1 array :rtype: Nx1 array
@ -186,7 +181,7 @@ class Gaussian(Likelihood):
hess = -(1.0/self.variance)*np.ones((N, 1)) hess = -(1.0/self.variance)*np.ones((N, 1))
return hess return hess
def d3logpdf_dlink3(self, link_f, y, extra_data=None): def d3logpdf_dlink3(self, link_f, y, Y_metadata=None):
""" """
Third order derivative log-likelihood function at y given link(f) w.r.t link(f) Third order derivative log-likelihood function at y given link(f) w.r.t link(f)
@ -197,7 +192,7 @@ class Gaussian(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in gaussian :param Y_metadata: Y_metadata not used in gaussian
:returns: third derivative of log likelihood evaluated at points link(f) :returns: third derivative of log likelihood evaluated at points link(f)
:rtype: Nx1 array :rtype: Nx1 array
""" """
@ -206,7 +201,7 @@ class Gaussian(Likelihood):
d3logpdf_dlink3 = np.zeros((N,1)) d3logpdf_dlink3 = np.zeros((N,1))
return d3logpdf_dlink3 return d3logpdf_dlink3
def dlogpdf_link_dvar(self, link_f, y, extra_data=None): def dlogpdf_link_dvar(self, link_f, y, Y_metadata=None):
""" """
Gradient of the log-likelihood function at y given link(f), w.r.t variance parameter (noise_variance) Gradient of the log-likelihood function at y given link(f), w.r.t variance parameter (noise_variance)
@ -217,7 +212,7 @@ class Gaussian(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in gaussian :param Y_metadata: Y_metadata not used in gaussian
:returns: derivative of log likelihood evaluated at points link(f) w.r.t variance parameter :returns: derivative of log likelihood evaluated at points link(f) w.r.t variance parameter
:rtype: float :rtype: float
""" """
@ -228,7 +223,7 @@ class Gaussian(Likelihood):
dlik_dsigma = -0.5*N/self.variance + 0.5*s_4*np.sum(np.square(e)) dlik_dsigma = -0.5*N/self.variance + 0.5*s_4*np.sum(np.square(e))
return np.sum(dlik_dsigma) # Sure about this sum? return np.sum(dlik_dsigma) # Sure about this sum?
def dlogpdf_dlink_dvar(self, link_f, y, extra_data=None): def dlogpdf_dlink_dvar(self, link_f, y, Y_metadata=None):
""" """
Derivative of the dlogpdf_dlink w.r.t variance parameter (noise_variance) Derivative of the dlogpdf_dlink w.r.t variance parameter (noise_variance)
@ -239,7 +234,7 @@ class Gaussian(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in gaussian :param Y_metadata: Y_metadata not used in gaussian
:returns: derivative of log likelihood evaluated at points link(f) w.r.t variance parameter :returns: derivative of log likelihood evaluated at points link(f) w.r.t variance parameter
:rtype: Nx1 array :rtype: Nx1 array
""" """
@ -248,7 +243,7 @@ class Gaussian(Likelihood):
dlik_grad_dsigma = -s_4*y + s_4*link_f dlik_grad_dsigma = -s_4*y + s_4*link_f
return dlik_grad_dsigma return dlik_grad_dsigma
def d2logpdf_dlink2_dvar(self, link_f, y, extra_data=None): def d2logpdf_dlink2_dvar(self, link_f, y, Y_metadata=None):
""" """
Gradient of the hessian (d2logpdf_dlink2) w.r.t variance parameter (noise_variance) Gradient of the hessian (d2logpdf_dlink2) w.r.t variance parameter (noise_variance)
@ -259,7 +254,7 @@ class Gaussian(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data not used in gaussian :param Y_metadata: Y_metadata not used in gaussian
:returns: derivative of log hessian evaluated at points link(f_i) and link(f_j) w.r.t variance parameter :returns: derivative of log hessian evaluated at points link(f_i) and link(f_j) w.r.t variance parameter
:rtype: Nx1 array :rtype: Nx1 array
""" """
@ -269,16 +264,16 @@ class Gaussian(Likelihood):
d2logpdf_dlink2_dvar = np.ones((N,1))*s_4 d2logpdf_dlink2_dvar = np.ones((N,1))*s_4
return d2logpdf_dlink2_dvar return d2logpdf_dlink2_dvar
def dlogpdf_link_dtheta(self, f, y, extra_data=None): def dlogpdf_link_dtheta(self, f, y, Y_metadata=None):
dlogpdf_dvar = self.dlogpdf_link_dvar(f, y, extra_data=extra_data) dlogpdf_dvar = self.dlogpdf_link_dvar(f, y, Y_metadata=Y_metadata)
return np.asarray([[dlogpdf_dvar]]) return np.asarray([[dlogpdf_dvar]])
def dlogpdf_dlink_dtheta(self, f, y, extra_data=None): def dlogpdf_dlink_dtheta(self, f, y, Y_metadata=None):
dlogpdf_dlink_dvar = self.dlogpdf_dlink_dvar(f, y, extra_data=extra_data) dlogpdf_dlink_dvar = self.dlogpdf_dlink_dvar(f, y, Y_metadata=Y_metadata)
return dlogpdf_dlink_dvar return dlogpdf_dlink_dvar
def d2logpdf_dlink2_dtheta(self, f, y, extra_data=None): def d2logpdf_dlink2_dtheta(self, f, y, Y_metadata=None):
d2logpdf_dlink2_dvar = self.d2logpdf_dlink2_dvar(f, y, extra_data=extra_data) d2logpdf_dlink2_dvar = self.d2logpdf_dlink2_dvar(f, y, Y_metadata=Y_metadata)
return d2logpdf_dlink2_dvar return d2logpdf_dlink2_dvar
def _mean(self, gp): def _mean(self, gp):

View file

@ -153,6 +153,10 @@ class Likelihood(Parameterized):
return mean 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, Y_metadata=None): def predictive_variance(self, mu,variance, predictive_mean=None, Y_metadata=None):
""" """
Numerical approximation to the predictive variance: V(Y_star) Numerical approximation to the predictive variance: V(Y_star)
@ -204,31 +208,31 @@ class Likelihood(Parameterized):
# V(Y_star) = E[ V(Y_star|f_star) ] + E(Y_star**2|f_star) - E[Y_star|f_star]**2 # V(Y_star) = E[ V(Y_star|f_star) ] + E(Y_star**2|f_star) - E[Y_star|f_star]**2
return exp_var + var_exp return exp_var + var_exp
def pdf_link(self, link_f, y, extra_data=None): def pdf_link(self, link_f, y, Y_metadata=None):
raise NotImplementedError raise NotImplementedError
def logpdf_link(self, link_f, y, extra_data=None): def logpdf_link(self, link_f, y, Y_metadata=None):
raise NotImplementedError raise NotImplementedError
def dlogpdf_dlink(self, link_f, y, extra_data=None): def dlogpdf_dlink(self, link_f, y, Y_metadata=None):
raise NotImplementedError raise NotImplementedError
def d2logpdf_dlink2(self, link_f, y, extra_data=None): def d2logpdf_dlink2(self, link_f, y, Y_metadata=None):
raise NotImplementedError raise NotImplementedError
def d3logpdf_dlink3(self, link_f, y, extra_data=None): def d3logpdf_dlink3(self, link_f, y, Y_metadata=None):
raise NotImplementedError raise NotImplementedError
def dlogpdf_link_dtheta(self, link_f, y, extra_data=None): def dlogpdf_link_dtheta(self, link_f, y, Y_metadata=None):
raise NotImplementedError raise NotImplementedError
def dlogpdf_dlink_dtheta(self, link_f, y, extra_data=None): def dlogpdf_dlink_dtheta(self, link_f, y, Y_metadata=None):
raise NotImplementedError raise NotImplementedError
def d2logpdf_dlink2_dtheta(self, link_f, y, extra_data=None): def d2logpdf_dlink2_dtheta(self, link_f, y, Y_metadata=None):
raise NotImplementedError raise NotImplementedError
def pdf(self, f, y, extra_data=None): def pdf(self, f, y, Y_metadata=None):
""" """
Evaluates the link function link(f) then computes the likelihood (pdf) using it Evaluates the link function link(f) then computes the likelihood (pdf) using it
@ -239,14 +243,14 @@ class Likelihood(Parameterized):
:type f: Nx1 array :type f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution - not used :param Y_metadata: Y_metadata which is not used in student t distribution - not used
:returns: likelihood evaluated for this point :returns: likelihood evaluated for this point
:rtype: float :rtype: float
""" """
link_f = self.gp_link.transf(f) link_f = self.gp_link.transf(f)
return self.pdf_link(link_f, y, extra_data=extra_data) return self.pdf_link(link_f, y, Y_metadata=Y_metadata)
def logpdf(self, f, y, extra_data=None): def logpdf(self, f, y, Y_metadata=None):
""" """
Evaluates the link function link(f) then computes the log likelihood (log pdf) using it Evaluates the link function link(f) then computes the log likelihood (log pdf) using it
@ -257,14 +261,14 @@ class Likelihood(Parameterized):
:type f: Nx1 array :type f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution - not used :param Y_metadata: Y_metadata which is not used in student t distribution - not used
:returns: log likelihood evaluated for this point :returns: log likelihood evaluated for this point
:rtype: float :rtype: float
""" """
link_f = self.gp_link.transf(f) link_f = self.gp_link.transf(f)
return self.logpdf_link(link_f, y, extra_data=extra_data) return self.logpdf_link(link_f, y, Y_metadata=Y_metadata)
def dlogpdf_df(self, f, y, extra_data=None): def dlogpdf_df(self, f, y, Y_metadata=None):
""" """
Evaluates the link function link(f) then computes the derivative of log likelihood using it Evaluates the link function link(f) then computes the derivative of log likelihood using it
Uses the Faa di Bruno's formula for the chain rule Uses the Faa di Bruno's formula for the chain rule
@ -276,16 +280,16 @@ class Likelihood(Parameterized):
:type f: Nx1 array :type f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution - not used :param Y_metadata: Y_metadata which is not used in student t distribution - not used
:returns: derivative of log likelihood evaluated for this point :returns: derivative of log likelihood evaluated for this point
:rtype: 1xN array :rtype: 1xN array
""" """
link_f = self.gp_link.transf(f) link_f = self.gp_link.transf(f)
dlogpdf_dlink = self.dlogpdf_dlink(link_f, y, extra_data=extra_data) dlogpdf_dlink = self.dlogpdf_dlink(link_f, y, Y_metadata=Y_metadata)
dlink_df = self.gp_link.dtransf_df(f) dlink_df = self.gp_link.dtransf_df(f)
return chain_1(dlogpdf_dlink, dlink_df) return chain_1(dlogpdf_dlink, dlink_df)
def d2logpdf_df2(self, f, y, extra_data=None): def d2logpdf_df2(self, f, y, Y_metadata=None):
""" """
Evaluates the link function link(f) then computes the second derivative of log likelihood using it Evaluates the link function link(f) then computes the second derivative of log likelihood using it
Uses the Faa di Bruno's formula for the chain rule Uses the Faa di Bruno's formula for the chain rule
@ -297,18 +301,18 @@ class Likelihood(Parameterized):
:type f: Nx1 array :type f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution - not used :param Y_metadata: Y_metadata which is not used in student t distribution - not used
:returns: second derivative of log likelihood evaluated for this point (diagonal only) :returns: second derivative of log likelihood evaluated for this point (diagonal only)
:rtype: 1xN array :rtype: 1xN array
""" """
link_f = self.gp_link.transf(f) link_f = self.gp_link.transf(f)
d2logpdf_dlink2 = self.d2logpdf_dlink2(link_f, y, extra_data=extra_data) d2logpdf_dlink2 = self.d2logpdf_dlink2(link_f, y, Y_metadata=Y_metadata)
dlink_df = self.gp_link.dtransf_df(f) dlink_df = self.gp_link.dtransf_df(f)
dlogpdf_dlink = self.dlogpdf_dlink(link_f, y, extra_data=extra_data) dlogpdf_dlink = self.dlogpdf_dlink(link_f, y, Y_metadata=Y_metadata)
d2link_df2 = self.gp_link.d2transf_df2(f) d2link_df2 = self.gp_link.d2transf_df2(f)
return chain_2(d2logpdf_dlink2, dlink_df, dlogpdf_dlink, d2link_df2) return chain_2(d2logpdf_dlink2, dlink_df, dlogpdf_dlink, d2link_df2)
def d3logpdf_df3(self, f, y, extra_data=None): def d3logpdf_df3(self, f, y, Y_metadata=None):
""" """
Evaluates the link function link(f) then computes the third derivative of log likelihood using it Evaluates the link function link(f) then computes the third derivative of log likelihood using it
Uses the Faa di Bruno's formula for the chain rule Uses the Faa di Bruno's formula for the chain rule
@ -320,44 +324,44 @@ class Likelihood(Parameterized):
:type f: Nx1 array :type f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution - not used :param Y_metadata: Y_metadata which is not used in student t distribution - not used
:returns: third derivative of log likelihood evaluated for this point :returns: third derivative of log likelihood evaluated for this point
:rtype: float :rtype: float
""" """
link_f = self.gp_link.transf(f) link_f = self.gp_link.transf(f)
d3logpdf_dlink3 = self.d3logpdf_dlink3(link_f, y, extra_data=extra_data) d3logpdf_dlink3 = self.d3logpdf_dlink3(link_f, y, Y_metadata=Y_metadata)
dlink_df = self.gp_link.dtransf_df(f) dlink_df = self.gp_link.dtransf_df(f)
d2logpdf_dlink2 = self.d2logpdf_dlink2(link_f, y, extra_data=extra_data) d2logpdf_dlink2 = self.d2logpdf_dlink2(link_f, y, Y_metadata=Y_metadata)
d2link_df2 = self.gp_link.d2transf_df2(f) d2link_df2 = self.gp_link.d2transf_df2(f)
dlogpdf_dlink = self.dlogpdf_dlink(link_f, y, extra_data=extra_data) dlogpdf_dlink = self.dlogpdf_dlink(link_f, y, Y_metadata=Y_metadata)
d3link_df3 = self.gp_link.d3transf_df3(f) d3link_df3 = self.gp_link.d3transf_df3(f)
return chain_3(d3logpdf_dlink3, dlink_df, d2logpdf_dlink2, d2link_df2, dlogpdf_dlink, d3link_df3) return chain_3(d3logpdf_dlink3, dlink_df, d2logpdf_dlink2, d2link_df2, dlogpdf_dlink, d3link_df3)
def dlogpdf_dtheta(self, f, y, extra_data=None): def dlogpdf_dtheta(self, f, y, Y_metadata=None):
""" """
TODO: Doc strings TODO: Doc strings
""" """
if self.size > 0: if self.size > 0:
link_f = self.gp_link.transf(f) link_f = self.gp_link.transf(f)
return self.dlogpdf_link_dtheta(link_f, y, extra_data=extra_data) return self.dlogpdf_link_dtheta(link_f, y, Y_metadata=Y_metadata)
else: else:
#Is no parameters so return an empty array for its derivatives #Is no parameters so return an empty array for its derivatives
return np.zeros([1, 0]) return np.zeros([1, 0])
def dlogpdf_df_dtheta(self, f, y, extra_data=None): def dlogpdf_df_dtheta(self, f, y, Y_metadata=None):
""" """
TODO: Doc strings TODO: Doc strings
""" """
if self.size > 0: if self.size > 0:
link_f = self.gp_link.transf(f) link_f = self.gp_link.transf(f)
dlink_df = self.gp_link.dtransf_df(f) dlink_df = self.gp_link.dtransf_df(f)
dlogpdf_dlink_dtheta = self.dlogpdf_dlink_dtheta(link_f, y, extra_data=extra_data) dlogpdf_dlink_dtheta = self.dlogpdf_dlink_dtheta(link_f, y, Y_metadata=Y_metadata)
return chain_1(dlogpdf_dlink_dtheta, dlink_df) return chain_1(dlogpdf_dlink_dtheta, dlink_df)
else: else:
#Is no parameters so return an empty array for its derivatives #Is no parameters so return an empty array for its derivatives
return np.zeros([f.shape[0], 0]) return np.zeros([f.shape[0], 0])
def d2logpdf_df2_dtheta(self, f, y, extra_data=None): def d2logpdf_df2_dtheta(self, f, y, Y_metadata=None):
""" """
TODO: Doc strings TODO: Doc strings
""" """
@ -365,17 +369,17 @@ class Likelihood(Parameterized):
link_f = self.gp_link.transf(f) link_f = self.gp_link.transf(f)
dlink_df = self.gp_link.dtransf_df(f) dlink_df = self.gp_link.dtransf_df(f)
d2link_df2 = self.gp_link.d2transf_df2(f) d2link_df2 = self.gp_link.d2transf_df2(f)
d2logpdf_dlink2_dtheta = self.d2logpdf_dlink2_dtheta(link_f, y, extra_data=extra_data) d2logpdf_dlink2_dtheta = self.d2logpdf_dlink2_dtheta(link_f, y, Y_metadata=Y_metadata)
dlogpdf_dlink_dtheta = self.dlogpdf_dlink_dtheta(link_f, y, extra_data=extra_data) dlogpdf_dlink_dtheta = self.dlogpdf_dlink_dtheta(link_f, y, Y_metadata=Y_metadata)
return chain_2(d2logpdf_dlink2_dtheta, dlink_df, dlogpdf_dlink_dtheta, d2link_df2) return chain_2(d2logpdf_dlink2_dtheta, dlink_df, dlogpdf_dlink_dtheta, d2link_df2)
else: else:
#Is no parameters so return an empty array for its derivatives #Is no parameters so return an empty array for its derivatives
return np.zeros([f.shape[0], 0]) return np.zeros([f.shape[0], 0])
def _laplace_gradients(self, f, y, extra_data=None): def _laplace_gradients(self, f, y, Y_metadata=None):
dlogpdf_dtheta = self.dlogpdf_dtheta(f, y, extra_data=extra_data) dlogpdf_dtheta = self.dlogpdf_dtheta(f, y, Y_metadata=Y_metadata)
dlogpdf_df_dtheta = self.dlogpdf_df_dtheta(f, y, extra_data=extra_data) dlogpdf_df_dtheta = self.dlogpdf_df_dtheta(f, y, Y_metadata=Y_metadata)
d2logpdf_df2_dtheta = self.d2logpdf_df2_dtheta(f, y, extra_data=extra_data) d2logpdf_df2_dtheta = self.d2logpdf_df2_dtheta(f, y, Y_metadata=Y_metadata)
#Parameters are stacked vertically. Must be listed in same order as 'get_param_names' #Parameters are stacked vertically. Must be listed in same order as 'get_param_names'
# ensure we have gradients for every parameter we want to optimize # ensure we have gradients for every parameter we want to optimize
@ -390,7 +394,7 @@ class Likelihood(Parameterized):
def predictive_values(self, mu, var, full_cov=False, Y_metadata=None): def predictive_values(self, mu, var, full_cov=False, Y_metadata=None):
""" """
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 mu: mean of the latent variable, f, of posterior
:param var: variance of the latent variable, f, of posterior :param var: variance of the latent variable, f, of posterior

View file

@ -24,10 +24,11 @@ class MixedNoise(Likelihood):
variance = np.zeros(ind.size) variance = np.zeros(ind.size)
for lik, j in zip(self.likelihoods_list, range(len(self.likelihoods_list))): for lik, j in zip(self.likelihoods_list, range(len(self.likelihoods_list))):
variance[ind==j] = lik.variance variance[ind==j] = lik.variance
return variance[:,None] return variance
def betaY(self,Y,Y_metadata): def betaY(self,Y,Y_metadata):
return Y/self.gaussian_variance(Y_metadata=Y_metadata) #TODO not here.
return Y/self.gaussian_variance(Y_metadata=Y_metadata)[:,None]
def update_gradients(self, gradients): def update_gradients(self, gradients):
self.gradient = gradients self.gradient = gradients
@ -60,10 +61,6 @@ class MixedNoise(Likelihood):
Q[ind==j,:] = np.hstack(q) Q[ind==j,:] = np.hstack(q)
return [q[:,None] for q in Q.T] return [q[:,None] for q in Q.T]
def covariance_matrix(self, Y, Y_metadata):
#TODO make more general, to allow non-gaussian likelihoods
return np.diag(self.gaussian_variance(Y_metadata).flatten())
def samples(self, gp, Y_metadata): def samples(self, gp, Y_metadata):
""" """
Returns a set of samples of observations based on a given value of the latent variable. Returns a set of samples of observations based on a given value of the latent variable.

View file

@ -25,10 +25,13 @@ class Poisson(Likelihood):
super(Poisson, self).__init__(gp_link, name='Poisson') super(Poisson, self).__init__(gp_link, name='Poisson')
def _preprocess_values(self,Y): def _conditional_mean(self, f):
return Y """
the expected value of y given a value of f
"""
return self.gp_link.transf(gp)
def pdf_link(self, link_f, y, extra_data=None): def pdf_link(self, link_f, y, Y_metadata=None):
""" """
Likelihood function given link(f) Likelihood function given link(f)
@ -39,14 +42,14 @@ class Poisson(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in poisson distribution :param Y_metadata: Y_metadata which is not used in poisson distribution
:returns: likelihood evaluated for this point :returns: likelihood evaluated for this point
:rtype: float :rtype: float
""" """
assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape
return np.prod(stats.poisson.pmf(y,link_f)) return np.prod(stats.poisson.pmf(y,link_f))
def logpdf_link(self, link_f, y, extra_data=None): def logpdf_link(self, link_f, y, Y_metadata=None):
""" """
Log Likelihood Function given link(f) Log Likelihood Function given link(f)
@ -57,7 +60,7 @@ class Poisson(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in poisson distribution :param Y_metadata: Y_metadata which is not used in poisson distribution
:returns: likelihood evaluated for this point :returns: likelihood evaluated for this point
:rtype: float :rtype: float
@ -65,7 +68,7 @@ class Poisson(Likelihood):
assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape
return np.sum(-link_f + y*np.log(link_f) - special.gammaln(y+1)) return np.sum(-link_f + y*np.log(link_f) - special.gammaln(y+1))
def dlogpdf_dlink(self, link_f, y, extra_data=None): def dlogpdf_dlink(self, link_f, y, Y_metadata=None):
""" """
Gradient of the log likelihood function at y, given link(f) w.r.t link(f) Gradient of the log likelihood function at y, given link(f) w.r.t link(f)
@ -76,7 +79,7 @@ class Poisson(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in poisson distribution :param Y_metadata: Y_metadata which is not used in poisson distribution
:returns: gradient of likelihood evaluated at points :returns: gradient of likelihood evaluated at points
:rtype: Nx1 array :rtype: Nx1 array
@ -84,7 +87,7 @@ class Poisson(Likelihood):
assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape assert np.atleast_1d(link_f).shape == np.atleast_1d(y).shape
return y/link_f - 1 return y/link_f - 1
def d2logpdf_dlink2(self, link_f, y, extra_data=None): def d2logpdf_dlink2(self, link_f, y, Y_metadata=None):
""" """
Hessian at y, given link(f), w.r.t link(f) Hessian at y, given link(f), w.r.t link(f)
i.e. second derivative logpdf at y given link(f_i) and link(f_j) w.r.t link(f_i) and link(f_j) i.e. second derivative logpdf at y given link(f_i) and link(f_j) w.r.t link(f_i) and link(f_j)
@ -97,7 +100,7 @@ class Poisson(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in poisson distribution :param Y_metadata: Y_metadata which is not used in poisson distribution
:returns: Diagonal of hessian matrix (second derivative of likelihood evaluated at points f) :returns: Diagonal of hessian matrix (second derivative of likelihood evaluated at points f)
:rtype: Nx1 array :rtype: Nx1 array
@ -112,7 +115,7 @@ class Poisson(Likelihood):
#transf = self.gp_link.transf(gp) #transf = self.gp_link.transf(gp)
#return obs * ((self.gp_link.dtransf_df(gp)/transf)**2 - d2_df/transf) + d2_df #return obs * ((self.gp_link.dtransf_df(gp)/transf)**2 - d2_df/transf) + d2_df
def d3logpdf_dlink3(self, link_f, y, extra_data=None): def d3logpdf_dlink3(self, link_f, y, Y_metadata=None):
""" """
Third order derivative log-likelihood function at y given link(f) w.r.t link(f) Third order derivative log-likelihood function at y given link(f) w.r.t link(f)
@ -123,7 +126,7 @@ class Poisson(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in poisson distribution :param Y_metadata: Y_metadata which is not used in poisson distribution
:returns: third derivative of likelihood evaluated at points f :returns: third derivative of likelihood evaluated at points f
:rtype: Nx1 array :rtype: Nx1 array
""" """

View file

@ -43,10 +43,10 @@ class StudentT(Likelihood):
Pull out the gradients, be careful as the order must match the order Pull out the gradients, be careful as the order must match the order
in which the parameters are added in which the parameters are added
""" """
self.sigma2.gradient = grads[0] self.sigma2.gradient = derivatives[0]
self.v.gradient = grads[1] self.v.gradient = derivatives[1]
def pdf_link(self, link_f, y, extra_data=None): def pdf_link(self, link_f, y, Y_metadata=None):
""" """
Likelihood function given link(f) Likelihood function given link(f)
@ -57,7 +57,7 @@ class StudentT(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution :param Y_metadata: Y_metadata which is not used in student t distribution
:returns: likelihood evaluated for this point :returns: likelihood evaluated for this point
:rtype: float :rtype: float
""" """
@ -70,7 +70,7 @@ class StudentT(Likelihood):
) )
return np.prod(objective) return np.prod(objective)
def logpdf_link(self, link_f, y, extra_data=None): def logpdf_link(self, link_f, y, Y_metadata=None):
""" """
Log Likelihood Function given link(f) Log Likelihood Function given link(f)
@ -81,7 +81,7 @@ class StudentT(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution :param Y_metadata: Y_metadata which is not used in student t distribution
:returns: likelihood evaluated for this point :returns: likelihood evaluated for this point
:rtype: float :rtype: float
@ -99,7 +99,7 @@ class StudentT(Likelihood):
) )
return np.sum(objective) return np.sum(objective)
def dlogpdf_dlink(self, link_f, y, extra_data=None): def dlogpdf_dlink(self, link_f, y, Y_metadata=None):
""" """
Gradient of the log likelihood function at y, given link(f) w.r.t link(f) Gradient of the log likelihood function at y, given link(f) w.r.t link(f)
@ -110,7 +110,7 @@ class StudentT(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution :param Y_metadata: Y_metadata which is not used in student t distribution
:returns: gradient of likelihood evaluated at points :returns: gradient of likelihood evaluated at points
:rtype: Nx1 array :rtype: Nx1 array
@ -120,7 +120,7 @@ class StudentT(Likelihood):
grad = ((self.v + 1) * e) / (self.v * self.sigma2 + (e**2)) grad = ((self.v + 1) * e) / (self.v * self.sigma2 + (e**2))
return grad return grad
def d2logpdf_dlink2(self, link_f, y, extra_data=None): def d2logpdf_dlink2(self, link_f, y, Y_metadata=None):
""" """
Hessian at y, given link(f), w.r.t link(f) Hessian at y, given link(f), w.r.t link(f)
i.e. second derivative logpdf at y given link(f_i) and link(f_j) w.r.t link(f_i) and link(f_j) i.e. second derivative logpdf at y given link(f_i) and link(f_j) w.r.t link(f_i) and link(f_j)
@ -133,7 +133,7 @@ class StudentT(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution :param Y_metadata: Y_metadata which is not used in student t distribution
:returns: Diagonal of hessian matrix (second derivative of likelihood evaluated at points f) :returns: Diagonal of hessian matrix (second derivative of likelihood evaluated at points f)
:rtype: Nx1 array :rtype: Nx1 array
@ -146,7 +146,7 @@ class StudentT(Likelihood):
hess = ((self.v + 1)*(e**2 - self.v*self.sigma2)) / ((self.sigma2*self.v + e**2)**2) hess = ((self.v + 1)*(e**2 - self.v*self.sigma2)) / ((self.sigma2*self.v + e**2)**2)
return hess return hess
def d3logpdf_dlink3(self, link_f, y, extra_data=None): def d3logpdf_dlink3(self, link_f, y, Y_metadata=None):
""" """
Third order derivative log-likelihood function at y given link(f) w.r.t link(f) Third order derivative log-likelihood function at y given link(f) w.r.t link(f)
@ -157,7 +157,7 @@ class StudentT(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution :param Y_metadata: Y_metadata which is not used in student t distribution
:returns: third derivative of likelihood evaluated at points f :returns: third derivative of likelihood evaluated at points f
:rtype: Nx1 array :rtype: Nx1 array
""" """
@ -168,7 +168,7 @@ class StudentT(Likelihood):
) )
return d3lik_dlink3 return d3lik_dlink3
def dlogpdf_link_dvar(self, link_f, y, extra_data=None): def dlogpdf_link_dvar(self, link_f, y, Y_metadata=None):
""" """
Gradient of the log-likelihood function at y given f, w.r.t variance parameter (t_noise) Gradient of the log-likelihood function at y given f, w.r.t variance parameter (t_noise)
@ -179,7 +179,7 @@ class StudentT(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution :param Y_metadata: Y_metadata which is not used in student t distribution
:returns: derivative of likelihood evaluated at points f w.r.t variance parameter :returns: derivative of likelihood evaluated at points f w.r.t variance parameter
:rtype: float :rtype: float
""" """
@ -188,7 +188,7 @@ class StudentT(Likelihood):
dlogpdf_dvar = self.v*(e**2 - self.sigma2)/(2*self.sigma2*(self.sigma2*self.v + e**2)) dlogpdf_dvar = self.v*(e**2 - self.sigma2)/(2*self.sigma2*(self.sigma2*self.v + e**2))
return np.sum(dlogpdf_dvar) return np.sum(dlogpdf_dvar)
def dlogpdf_dlink_dvar(self, link_f, y, extra_data=None): def dlogpdf_dlink_dvar(self, link_f, y, Y_metadata=None):
""" """
Derivative of the dlogpdf_dlink w.r.t variance parameter (t_noise) Derivative of the dlogpdf_dlink w.r.t variance parameter (t_noise)
@ -199,7 +199,7 @@ class StudentT(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution :param Y_metadata: Y_metadata which is not used in student t distribution
:returns: derivative of likelihood evaluated at points f w.r.t variance parameter :returns: derivative of likelihood evaluated at points f w.r.t variance parameter
:rtype: Nx1 array :rtype: Nx1 array
""" """
@ -208,7 +208,7 @@ class StudentT(Likelihood):
dlogpdf_dlink_dvar = (self.v*(self.v+1)*(-e))/((self.sigma2*self.v + e**2)**2) dlogpdf_dlink_dvar = (self.v*(self.v+1)*(-e))/((self.sigma2*self.v + e**2)**2)
return dlogpdf_dlink_dvar return dlogpdf_dlink_dvar
def d2logpdf_dlink2_dvar(self, link_f, y, extra_data=None): def d2logpdf_dlink2_dvar(self, link_f, y, Y_metadata=None):
""" """
Gradient of the hessian (d2logpdf_dlink2) w.r.t variance parameter (t_noise) Gradient of the hessian (d2logpdf_dlink2) w.r.t variance parameter (t_noise)
@ -219,7 +219,7 @@ class StudentT(Likelihood):
:type link_f: Nx1 array :type link_f: Nx1 array
:param y: data :param y: data
:type y: Nx1 array :type y: Nx1 array
:param extra_data: extra_data which is not used in student t distribution :param Y_metadata: Y_metadata which is not used in student t distribution
:returns: derivative of hessian evaluated at points f and f_j w.r.t variance parameter :returns: derivative of hessian evaluated at points f and f_j w.r.t variance parameter
:rtype: Nx1 array :rtype: Nx1 array
""" """
@ -230,25 +230,22 @@ class StudentT(Likelihood):
) )
return d2logpdf_dlink2_dvar return d2logpdf_dlink2_dvar
def dlogpdf_link_dtheta(self, f, y, extra_data=None): def dlogpdf_link_dtheta(self, f, y, Y_metadata=None):
dlogpdf_dvar = self.dlogpdf_link_dvar(f, y, extra_data=extra_data) dlogpdf_dvar = self.dlogpdf_link_dvar(f, y, Y_metadata=Y_metadata)
dlogpdf_dv = np.zeros_like(dlogpdf_dvar) #FIXME: Not done yet dlogpdf_dv = np.zeros_like(dlogpdf_dvar) #FIXME: Not done yet
return np.hstack((dlogpdf_dvar, dlogpdf_dv)) return np.hstack((dlogpdf_dvar, dlogpdf_dv))
def dlogpdf_dlink_dtheta(self, f, y, extra_data=None): def dlogpdf_dlink_dtheta(self, f, y, Y_metadata=None):
dlogpdf_dlink_dvar = self.dlogpdf_dlink_dvar(f, y, extra_data=extra_data) dlogpdf_dlink_dvar = self.dlogpdf_dlink_dvar(f, y, Y_metadata=Y_metadata)
dlogpdf_dlink_dv = np.zeros_like(dlogpdf_dlink_dvar) #FIXME: Not done yet dlogpdf_dlink_dv = np.zeros_like(dlogpdf_dlink_dvar) #FIXME: Not done yet
return np.hstack((dlogpdf_dlink_dvar, dlogpdf_dlink_dv)) return np.hstack((dlogpdf_dlink_dvar, dlogpdf_dlink_dv))
def d2logpdf_dlink2_dtheta(self, f, y, extra_data=None): def d2logpdf_dlink2_dtheta(self, f, y, Y_metadata=None):
d2logpdf_dlink2_dvar = self.d2logpdf_dlink2_dvar(f, y, extra_data=extra_data) d2logpdf_dlink2_dvar = self.d2logpdf_dlink2_dvar(f, y, Y_metadata=Y_metadata)
d2logpdf_dlink2_dv = np.zeros_like(d2logpdf_dlink2_dvar) #FIXME: Not done yet d2logpdf_dlink2_dv = np.zeros_like(d2logpdf_dlink2_dvar) #FIXME: Not done yet
return np.hstack((d2logpdf_dlink2_dvar, d2logpdf_dlink2_dv)) return np.hstack((d2logpdf_dlink2_dvar, d2logpdf_dlink2_dv))
def predictive_mean(self, mu, sigma, Y_metadata=None): def predictive_mean(self, mu, sigma, Y_metadata=None):
"""
Compute mean of the prediction
"""
return self.gp_link.transf(mu) # only true in link is monotoci, which it is. return self.gp_link.transf(mu) # only true in link is monotoci, which it is.
def predictive_variance(self, mu,variance, predictive_mean=None, Y_metadata=None): def predictive_variance(self, mu,variance, predictive_mean=None, Y_metadata=None):

View file

@ -50,7 +50,7 @@ class SparseGPCoregionalizedRegression(SparseGP):
#Inducing inputs list #Inducing inputs list
if len(Z_list): if len(Z_list):
assert len(Z_list) == self.output_dim, 'Number of outputs do not match length of inducing inputs list.' assert len(Z_list) == Ny, 'Number of outputs do not match length of inducing inputs list.'
else: else:
if isinstance(num_inducing,np.int): if isinstance(num_inducing,np.int):
num_inducing = [num_inducing] * Ny num_inducing = [num_inducing] * Ny

View file

@ -255,25 +255,25 @@ class TestNoiseModels(object):
"Y": self.binary_Y, "Y": self.binary_Y,
"ep": False # FIXME: Should be True when we have it working again "ep": False # FIXME: Should be True when we have it working again
}, },
#"Exponential_default": { "Exponential_default": {
#"model": GPy.likelihoods.exponential(), "model": GPy.likelihoods.Exponential(),
#"link_f_constraints": [constrain_positive], "link_f_constraints": [constrain_positive],
#"Y": self.positive_Y, "Y": self.positive_Y,
#"laplace": True, "laplace": True,
#}, },
#"Poisson_default": { "Poisson_default": {
#"model": GPy.likelihoods.poisson(), "model": GPy.likelihoods.Poisson(),
#"link_f_constraints": [constrain_positive], "link_f_constraints": [constrain_positive],
#"Y": self.integer_Y, "Y": self.integer_Y,
#"laplace": True, "laplace": True,
#"ep": False #Should work though... "ep": False #Should work though...
#}, },
#"Gamma_default": { "Gamma_default": {
#"model": GPy.likelihoods.gamma(), "model": GPy.likelihoods.Gamma(),
#"link_f_constraints": [constrain_positive], "link_f_constraints": [constrain_positive],
#"Y": self.positive_Y, "Y": self.positive_Y,
#"laplace": True "laplace": True
#} }
} }
for name, attributes in noise_models.iteritems(): for name, attributes in noise_models.iteritems():

View file

@ -56,8 +56,6 @@ def ICM(input_dim, num_outputs, kernel, W_rank=1,W=None,kappa=None,name='X'):
warnings.warn("kernel's input dimension overwritten to fit input_dim parameter.") warnings.warn("kernel's input dimension overwritten to fit input_dim parameter.")
K = kernel.prod(GPy.kern.Coregionalize(1, num_outputs, active_dims=[input_dim], rank=W_rank,W=W,kappa=kappa,name='B'),name=name) K = kernel.prod(GPy.kern.Coregionalize(1, num_outputs, active_dims=[input_dim], rank=W_rank,W=W,kappa=kappa,name='B'),name=name)
#K = kernel * GPy.kern.Coregionalize(1, num_outputs, active_dims=[input_dim], rank=W_rank,W=W,kappa=kappa,name='B')
#K = kernel ** GPy.kern.Coregionalize(input_dim, num_outputs,W_rank,W,kappa, name= 'B')
K['.*variance'] = 1. K['.*variance'] = 1.
K['.*variance'].fix() K['.*variance'].fix()
return K return K