mirror of
https://github.com/SheffieldML/GPy.git
synced 2026-06-02 14:45:15 +02:00
Adding gradients, shapes starting to make sense
This commit is contained in:
parent
1420aa532c
commit
267a8e427c
4 changed files with 60 additions and 29 deletions
|
|
@ -1,11 +1,12 @@
|
|||
import numpy as np
|
||||
import scipy as sp
|
||||
import GPy
|
||||
from scipy.linalg import cholesky, eig, inv, cho_solve, det
|
||||
from scipy.linalg import inv, cho_solve, det
|
||||
from numpy.linalg import cond
|
||||
from GPy.likelihoods.likelihood import likelihood
|
||||
from GPy.util.linalg import pdinv, mdot, jitchol, chol_inv, det_ln_diag, pddet
|
||||
from scipy.linalg.lapack import dtrtrs
|
||||
import pylab as plt
|
||||
|
||||
|
||||
class Laplace(likelihood):
|
||||
|
|
@ -62,7 +63,7 @@ class Laplace(likelihood):
|
|||
return self.likelihood_function._get_param_names()
|
||||
|
||||
def _set_params(self, p):
|
||||
return self.likelihood_function._set_params()
|
||||
return self.likelihood_function._set_params(p)
|
||||
|
||||
def both_gradients(self, dL_d_K_Sigma, dK_dthetaK):
|
||||
"""
|
||||
|
|
@ -77,8 +78,8 @@ class Laplace(likelihood):
|
|||
return (self._Kgradients(dL_d_K_Sigma, dK_dthetaK), self._gradients(dL_d_K_Sigma))
|
||||
|
||||
def _shared_gradients_components(self):
|
||||
dL_dytil = -np.dot((self.K+self.Sigma_tilde), self.Y)
|
||||
dytil_dfhat = np.dot(self.Sigma_tilde, self.Ki) + np.eye(self.N) # or self.Wi__Ki_W?
|
||||
dL_dytil = -np.dot(self.Y.T, (self.K+self.Sigma_tilde))
|
||||
dytil_dfhat = self.Wi__Ki_W # np.dot(self.Sigma_tilde, self.Ki) + np.eye(self.N) # or self.Wi__Ki_W?
|
||||
return dL_dytil, dytil_dfhat
|
||||
|
||||
def _Kgradients(self, dL_d_K_Sigma, dK_dthetaK):
|
||||
|
|
@ -91,12 +92,18 @@ class Laplace(likelihood):
|
|||
"""
|
||||
dL_dytil, dytil_dfhat = self._shared_gradients_components()
|
||||
|
||||
I_KW_i, _, _, _ = pdinv(np.eye(self.N) + np.dot(self.K, self.W))
|
||||
A = np.eye(self.N) + np.dot(self.K, self.W)
|
||||
plt.imshow(A)
|
||||
plt.show()
|
||||
I_KW_i, _, _, _ = pdinv(A)
|
||||
|
||||
#FIXME: Careful dK_dthetaK is not the derivative with respect to the marginal just prior K!
|
||||
dfhat_dthetaK = I_KW_i*dK_dthetaK*self.likelihood_function.link_grad(self.data, self.f_hat, self.extra_data)
|
||||
|
||||
dytil_dthetaK = dytil_dfhat*dfhat_dthetaK
|
||||
#Derivative for each f dimension, for each of K's hyper parameters
|
||||
dfhat_dthetaK = np.zeros((self.f_hat.shape[0], dK_dthetaK.shape[0]))
|
||||
for ind_j, thetaj in enumerate(dK_dthetaK):
|
||||
dfhat_dthetaK[:, ind_j] = mdot(I_KW_i, thetaj, self.likelihood_function.link_grad(self.data, self.f_hat, self.extra_data))
|
||||
|
||||
dytil_dthetaK = np.dot(dytil_dfhat, dfhat_dthetaK) # should be (D,thetaK)
|
||||
#FIXME: Careful dL_dK = dL_d_K_Sigma
|
||||
#FIXME: Careful the -D*0.5 in dL_d_K_sigma might need to be -0.5?
|
||||
dL_dSigma = dL_d_K_Sigma
|
||||
|
|
@ -105,8 +112,9 @@ class Laplace(likelihood):
|
|||
dSigmai_dthetaK = 0 #+ np.sum(d3phi_d3fhat*dfhat_dthetaK) #FIXME: CAREFUL OF THIS SUM! SHOULD SUM OVER FHAT NOT THETAS
|
||||
dSigma_dthetaK = -mdot(self.Sigma_tilde, dSigmai_dthetaK, self.Sigma_tilde)
|
||||
|
||||
dL_dthetaK_implicit = dL_dytil*dytil_dthetaK + dL_dSigma*dSigma_dthetaK
|
||||
return dL_dthetaK_implicit
|
||||
dL_dthetaK_implicit = np.sum(np.dot(dL_dytil, dytil_dthetaK), axis=0)# + np.dot(dL_dSigma, dSigma_dthetaK)
|
||||
#dL_dthetaK_implicit = np.dot(dL_dytil.T, dytil_dthetaK.T)
|
||||
return np.squeeze(dL_dthetaK_implicit)
|
||||
|
||||
def _gradients(self, partial):
|
||||
"""
|
||||
|
|
@ -132,16 +140,25 @@ class Laplace(likelihood):
|
|||
partial = dL_dK
|
||||
"""
|
||||
dL_dytil, dytil_dfhat = self._shared_gradients_components()
|
||||
dfhat_dthetaL = self.likelihood_function.df_dtheta()
|
||||
dfhat_dthetaL, dSigmai_dthetaL = self.likelihood_function._gradients(self.data, self.f_hat, self.extra_data) #FIXME: Shouldn't this have a implicit component aswell?
|
||||
|
||||
dSigmai_dthetaL = self.likelihood_function._gradients(self.data, self.f_hat, self.extra_data) #FIXME: Shouldn't this have a implicit component aswell?
|
||||
dSigma_dthetaL = -mdot(self.Sigma_tilde, dSigmai_dthetaL, self.Sigma_tilde)
|
||||
#dSigmai_dthetaL = self.likelihood_function._gradients(self.data, self.f_hat, self.extra_data) #FIXME: Shouldn't this have a implicit component aswell?
|
||||
#Derivative for each f dimension, for each of K's hyper parameters
|
||||
dSigma_dthetaL = np.empty((self.N, len(self.likelihood_function._get_param_names())))
|
||||
for ind_l, dSigmai_dtheta_l in enumerate(dSigmai_dthetaL.T):
|
||||
dSigma_dthetaL[:, ind_l] = -mdot(self.Sigma_tilde,
|
||||
dSigmai_dtheta_l, # Careful, shouldn't this be (N, 1)?
|
||||
self.Sigma_tilde
|
||||
)
|
||||
|
||||
#TODO: This is Wi*A*Wi, can be more numerically stable with a trick
|
||||
#dSigma_dthetaL = -mdot(self.Sigma_tilde, dSigmai_dthetaL, self.Sigma_tilde)
|
||||
dL_dSigma = partial # partial is dL_dK but K here is K+Sigma_tilde.... which is fine in this case
|
||||
|
||||
dytil_dthetaL = dytil_dfhat*dfhat_dthetaL
|
||||
dL_dthetaL = 0 + dL_dytil*dytil_dthetaL + dL_dSigma*dSigma_dthetaL
|
||||
return dL_dthetaL
|
||||
#return np.zeros(0) # TODO: Laplace likelihood might want to take some parameters...
|
||||
#dytil_dthetaL = dytil_dfhat*dfhat_dthetaL
|
||||
dytil_dthetaL = np.dot(dytil_dfhat, dfhat_dthetaL)
|
||||
dL_dthetaL = 0 + np.dot(dL_dytil, dytil_dthetaL)# + np.dot(dL_dSigma, dSigma_dthetaL)
|
||||
return np.squeeze(dL_dthetaL) #should be array of length *params-being optimized*, for student t just optimising 1 parameter, this is (1,)
|
||||
|
||||
def _compute_GP_variables(self):
|
||||
"""
|
||||
|
|
@ -335,7 +352,7 @@ class Laplace(likelihood):
|
|||
rs = 0
|
||||
i = 0
|
||||
while difference > epsilon and i < MAX_ITER and rs < MAX_RESTART:
|
||||
f_old = f.copy()
|
||||
#f_old = f.copy()
|
||||
W = -np.diag(self.likelihood_function.link_hess(self.data, f, extra_data=self.extra_data))
|
||||
if not self.likelihood_function.log_concave:
|
||||
W[W < 0] = 1e-6 # FIXME-HACK: This is a hack since GPy can't handle negative variances which can occur
|
||||
|
|
|
|||
|
|
@ -159,10 +159,10 @@ class student_t(likelihood_function):
|
|||
d2ln p(yi|fi)_d2fifj
|
||||
"""
|
||||
def __init__(self, deg_free, sigma=2):
|
||||
super(student_t, self).__init__()
|
||||
self.v = deg_free
|
||||
self.sigma = sigma
|
||||
self.log_concave = False
|
||||
#super(student_t, self).__init__()
|
||||
|
||||
def _get_params(self):
|
||||
return np.asarray(self.sigma)
|
||||
|
|
@ -258,9 +258,9 @@ class student_t(likelihood_function):
|
|||
)
|
||||
return d3link_d3f
|
||||
|
||||
def link_hess_grad_sigma(self, y, f, extra_data=None):
|
||||
def link_hess_grad_std(self, y, f, extra_data=None):
|
||||
"""
|
||||
Gradient of the hessian w.r.t sigma parameter
|
||||
Gradient of the hessian w.r.t sigma parameter (standard deviation)
|
||||
|
||||
$$\frac{2\sigma v(v+1)(\sigma^{2}v - 3(f-y)^2)}{((f-y)^{2} + \sigma^{2}v)^{3}}
|
||||
"""
|
||||
|
|
@ -273,8 +273,24 @@ class student_t(likelihood_function):
|
|||
)
|
||||
return hess_grad_sigma
|
||||
|
||||
def link_grad_std(self, y, f, extra_data=None):
|
||||
"""
|
||||
Gradient of the likelihood w.r.t sigma parameter (standard deviation)
|
||||
|
||||
$$\frac{-2\sigma(v+1)(y-f)}{(v\sigma^{2} + (y-f)^{2})^{2}}$$
|
||||
"""
|
||||
y = np.squeeze(y)
|
||||
f = np.squeeze(f)
|
||||
assert y.shape == f.shape
|
||||
e = y - f
|
||||
grad_sigma = ( (-2*self.sigma*self.v*(self.v + 1)*e)
|
||||
/ ((self.v*(self.sigma**2) + e**2)**2)
|
||||
)
|
||||
return grad_sigma
|
||||
|
||||
def _gradients(self, y, f, extra_data=None):
|
||||
return [self.link_hess_grad_sigma] # list as we might learn many parameters
|
||||
return [self.link_grad_std(y, f, extra_data=extra_data)[:, None],
|
||||
self.link_hess_grad_std(y, f, extra_data=extra_data)[:, None]] # list as we might learn many parameters
|
||||
|
||||
def predictive_values(self, mu, var):
|
||||
"""
|
||||
|
|
@ -372,9 +388,7 @@ class weibull_survival(likelihood_function):
|
|||
def __init__(self, shape, scale):
|
||||
self.shape = shape
|
||||
self.scale = scale
|
||||
|
||||
#FIXME: This should be in the superclass
|
||||
self.log_concave = True
|
||||
self.log_concave = True # Or false?
|
||||
|
||||
def link_function(self, y, f, extra_data=None):
|
||||
"""
|
||||
|
|
|
|||
|
|
@ -128,17 +128,17 @@ class GP(model):
|
|||
|
||||
For the likelihood parameters, pass in alpha = K^-1 y
|
||||
"""
|
||||
dL_dthetaK = self.kern.dK_dtheta(dL_dK=self.dL_dK, X=self.X, slices1=self.Xslices, slices2=self.Xslices)
|
||||
if isinstance(self.likelihood, Laplace):
|
||||
dL_dthetaK_explicit = self.kern.dK_dtheta(dL_dK=self.dL_dK, X=self.X, slices1=self.Xslices, slices2=self.Xslices)
|
||||
dL_dthetaK_explicit = dL_dthetaK
|
||||
#Need to pass in a matrix of ones to get access to raw dK_dthetaK values without being chained
|
||||
fake_dL_dKs = np.ones(self.dL_dK.shape)
|
||||
dK_dthetaK = self.kern.dK_dtheta(dL_dK=fake_dL_dKs, X=self.X, slices1=self.Xslices, slices2=self.Xslices)
|
||||
|
||||
dL_dthetaK_implicit = self.likelihood._Kgradients(self.dL_dK, dK_dthetaK)
|
||||
dL_dthetaK = dL_dthetaK_explicit + dL_dthetaK_implicit
|
||||
dL_dthetaL = self.likelihood._gradients(partial=np.diag(self.dL_dK))
|
||||
dL_dthetaL = self.likelihood._gradients(partial=self.dL_dK)
|
||||
else:
|
||||
dL_dthetaK = self.kern.dK_dtheta(dL_dK=self.dL_dK, X=self.X, slices1=self.Xslices, slices2=self.Xslices)
|
||||
dL_dthetaL = self.likelihood._gradients(partial=np.diag(self.dL_dK))
|
||||
return np.hstack((dL_dthetaK, dL_dthetaL))
|
||||
|
||||
|
|
|
|||
|
|
@ -25,7 +25,7 @@ def pddet(A):
|
|||
"""
|
||||
Determinant of a positive definite matrix
|
||||
"""
|
||||
L = cholesky(A)
|
||||
L = jitchol(A)
|
||||
logdetA = 2*sum(np.log(np.diag(L)))
|
||||
return logdetA
|
||||
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue