Example script showing the ctps / ica bug in mne-python

plot_show_ica_bug.py

import mnedev as mne
from mnedev.preprocessing import ICA, create_ecg_epochs
from mnedev.datasets import sample

print(__doc__)
print mne.__version__

# Read and preprocess the data. Preprocessing consists of:
#
# - MEG channel selection
#
# - 1-30 Hz band-pass filter

data_path = sample.data_path()
raw_fname = data_path + '/MEG/sample/sample_audvis_filt-0-40_raw.fif'

raw = mne.io.read_raw_fif(raw_fname, preload=True)

# include annotations for testing
annotations = mne.Annotations([2.], [30.], 'test reject annot')
raw.annotations = annotations

raw.pick_types(meg=True, eeg=False, exclude='bads', stim=True)
raw.filter(1, 30, fir_design='firwin')
chop = raw.copy().crop(tmin=0., tmax=60.)

ica = ICA(n_components=0.95, method='fastica').fit(chop)

print '\nICA on raw between 0. to 60. with CTPS \n'
# apply the ica on raw (between start and stop) and identify ECG related components
# using ctps
ecg_inds, scores = ica.find_bads_ecg(raw, ch_name='MEG 0143', threshold=None,
                                     start=0., stop=60., l_freq=8, h_freq=16,
                                     method='ctps', reject_by_annotation=True,
                                     verbose=None)
print ecg_inds

print '\nICA on raw between 0. to 60 with Correlation. \n'
# apply the ica on raw (between start and stop) and identify ECG related components
# using correlation
ecg_inds, scores = ica.find_bads_ecg(raw, ch_name='MEG 0143', threshold=None,
                                     start=0., stop=60., l_freq=8, h_freq=16,
                                     method='correlation', reject_by_annotation=True,
                                     verbose=None)
print ecg_inds


print '\nICA on previously chopped data  with ctps\n'
# apply the ica on chop and identify the ECG related components using ctps
ecg_inds, scores = ica.find_bads_ecg(chop, ch_name='MEG 0143', threshold=None,
                                     start=None, stop=None, l_freq=8, h_freq=16,
                                     method='ctps', reject_by_annotation=True,
                                     verbose=None)
print ecg_inds