Samples in plot_f fixed for sparse_models

GPy/core/gp_base.py

@@ -89,7 +89,9 @@ class GPBase(Model):
             fig = pb.figure(num=fignum)
             ax = fig.add_subplot(111)
 
-        if self.X.shape[1] == 1 and not hasattr(self,'multioutput'):
+        if not hasattr(self,'multioutput'):
+
+            if self.X.shape[1] == 1:
                 Xnew, xmin, xmax = x_frame1D(self.X, plot_limits=plot_limits)
                 if samples == 0:
                     m, v = self._raw_predict(Xnew, which_parts=which_parts)
@@ -97,17 +99,23 @@ class GPBase(Model):
                     ax.plot(self.X[which_data], self.likelihood.Y[which_data], 'kx', mew=1.5)
                 else:
                     m, v = self._raw_predict(Xnew, which_parts=which_parts, full_cov=True)
+                    v = v.reshape(m.size,-1) if len(v.shape)==3 else v
                     Ysim = np.random.multivariate_normal(m.flatten(), v, samples)
                     gpplot(Xnew, m, m - 2 * np.sqrt(np.diag(v)[:, None]), m + 2 * np.sqrt(np.diag(v))[:, None, ], axes=ax)
                     for i in range(samples):
                         ax.plot(Xnew, Ysim[i, :], Tango.colorsHex['darkBlue'], linewidth=0.25)
+
                 ax.plot(self.X[which_data], self.likelihood.Y[which_data], 'kx', mew=1.5)
                 ax.set_xlim(xmin, xmax)
                 ymin, ymax = min(np.append(self.likelihood.Y, m - 2 * np.sqrt(np.diag(v)[:, None]))), max(np.append(self.likelihood.Y, m + 2 * np.sqrt(np.diag(v)[:, None])))
                 ymin, ymax = ymin - 0.1 * (ymax - ymin), ymax + 0.1 * (ymax - ymin)
                 ax.set_ylim(ymin, ymax)
 
-        elif self.X.shape[1] == 2 and not hasattr(self,'multioutput'):
+                if hasattr(self,'Z'):
+                    Zu = self.Z * self._Xscale + self._Xoffset
+                    ax.plot(Zu, np.zeros_like(Zu) + ax.get_ylim()[0], 'r|', mew=1.5, markersize=12)
+
+            elif self.X.shape[1] == 2:
                 resolution = resolution or 50
                 Xnew, xmin, xmax, xx, yy = x_frame2D(self.X, plot_limits, resolution)
                 m, v = self._raw_predict(Xnew, which_parts=which_parts)
@@ -117,9 +125,12 @@ class GPBase(Model):
                 ax.set_xlim(xmin[0], xmax[0])
                 ax.set_ylim(xmin[1], xmax[1])
 
+            else:
+                raise NotImplementedError, "Cannot define a frame with more than two input dimensions"
+        else:
+            assert self.num_outputs > output, 'The model has only %s outputs.' %self.num_outputs
 
-        elif self.X.shape[1] == 2 and hasattr(self,'multioutput'):
+            if self.X.shape[1] == 2:
-            output -= 1
                 assert self.num_outputs >= output, 'The model has only %s outputs.' %self.num_outputs
                 Xu = self.X[self.X[:,-1]==output ,0:1]
                 Xnew, xmin, xmax = x_frame1D(Xu, plot_limits=plot_limits)
@@ -130,6 +141,7 @@ class GPBase(Model):
                     ax.plot(Xu[which_data], self.likelihood.Y[self.likelihood.index==output][:,None], 'kx', mew=1.5)
                 else:
                     m, v = self._raw_predict_single_output(Xnew, output=output, which_parts=which_parts, full_cov=True)
+                    v = v.reshape(m.size,-1) if len(v.shape)==3 else v
                     Ysim = np.random.multivariate_normal(m.flatten(), v, samples)
                     gpplot(Xnew, m, m - 2 * np.sqrt(np.diag(v)[:, None]), m + 2 * np.sqrt(np.diag(v))[:, None, ], axes=ax)
                     for i in range(samples):
@@ -139,19 +151,19 @@ class GPBase(Model):
                 ymin, ymax = ymin - 0.1 * (ymax - ymin), ymax + 0.1 * (ymax - ymin)
                 ax.set_ylim(ymin, ymax)
 
+            elif self.X.shape[1] == 3:
+                raise NotImplementedError, "Plots not implemented for multioutput models with 2D inputs...yet"
+                assert self.num_outputs >= output, 'The model has only %s outputs.' %self.num_outputs
+
+            else:
+                raise NotImplementedError, "Cannot define a frame with more than two input dimensions"
+
         if hasattr(self,'Z'):
             Zu = self.Z[self.Z[:,-1]==output,:]
             Zu = self.Z * self._Xscale + self._Xoffset
             Zu = self.Z[self.Z[:,-1]==output ,0:1] #??
             ax.plot(Zu, np.zeros_like(Zu) + ax.get_ylim()[0], 'r|', mew=1.5, markersize=12)
 
-        elif self.X.shape[1] == 3 and hasattr(self,'multioutput'):
-            raise NotImplementedError, "Plots not implemented for multioutput models with 2D inputs...yet"
-            output -= 1
-            assert self.num_outputs >= output, 'The model has only %s outputs.' %self.num_outputs
-
-        else:
-            raise NotImplementedError, "Cannot define a frame with more than two input dimensions"
 
     def plot(self, plot_limits=None, which_data='all', which_parts='all', resolution=None, levels=20, samples=0, fignum=None, ax=None, output=None, fixed_inputs=[], linecol=Tango.colorsHex['darkBlue'],fillcol=Tango.colorsHex['lightBlue']):
         """
@@ -203,7 +215,7 @@ class GPBase(Model):
             if plotdims == 1:
                 resolution = resolution or 200
 
-                Xu = self.X * self._Xscale + self._Xoffset # NOTE self.X are the normalized values now
+                Xu = self.X * self._Xscale + self._Xoffset #NOTE self.X are the normalized values now
 
                 fixed_dims = np.array([i for i,v in fixed_inputs])
                 freedim = np.setdiff1d(np.arange(self.input_dim),fixed_dims)