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

GPy/core/transformations.py

@@ -39,8 +39,8 @@ class logexp(transformation):
         return '(+ve)'
 
 class logexp_clipped(transformation):
-    max_bound = 1e300
-    min_bound = 1e-10
+    max_bound = 1e250
+    min_bound = 1e-9
     log_max_bound = np.log(max_bound)
     log_min_bound = np.log(min_bound)
     def __init__(self, lower=1e-6):
@@ -49,11 +49,13 @@ class logexp_clipped(transformation):
     def f(self, x):
         exp = np.exp(np.clip(x, self.log_min_bound, self.log_max_bound))
         f = np.log(1. + exp)
+        if np.isnan(f).any():
+            import ipdb;ipdb.set_trace()
         return f
     def finv(self, f):
         return np.log(np.exp(np.clip(f, self.min_bound, self.max_bound)) - 1.)
     def gradfactor(self, f):
-        ef = np.exp(f)
+        ef = np.exp(f) # np.clip(f, self.min_bound, self.max_bound))
         gf = (ef - 1.) / ef
         return np.where(f < self.lower, 0, gf)
     def initialize(self, f):

GPy/examples/dimensionality_reduction.py

@@ -273,8 +273,8 @@ def bgplvm_simulation(optimize='scg',
         pylab.figure(); pylab.axis(); m.kern.plot_ARD()
     return m
 
-def mrd_simulation(plot_sim=False):
-    D1, D2, D3, N, M, Q = 150, 250, 300, 700, 3, 7
+def mrd_simulation(optimize=True, plot_sim=False):
+    D1, D2, D3, N, M, Q = 150, 250, 30, 300, 3, 7
     slist, Slist, Ylist = _simulate_sincos(D1, D2, D3, N, M, Q, plot_sim)
 
     from GPy.models import mrd
@@ -292,6 +292,13 @@ def mrd_simulation(plot_sim=False):
     m.constrain('variance|noise', logexp_clipped())
     m.ensure_default_constraints()
 
+    # DEBUG
+    np.seterr("raise")
+
+    if optimize:
+        print "Optimizing Model:"
+        m.optimize('scg', messages=1, max_iters=3e3)
+
     return m
 
 def brendan_faces():

GPy/inference/SCG.py

@@ -85,8 +85,6 @@ def SCG(f, gradf, x, optargs=(), maxiters=500, max_f_eval=500, display=True, xto
         # Increase effective curvature and evaluate step size alpha.
         delta = theta + beta * kappa
         if delta <= 0:
-            if display:
-                print ""
             delta = beta * kappa
             beta = beta - theta / kappa
 

GPy/models/Bayesian_GPLVM.py

@@ -171,9 +171,6 @@ class Bayesian_GPLVM(sparse_GP, GPLVM):
         self.dbound_dZtheta = sparse_GP._log_likelihood_gradients(self)
         return np.hstack((self.dbound_dmuS.flatten(), self.dbound_dZtheta))
 
-    def _log_likelihood_normal_gradients(self):
-        Si, _, _, _ = pdinv(self.X_variance)
-
     def plot_latent(self, which_indices=None, *args, **kwargs):
 
         if which_indices is None:

GPy/models/sparse_GP.py

@@ -3,7 +3,7 @@
 
 import numpy as np
 import pylab as pb
-from ..util.linalg import mdot, jitchol, tdot, symmetrify, backsub_both_sides
+from ..util.linalg import mdot, jitchol, tdot, symmetrify, backsub_both_sides,chol_inv
 from ..util.plot import gpplot
 from .. import kern
 from GP import GP
@@ -111,7 +111,7 @@ class sparse_GP(GP):
 
         if self.likelihood.is_heteroscedastic:
             if self.has_uncertain_inputs:
-                self.dL_dpsi2 = self.likelihood.precision[:, None, None] * dL_dpsi2_beta[None, :, :]
+                self.dL_dpsi2 = self.likelihood.precision.flatten()[:, None, None] * dL_dpsi2_beta[None, :, :]
             else:
                 self.dL_dpsi1 += 2.*np.dot(dL_dpsi2_beta, self.psi1 * self.likelihood.precision.reshape(1, self.N))
                 self.dL_dpsi2 = None
@@ -173,7 +173,13 @@ class sparse_GP(GP):
         this function does nothing
         """
         if self.has_uncertain_inputs:
-            raise NotImplementedError, "EP approximation not implemented for uncertain inputs"
+
+            Lmi = chol_inv(self.Lm)
+            Kmmi = tdot(Lmi.T)
+            diag_tr_psi2Kmmi = np.array([np.trace(psi2_Kmmi) for psi2_Kmmi in np.dot(self.psi2,Kmmi)])
+
+            self.likelihood.fit_FITC(self.Kmm,self.psi1,diag_tr_psi2Kmmi) #This uses the fit_FITC code, but does not perfomr a FITC-EP.#TODO solve potential confusion
+            #raise NotImplementedError, "EP approximation not implemented for uncertain inputs"
         else:
             self.likelihood.fit_DTC(self.Kmm, self.psi1)
             # self.likelihood.fit_FITC(self.Kmm,self.psi1,self.psi0)

GPy/util/datasets.py

@@ -5,6 +5,8 @@ import GPy
 import scipy.sparse
 import scipy.io
 import cPickle as pickle
+import urllib2 as url
+
 data_path = os.path.join(os.path.dirname(__file__), 'datasets')
 default_seed = 10000
 
@@ -15,6 +17,18 @@ def sample_class(f):
     c = np.where(c, 1, -1)
     return c
 
+def fetch_dataset(resource, file_name, messages = True):
+    if messages:
+        print "Downloading resource: " , resource, " ... "
+    response = url.urlopen(resource)
+    # TODO: Some error checking...
+    html = response.read()
+    response.close()
+    with open(file_name, "w") as text_file:
+        text_file.write("%s"%html)
+    if messages:
+        print "Done!"
+
 def della_gatta_TRP63_gene_expression(gene_number=None):
     mat_data = scipy.io.loadmat(os.path.join(data_path, 'DellaGattadata.mat'))
     X = np.double(mat_data['timepoints'])

GPy/util/mocap_fetch.py

@@ -0,0 +1,13 @@
+import GPy
+import urllib2
+
+# TODO...
+class mocap_fetch(base_url = 'http://mocap.cs.cmu.edu:8080/subjects/', skel_store_dir = './', motion_store_dir = './'):
+    def __init__(self):
+        self.base_url = base_url
+        self.store_dir = store_dir
+        self.motion_dict = []
+    
+    def fetch_motions(self, motion_dict = None):
+        response = urllib2.urlopen(...)
+        html = response.read()