[med-svn] [python-mne] 138/376: computing single_trial_power in parallel

[Yaroslav Halchenko]

This is an automated email from the git hooks/post-receive script.

yoh pushed a commit to annotated tag v0.1
in repository python-mne.

commit 576186356967bb3f411feab3dbf56df733d07b89
Author: Alexandre Gramfort <alexandre.gramfort at inria.fr>
Date:   Tue Mar 15 15:37:33 2011 -0400

    computing single_trial_power in parallel
---
 .../plot_cluster_1samp_test_time_frequency.py      |  2 +-
 mne/time_frequency/tfr.py                          | 86 +++++++++++++++++++---
 2 files changed, 77 insertions(+), 11 deletions(-)

diff --git a/examples/stats/plot_cluster_1samp_test_time_frequency.py b/examples/stats/plot_cluster_1samp_test_time_frequency.py
index 34cb361..b562dc8 100644
--- a/examples/stats/plot_cluster_1samp_test_time_frequency.py
+++ b/examples/stats/plot_cluster_1samp_test_time_frequency.py
@@ -66,7 +66,7 @@ frequencies = np.arange(7, 30, 3) # define frequencies of interest
 Fs = raw.info['sfreq'] # sampling in Hz
 epochs_power = single_trial_power(data, Fs=Fs,
                                    frequencies=frequencies,
-                                   n_cycles=3, use_fft=False)
+                                   n_cycles=3, use_fft=False, n_jobs=1)
 
 epochs_power = epochs_power[:,0,:,:] # only 1 channel to get a 3D matrix
 # do ratio with baseline power:
diff --git a/mne/time_frequency/tfr.py b/mne/time_frequency/tfr.py
index 7f2593b..1678211 100644
--- a/mne/time_frequency/tfr.py
+++ b/mne/time_frequency/tfr.py
@@ -168,6 +168,42 @@ def cwt_morlet(X, Fs, freqs, use_fft=True, n_cycles=7.0):
 
     return tfrs
 
+def cwt(X, Ws, use_fft=True, mode='same'):
+    """Compute time freq decomposition with continuous wavelet transform
+
+    Parameters
+    ----------
+    X : array of shape [n_signals, n_times]
+        signals (one per line)
+
+    Ws : list of array
+        Wavelets time series
+
+    use_fft : bool
+        Use FFT for convolutions
+
+    mode : 'same' | 'valid' | 'full'
+        Convention for convolution
+
+    Returns
+    -------
+    tfr : 3D array
+        Time Frequency Decompositions (n_signals x n_frequencies x n_times)
+    """
+    n_signals, n_times = X.shape
+    n_frequencies = len(Ws)
+
+    if use_fft:
+        coefs = _cwt_fft(X, Ws, mode)
+    else:
+        coefs = _cwt_convolve(X, Ws, mode)
+
+    tfrs = np.empty((n_signals, n_frequencies, n_times), dtype=np.complex)
+    for k, tfr in enumerate(coefs):
+        tfrs[k] = tfr
+
+    return tfrs
+
 
 def _time_frequency(X, Ws, use_fft):
     """Aux of time_frequency for parallel computing over channels
@@ -191,27 +227,57 @@ def _time_frequency(X, Ws, use_fft):
     return psd, plf
 
 
-def single_trial_power(epochs, Fs, frequencies, use_fft=True, n_cycles=7):
+def single_trial_power(epochs, Fs, frequencies, use_fft=True, n_cycles=7,
+                       n_jobs=1):
     """Compute time-frequency power on single epochs
+
+    Parameters
+    ----------
+    epochs : instance of Epochs | 3D array
+        The epochs
+    Fs : float
+        Sampling rate
+    frequencies : array-like
+        The frequencies
+    use_fft : bool
+        Use the FFT for convolutions or not.
+    n_cycles : float
+        The number of cycles in the Morlet wavelet
+    n_jobs : int
+        The number of epochs to process at the same time
+
+    Returns
+    -------
+    power : list of 2D array
+        Each element of the list the the power estimate for an epoch.
     """
+    mode = 'same'
     n_frequencies = len(frequencies)
     n_epochs, n_channels, n_times = epochs.shape
 
     # Precompute wavelets for given frequency range to save time
     Ws = morlet(Fs, frequencies, n_cycles=n_cycles)
 
+    try:
+        from scikits.learn.externals.joblib import Parallel, delayed
+        parallel = Parallel(n_jobs)
+        my_cwt = delayed(cwt)
+    except ImportError:
+        print "joblib not installed. Cannot run in parallel."
+        n_jobs = 1
+        my_cwt = cwt
+        parallel = list
+
     power = np.empty((n_epochs, n_channels, n_frequencies, n_times),
                      dtype=np.float)
-
-    mode = 'same'
-    if use_fft:
-        _cwt = _cwt_fft
+    if n_jobs == 1:
+        for k, e in enumerate(epochs):
+            power[k] = np.abs(cwt(e, Ws, mode))**2
     else:
-        _cwt = _cwt_convolve
-
-    for k, e in enumerate(epochs):
-        mode = 'same'
-        power[k] = np.abs(list(_cwt(e, Ws, mode)))**2
+        # Precompute tf decompositions in parallel
+        tfrs = parallel(my_cwt(e, Ws, use_fft, mode) for e in epochs)
+        for k, tfr in enumerate(tfrs):
+            power[k] = np.abs(tfr)**2
 
     return power
 

-- 
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-med/python-mne.git