diff --git a/GPy/core/gp.py b/GPy/core/gp.py index ca2583e0..f4b02128 100644 --- a/GPy/core/gp.py +++ b/GPy/core/gp.py @@ -40,7 +40,7 @@ class GP(Model): assert Y.ndim == 2 self.Y = ObsAr(Y) - assert Y.shape[0] == self.num_data +# assert Y.shape[0] == self.num_data _, self.output_dim = self.Y.shape #TODO: check the type of this is okay? diff --git a/GPy/inference/latent_function_inference/posterior.py b/GPy/inference/latent_function_inference/posterior.py index 70b0e0ea..66c68261 100644 --- a/GPy/inference/latent_function_inference/posterior.py +++ b/GPy/inference/latent_function_inference/posterior.py @@ -153,9 +153,14 @@ class Posterior(object): $$ """ if self._woodbury_inv is None: - self._woodbury_inv, _ = dpotri(self.woodbury_chol, lower=1) - #self._woodbury_inv, _ = dpotrs(self.woodbury_chol, np.eye(self.woodbury_chol.shape[0]), lower=1) - symmetrify(self._woodbury_inv) + if self._woodbury_chol is not None: + self._woodbury_inv, _ = dpotri(self._woodbury_chol, lower=1) + #self._woodbury_inv, _ = dpotrs(self.woodbury_chol, np.eye(self.woodbury_chol.shape[0]), lower=1) + symmetrify(self._woodbury_inv) + elif self._covariance is not None: + B = self._K - self._covariance + tmp, _ = dpotrs(self.K_chol, B) + self._woodbury_inv, _ = dpotrs(self.K_chol, tmp.T) return self._woodbury_inv @property diff --git a/GPy/kern/_src/linear.py b/GPy/kern/_src/linear.py index 649c1222..9a1388c5 100644 --- a/GPy/kern/_src/linear.py +++ b/GPy/kern/_src/linear.py @@ -11,8 +11,8 @@ from ...core.parameterization import Param from ...core.parameterization.transformations import Logexp from ...util.caching import Cache_this from ...core.parameterization import variational -from psi_comp import linear_psi_comp from ...util.config import * +from .psi_comp import PSICOMP_Linear class Linear(Kern): """ @@ -53,10 +53,8 @@ class Linear(Kern): self.variances = Param('variances', variances, Logexp()) self.add_parameter(self.variances) + self.psicomp = PSICOMP_Linear() - def set_for_SpikeAndSlab(self): - self.psicomp = linear_psi_comp.PSICOMP_SSLinear() - @Cache_this(limit=2) def K(self, X, X2=None): if self.ARD: diff --git a/GPy/kern/_src/psi_comp/__init__.py b/GPy/kern/_src/psi_comp/__init__.py index 77668e7f..38094394 100644 --- a/GPy/kern/_src/psi_comp/__init__.py +++ b/GPy/kern/_src/psi_comp/__init__.py @@ -6,20 +6,7 @@ from GPy.util.caching import Cache_this from ....core.parameterization import variational import rbf_psi_comp import ssrbf_psi_comp - -class PSICOMP(Pickleable): - - def psicomputations(self, variance, Z, variational_posterior): - """ - Compute psi-statistics - """ - pass - - def psiDerivativecomputations(self, dL_dpsi0, dL_dpsi1, dL_dpsi2, variance, Z, variational_posterior): - """ - Compute the derivatives of parameters by combing dL_dpsi and dpsi_dparam - """ - pass +import sslinear_psi_comp class PSICOMP_RBF(Pickleable): @@ -38,5 +25,25 @@ class PSICOMP_RBF(Pickleable): return rbf_psi_comp.psiDerivativecomputations(dL_dpsi0, dL_dpsi1, dL_dpsi2, variance, lengthscale, Z, variational_posterior) elif isinstance(variational_posterior, variational.SpikeAndSlabPosterior): return ssrbf_psi_comp.psiDerivativecomputations(dL_dpsi0, dL_dpsi1, dL_dpsi2, variance, lengthscale, Z, variational_posterior) + else: + raise ValueError, "unknown distriubtion received for psi-statistics" + +class PSICOMP_Linear(Pickleable): + + @Cache_this(limit=1, ignore_args=(0,)) + def psicomputations(self, variance, Z, variational_posterior): + if isinstance(variational_posterior, variational.NormalPosterior): + raise NotImplementedError + elif isinstance(variational_posterior, variational.SpikeAndSlabPosterior): + return sslinear_psi_comp.psicomputations(variance, Z, variational_posterior) + else: + raise ValueError, "unknown distriubtion received for psi-statistics" + + @Cache_this(limit=1, ignore_args=(0,1,2,3)) + def psiDerivativecomputations(self, dL_dpsi0, dL_dpsi1, dL_dpsi2, variance, Z, variational_posterior): + if isinstance(variational_posterior, variational.NormalPosterior): + raise NotImplementedError + elif isinstance(variational_posterior, variational.SpikeAndSlabPosterior): + return sslinear_psi_comp.psiDerivativecomputations(dL_dpsi0, dL_dpsi1, dL_dpsi2, variance, Z, variational_posterior) else: raise ValueError, "unknown distriubtion received for psi-statistics" \ No newline at end of file diff --git a/GPy/kern/_src/psi_comp/linear_psi_comp.py b/GPy/kern/_src/psi_comp/linear_psi_comp.py deleted file mode 100644 index 3404c987..00000000 --- a/GPy/kern/_src/psi_comp/linear_psi_comp.py +++ /dev/null @@ -1,93 +0,0 @@ -# Copyright (c) 2012, GPy authors (see AUTHORS.txt). -# Licensed under the BSD 3-clause license (see LICENSE.txt) - -""" -The package for the Psi statistics computation of the linear kernel for SSGPLVM -""" - -import numpy as np -from . import PSICOMP -from GPy.util.caching import Cache_this - -class PSICOMP_SSLinear(PSICOMP): - @Cache_this(limit=1, ignore_args=(0,)) - def psicomputations(self, variance, Z, variational_posterior): - """ - Compute psi-statistics for ss-linear kernel - """ - # here are the "statistics" for psi0, psi1 and psi2 - # Produced intermediate results: - # psi0 N - # psi1 NxM - # psi2 MxM - mu = variational_posterior.mean - S = variational_posterior.variance - gamma = variational_posterior.binary_prob - - psi0 = np.einsum('q,nq,nq->n',variance,gamma,np.square(mu)+S) - psi1 = np.einsum('nq,q,mq,nq->nm',gamma,variance,Z,mu) - mu2 = np.square(mu) - variances2 = np.square(variance) - tmp = np.einsum('nq,q,mq,nq->nm',gamma,variance,Z,mu) - psi2 = np.einsum('nq,q,mq,oq,nq->mo',gamma,variances2,Z,Z,mu2+S)+\ - np.einsum('nm,no->mo',tmp,tmp) - np.einsum('nq,q,mq,oq,nq->mo',np.square(gamma),variances2,Z,Z,mu2) - - return psi0, psi1, psi2 - - @Cache_this(limit=1, ignore_args=(0,1,2,3)) - def psiDerivativecomputations(self, dL_dpsi0, dL_dpsi1, dL_dpsi2, variance, Z, variational_posterior): - mu = variational_posterior.mean - S = variational_posterior.variance - gamma = variational_posterior.binary_prob - - dL_dvar, dL_dgamma, dL_dmu, dL_dS, dL_dZ = self._psi2computations(dL_dpsi2, variance, Z, mu, S, gamma) - - # Compute for psi0 and psi1 - mu2S = np.square(mu)+S - dL_dvar += np.einsum('n,nq,nq->q',dL_dpsi0,gamma,mu2S) + np.einsum('nm,nq,mq,nq->q',dL_dpsi1,gamma,Z,mu) - dL_dgamma += np.einsum('n,q,nq->nq',dL_dpsi0,variance,mu2S) + np.einsum('nm,q,mq,nq->nq',dL_dpsi1,variance,Z,mu) - dL_dmu += np.einsum('n,nq,q,nq->nq',dL_dpsi0,gamma,2.*variance,mu) + np.einsum('nm,nq,q,mq->nq',dL_dpsi1,gamma,variance,Z) - dL_dS += np.einsum('n,nq,q->nq',dL_dpsi0,gamma,variance) - dL_dZ += np.einsum('nm,nq,q,nq->mq',dL_dpsi1,gamma, variance,mu) - - return dL_dvar, dL_dZ, dL_dmu, dL_dS, dL_dgamma - - def _psi2computations(self, dL_dpsi2, variance, Z, mu, S, gamma): - """ - Z - MxQ - mu - NxQ - S - NxQ - gamma - NxQ - """ - # here are the "statistics" for psi1 and psi2 - # Produced intermediate results: - # _psi2_dvariance Q - # _psi2_dZ MxQ - # _psi2_dgamma NxQ - # _psi2_dmu NxQ - # _psi2_dS NxQ - - mu2 = np.square(mu) - gamma2 = np.square(gamma) - variance2 = np.square(variance) - mu2S = mu2+S # NxQ - common_sum = np.einsum('nq,q,mq,nq->nm',gamma,variance,Z,mu) # NxM - - dL_dvar = np.einsum('mo,nq,q,mq,oq->q',dL_dpsi2,2.*(gamma*mu2S-gamma2*mu2),variance,Z,Z)+\ - np.einsum('mo,nq,mq,nq,no->q',dL_dpsi2,gamma,Z,mu,common_sum)+\ - np.einsum('mo,nq,oq,nq,nm->q',dL_dpsi2,gamma,Z,mu,common_sum) - - dL_dgamma = np.einsum('mo,q,mq,oq,nq->nq',dL_dpsi2,variance2,Z,Z,(mu2S-2.*gamma*mu2))+\ - np.einsum('mo,q,mq,nq,no->nq',dL_dpsi2,variance,Z,mu,common_sum)+\ - np.einsum('mo,q,oq,nq,nm->nq',dL_dpsi2,variance,Z,mu,common_sum) - - dL_dmu = np.einsum('mo,q,mq,oq,nq,nq->nq',dL_dpsi2,variance2,Z,Z,mu,2.*(gamma-gamma2))+\ - np.einsum('mo,nq,q,mq,no->nq',dL_dpsi2,gamma,variance,Z,common_sum)+\ - np.einsum('mo,nq,q,oq,nm->nq',dL_dpsi2,gamma,variance,Z,common_sum) - - dL_dS = np.einsum('mo,nq,q,mq,oq->nq',dL_dpsi2,gamma,variance2,Z,Z) - - dL_dZ = 2.*(np.einsum('om,nq,q,mq,nq->oq',dL_dpsi2,gamma,variance2,Z,mu2S)+np.einsum('om,nq,q,nq,nm->oq',dL_dpsi2,gamma,variance,mu,common_sum) - -np.einsum('om,nq,q,mq,nq->oq',dL_dpsi2,gamma2,variance2,Z,mu2)) - - return dL_dvar, dL_dgamma, dL_dmu, dL_dS, dL_dZ diff --git a/GPy/util/linalg.py b/GPy/util/linalg.py index 661a2b47..8f79ae91 100644 --- a/GPy/util/linalg.py +++ b/GPy/util/linalg.py @@ -16,8 +16,7 @@ import warnings import os from config import * -if np.all(np.float64((scipy.__version__).split('.')[:2]) >= np.array([0, 12])): - #import scipy.linalg.lapack.clapack as lapack +if float('.'.join((scipy.__version__).split('.')[:2])) >= 0.12: from scipy.linalg import lapack else: from scipy.linalg.lapack import flapack as lapack