mirror of
https://github.com/SheffieldML/GPy.git
synced 2026-06-11 15:15:15 +02:00
Took out all the asserts and using pure broadcasting method of diagonal now
This commit is contained in:
Alan Saul
parent
e842f6e687
commit
6c29750795
3 changed files with 20 additions and 57 deletions
|
|
@ -39,8 +39,8 @@ def debug_student_t_noise_approx():
|
|||
plot = False
|
||||
real_var = 0.1
|
||||
#Start a function, any function
|
||||
X = np.linspace(0.0, 10.0, 100)[:, None]
|
||||
#X = np.array([0.5])[:, None]
|
||||
#X = np.linspace(0.0, 10.0, 100)[:, None]
|
||||
X = np.array([0.5])[:, None]
|
||||
Y = np.sin(X) + np.random.randn(*X.shape)*real_var
|
||||
|
||||
X_full = np.linspace(0.0, 10.0, 500)[:, None]
|
||||
|
|
|
|||
|
|
@ -69,9 +69,7 @@ class Laplace(likelihood):
|
|||
#FIXME: Careful of side effects! And make sure W and K are up to date!
|
||||
d3lik_d3fhat = self.likelihood_function.d3lik_d3f(self.data, self.f_hat)
|
||||
dL_dfhat = -0.5*(np.diag(self.Ki_W_i)[:, None]*d3lik_d3fhat)
|
||||
Wi_K_i = mdot(np.diagflat(self.W_12), self.Bi, np.diagflat(self.W_12)) #same as rasms R
|
||||
Wi_K_inew = self.W_12*self.Bi*self.W_12.T #same as rasms R
|
||||
assert np.all(Wi_K_i == Wi_K_inew)
|
||||
Wi_K_i = self.W_12*self.Bi*self.W_12.T #same as rasms R
|
||||
|
||||
I_KW_i = np.eye(self.N) - np.dot(self.K, Wi_K_i)
|
||||
return dL_dfhat, I_KW_i, Wi_K_i
|
||||
|
|
@ -150,9 +148,7 @@ class Laplace(likelihood):
|
|||
#((L.T*w)_i + I)f_hat = y_tilde
|
||||
L = jitchol(self.K)
|
||||
Li = chol_inv(L)
|
||||
Lt_W = np.dot(L.T, np.diagflat(self.W)) #FIXME: Can make Faster
|
||||
Lt_Wnew = L.T*self.W.T
|
||||
assert np.all(Lt_Wnew == Lt_W)
|
||||
Lt_W = L.T*self.W.T
|
||||
|
||||
##Check it isn't singular!
|
||||
if cond(Lt_W) > epsilon:
|
||||
|
|
@ -164,25 +160,15 @@ class Laplace(likelihood):
|
|||
|
||||
#f.T(Ki + W)f
|
||||
f_Ki_W_f = (np.dot(self.f_hat.T, cho_solve((L, True), self.f_hat))
|
||||
+ mdot(self.f_hat.T, np.diagflat(self.W), self.f_hat)
|
||||
)
|
||||
f_Ki_W_fnew = (np.dot(self.f_hat.T, cho_solve((L, True), self.f_hat))
|
||||
+ mdot(self.f_hat.T, self.W*self.f_hat)
|
||||
)
|
||||
assert np.all(f_Ki_W_f == f_Ki_W_fnew)
|
||||
|
||||
y_W_f = mdot((Y_tilde.T, np.diagflat(self.W)), self.f_hat)
|
||||
y_W_fnew = mdot(Y_tilde.T*self.W.T, self.f_hat)
|
||||
assert np.all(y_W_f == y_W_fnew)
|
||||
y_W_f = mdot(Y_tilde.T*self.W.T, self.f_hat)
|
||||
|
||||
|
||||
y_W_y = mdot((Y_tilde.T, np.diagflat(self.W)), Y_tilde)
|
||||
y_W_ynew = mdot(Y_tilde.T, self.W*Y_tilde)
|
||||
assert np.all(y_W_y == y_W_ynew)
|
||||
y_W_y = mdot(Y_tilde.T, self.W*Y_tilde)
|
||||
|
||||
ln_W_det = det_ln_diag(np.diagflat(self.W))
|
||||
ln_W_detnew = np.log(self.W).sum()
|
||||
assert np.all(ln_W_det == ln_W_detnew)
|
||||
ln_W_det = np.log(self.W).sum()
|
||||
|
||||
#FIXME: Revisit this
|
||||
Z_tilde = (- self.NORMAL_CONST
|
||||
|
|
@ -203,15 +189,13 @@ class Laplace(likelihood):
|
|||
#+ y_W_f
|
||||
#+ self.ln_z_hat
|
||||
#)
|
||||
self.Z_tilde = 0
|
||||
#self.Z_tilde = 0
|
||||
|
||||
##Check it isn't singular!
|
||||
if cond(self.W) > epsilon:
|
||||
print "WARNING: Transformed covariance matrix is singular,\nnumerical stability may be a problem"
|
||||
|
||||
self.Sigma_tilde = inv(np.diagflat(self.W)) # Damn
|
||||
Sigma_tildenew = np.diagflat(1.0/self.W)
|
||||
assert np.all(self.Sigma_tilde == Sigma_tildenew)
|
||||
self.Sigma_tilde = np.diagflat(1.0/self.W)
|
||||
|
||||
#Convert to float as its (1, 1) and Z must be a scalar
|
||||
self.Z = np.float64(Z_tilde)
|
||||
|
|
@ -251,23 +235,15 @@ class Laplace(likelihood):
|
|||
self.B, self.B_chol, self.W_12 = self._compute_B_statistics(self.K, self.W)
|
||||
self.Bi, _, _, B_det = pdinv(self.B)
|
||||
|
||||
self.Ki_W_i = self.K - mdot(self.K, (np.diagflat(self.W_12), self.Bi, np.diagflat(self.W_12)), self.K) # Funky, order matters on stability!
|
||||
Ki_W_inew = self.K - mdot(self.K, self.W_12*self.Bi*self.W_12.T, self.K)
|
||||
assert np.all(self.Ki_W_i == Ki_W_inew)
|
||||
self.Ki_W_i = self.K - mdot(self.K, self.W_12*self.Bi*self.W_12.T, self.K)
|
||||
|
||||
self.ln_Ki_W_i_det = np.linalg.det(self.Ki_W_i)
|
||||
|
||||
b = np.dot(np.diagflat(self.W), self.f_hat) + self.likelihood_function.dlik_df(self.data, self.f_hat, extra_data=self.extra_data)
|
||||
bnew = self.W*self.f_hat + self.likelihood_function.dlik_df(self.data, self.f_hat, extra_data=self.extra_data)
|
||||
assert np.all(b == bnew)
|
||||
b = self.W*self.f_hat + self.likelihood_function.dlik_df(self.data, self.f_hat, extra_data=self.extra_data)
|
||||
|
||||
solve_chol = cho_solve((self.B_chol, True), mdot((np.diagflat(self.W_12), self.K), b))
|
||||
solve_cholnew = cho_solve((self.B_chol, True), np.dot(self.W_12*self.K, b))
|
||||
assert np.all(solve_chol == solve_cholnew)
|
||||
solve_chol = cho_solve((self.B_chol, True), np.dot(self.W_12*self.K, b))
|
||||
|
||||
a = b - mdot(np.diagflat(self.W_12), solve_chol)
|
||||
anew = b - self.W_12*solve_chol
|
||||
assert np.all(a == anew)
|
||||
a = b - self.W_12*solve_chol
|
||||
|
||||
self.Ki_f = a
|
||||
self.f_Ki_f = np.dot(self.f_hat.T, self.Ki_f)
|
||||
|
|
@ -291,10 +267,6 @@ class Laplace(likelihood):
|
|||
"""
|
||||
#W is diagonal so its sqrt is just the sqrt of the diagonal elements
|
||||
W_12 = np.sqrt(W)
|
||||
# FIXME Take this out when you've done multiinput, Weirdly this is
|
||||
# better when its W_12.T*K*W_12 which shouldnt make a difference
|
||||
# because K is symmetrical
|
||||
assert np.allclose(W_12*K*W_12.T, np.dot(np.diagflat(W_12), np.dot(K, np.diagflat(W_12))))
|
||||
B = np.eye(self.N) + W_12*K*W_12.T
|
||||
L = jitchol(B)
|
||||
return (B, L, W_12)
|
||||
|
|
@ -360,9 +332,7 @@ class Laplace(likelihood):
|
|||
# This is a property only held by non-log-concave likelihoods
|
||||
B, L, W_12 = self._compute_B_statistics(K, W)
|
||||
|
||||
W_f = np.dot(np.diagflat(W), f)
|
||||
W_fnew = W*f
|
||||
assert np.all(W_f == W_fnew)
|
||||
W_f = W*f
|
||||
grad = self.likelihood_function.dlik_df(self.data, f, extra_data=self.extra_data)
|
||||
#Find K_i_f
|
||||
b = W_f + grad
|
||||
|
|
@ -370,21 +340,13 @@ class Laplace(likelihood):
|
|||
#a should be equal to Ki*f now so should be able to use it
|
||||
c = np.dot(K, W_f) + f*(1-step_size) + step_size*np.dot(K, grad)
|
||||
|
||||
solve_L = cho_solve((L, True), np.dot(np.diagflat(W_12), c))
|
||||
solve_Lnew = cho_solve((L, True), W_12*c)
|
||||
assert np.all(solve_L == solve_Lnew)
|
||||
solve_L = cho_solve((L, True), W_12*c)
|
||||
|
||||
f = c - np.dot(K, np.dot(np.diagflat(W_12), solve_L))
|
||||
fnew = c - np.dot(K, W_12*solve_L)
|
||||
assert np.all(f == fnew)
|
||||
f = c - np.dot(K, W_12*solve_L)
|
||||
|
||||
solve_L = cho_solve((L, True), np.dot(np.diagflat(W_12), np.dot(K, b)))
|
||||
solve_Lnew = cho_solve((L, True), W_12*np.dot(K, b))
|
||||
assert np.all(solve_L == solve_Lnew)
|
||||
solve_L = cho_solve((L, True), W_12*np.dot(K, b))
|
||||
|
||||
a = b - np.dot(np.diagflat(W_12), solve_L)
|
||||
anew = b - W_12*solve_L
|
||||
assert np.all(a == anew)
|
||||
a = b - W_12*solve_L
|
||||
|
||||
tmp_old_obj = old_obj
|
||||
old_obj = new_obj
|
||||
|
|
|
|||
|
|
@ -152,8 +152,9 @@ class GP(model):
|
|||
#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) #FIXME: Check this is right...
|
||||
#fake_dL_dKs = np.eye(self.dL_dK.shape[0]) #FIXME: Check this is right...
|
||||
|
||||
#BUG: THIS SHOULD NOT BE (1,num_k_params) matrix it should be (N,N,num_k_params)
|
||||
dK_dthetaK = self.kern.dK_dtheta(dL_dK=fake_dL_dKs, X=self.X)
|
||||
#THIS SHOULD NOT BE (1,num_k_params) matrix it should be (N,N,num_k_params)
|
||||
|
||||
dL_dthetaK = self.likelihood._Kgradients(dK_dthetaK=dK_dthetaK)
|
||||
dL_dthetaL = self.likelihood._gradients(partial=np.diag(self.dL_dK))
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue