Contrast

nipy/modalities/fmri/design.py

@@ -79,7 +79,6 @@ def natural_spline(tvals, knots=None, order=3, intercept=True):
     f = formulae.natural_spline(t, knots=knots, order=order, intercept=intercept)
     return f.design(tvals, return_float=True)
 
-
 def event_design(event_spec, t, order=2, hrfs=(glover,),
                  level_contrasts=False):
     """ Create design matrix from event specification `event_spec`
@@ -132,6 +131,11 @@ def event_design(event_spec, t, order=2, hrfs=(glover,),
                 names = [c[0] for c in comb]
                 fs = [c[1].main_effect for c in comb]
                 e_contrasts[":".join(names)] = np.product(fs).design(event_spec)
+    else:
+        # only one factor, produce the main effect
+        field = fields[0]
+        factor = e_factors[0]
+        e_contrasts[field] = factor.main_effect.design(event_spec)
 
     e_contrasts['constant'] = formulae.I.design(event_spec)
 
@@ -161,7 +165,6 @@ def event_design(event_spec, t, order=2, hrfs=(glover,),
     X_t, c_t = t_formula.design(tval, contrasts=t_contrasts)
     return X_t, c_t
 
-
 def block_design(block_spec, t, order=2, hrfs=(glover,),
                  convolution_padding=5.,
                  convolution_dt=0.02,
@@ -221,12 +224,18 @@ def block_design(block_spec, t, order=2, hrfs=(glover,),
     e_factors = [Factor(n, np.unique(block_spec[n])) for n in fields]
     e_formula = np.product(e_factors)
     e_contrasts = {}
+
     if len(e_factors) > 1:
         for i in range(1, order+1):
             for comb in combinations(zip(fields, e_factors), i):
                 names = [c[0] for c in comb]
                 fs = [c[1].main_effect for c in comb]
                 e_contrasts[":".join(names)] = np.product(fs).design(block_spec)
+    else:
+        # only one factor, produce the main effect
+        field = fields[0]
+        factor = e_factors[0]
+        e_contrasts[field + '_effect'] = factor.main_effect.design(event_spec)
 
     e_contrasts['constant'] = formulae.I.design(block_spec)
 
@@ -272,7 +281,6 @@ def block_design(block_spec, t, order=2, hrfs=(glover,),
     X_t, c_t = t_formula.design(tval, contrasts=t_contrasts)
     return X_t, c_t
 
-
 def stack2designs(old_X, new_X, old_contrasts={}, new_contrasts={}):
     """
     Add some columns to a design matrix that has contrasts matrices

nipy/modalities/fmri/design_matrix.py

@@ -275,15 +275,17 @@ def write_csv(self, path):
             writer.writerow(self.names)
             writer.writerows(self.matrix)
 
-    def show(self, rescale=True, ax=None):
+    def show(self, rescale=True, ax=None, cmap=None):
         """Visualization of a design matrix
 
         Parameters
         ----------
         rescale: bool, optional
-                 rescale columns magnitude for visualization or not
+                 rescale columns magnitude for visualization or not.
         ax: axis handle, optional
-            Handle to axis onto which we will draw design matrix
+            Handle to axis onto which we will draw design matrix.
+        cmap: colormap, optional
+            Matplotlib colormap to use, passed to `imshow`.
 
         Returns
         -------
@@ -299,7 +301,8 @@ def show(self, rescale=True, ax=None):
             plt.figure()
             ax = plt.subplot(1, 1, 1)
 
-        ax.imshow(x, interpolation='Nearest', aspect='auto')
+        ax.imshow(x, interpolation='Nearest', aspect='auto',
+                  cmap=cmap)
         ax.set_label('conditions')
         ax.set_ylabel('scan number')
 
@@ -308,6 +311,44 @@ def show(self, rescale=True, ax=None):
             ax.set_xticklabels(self.names, rotation=60, ha='right')
         return ax
 
+    def show_contrast(self, contrast, ax=None, cmap=None):
+        """
+        Plot a contrast for a design matrix.
+
+        Parameters
+        ----------
+        contrast : np.float
+            Array forming contrast with respect to the design matrix.
+        ax: axis handle, optional
+            Handle to axis onto which we will draw design matrix.
+        cmap: colormap, optional
+            Matplotlib colormap to use, passed to `imshow`.
+
+        Returns
+        -------
+        ax: axis handle
+        """
+        import matplotlib.pyplot as plt
+
+        contrast = np.atleast_2d(contrast)
+
+        # normalize the values per column for better visualization
+        if ax is None:
+            plt.figure()
+            ax = plt.subplot(1, 1, 1)
+
+        ax.imshow(contrast, interpolation='Nearest', 
+                  aspect='auto',
+                  cmap=cmap)
+        ax.set_label('conditions')
+        ax.set_yticks(range(contrast.shape[0]))
+        ax.set_yticklabels([])
+
+        if self.names is not None:
+            ax.set_xticks(range(len(self.names)))
+            ax.set_xticklabels(self.names, rotation=60, ha='right')
+        return ax
+
 
 def make_dmtx(frametimes, paradigm=None, hrf_model='canonical',
               drift_model='cosine', hfcut=128, drift_order=1, fir_delays=[0],

nipy/modalities/fmri/tests/test_dmtx.py

@@ -71,6 +71,23 @@ def test_show_dmtx():
     ax = DM.show()
     assert (ax is not None)
 
+    # test the colormap
+    ax = DM.show(cmap=matplotlib.pyplot.cm.gray)
+    assert (ax is not None)
+
+@dec.skipif(not have_mpl)
+def test_show_constrast():
+    # test that the show code indeed (formally) runs
+    frametimes = np.linspace(0, 127 * 1.,128)
+    DM = make_dmtx(frametimes, drift_model='polynomial', drift_order=3)
+    contrast = np.random.standard_normal((3, DM.matrix.shape[1]))
+    ax = DM.show_contrast(contrast)
+    assert (ax is not None)
+
+    # test the colormap
+    ax = DM.show_contrast(contrast, cmap=matplotlib.pyplot.cm.gray)
+    assert (ax is not None)
+
 
 def test_dmtx0():
     # Test design matrix creation when no paradigm is provided