Python mne.pick_types() Examples

def interpolate_bads(inst, picks, dots=None, reset_bads=True, mode='accurate'):
    """Interpolate bad MEG and EEG channels."""
    import mne
    # to prevent cobyla printf error
    # XXX putting to critical for now unless better solution
    # emerges
    verbose = mne.set_log_level('CRITICAL', return_old_level=True)

    eeg_picks = set(pick_types(inst.info, meg=False, eeg=True, exclude=[]))
    eeg_picks_interp = [p for p in picks if p in eeg_picks]
    if len(eeg_picks_interp) > 0:
        _interpolate_bads_eeg(inst, picks=eeg_picks_interp)

    meg_picks = set(pick_types(inst.info, meg=True, eeg=False, exclude=[]))
    meg_picks_interp = [p for p in picks if p in meg_picks]
    if len(meg_picks_interp) > 0:
        _interpolate_bads_meg_fast(inst, picks=meg_picks_interp,
                                   dots=dots, mode=mode)

    if reset_bads is True:
        inst.info['bads'] = []

    mne.set_log_level(verbose)

    return inst 

def calc_noise_epoches_from_empty_room(events_id, data_raw_fname, empty_room_raw_fname, from_t, to_t,
        overwrite_epochs=False):
    from mne.event import make_fixed_length_events
    from mne.io import Raw

    epochs_noise_dic = {}
    epochs_noise_fnames = [get_cond_fname(EPO_NOISE, event) for event in events_id.keys()]
    if np.all([os.path.isfile(fname) for fname in epochs_noise_fnames]) and not overwrite_epochs:
        for event in events_id.keys():
            epochs_noise_dic[event] = mne.read_epochs(get_cond_fname(EPO_NOISE, event))
    else:
        raw = Raw(data_raw_fname)
        raw_noise = Raw(empty_room_raw_fname)
        # raw_noise.info['bads'] = ['MEG0321']  # 1 bad MEG channel
        picks = mne.pick_types(raw.info, meg=True)#, exclude='bads')
        events_noise = make_fixed_length_events(raw_noise, 1)
        epochs_noise = mne.Epochs(raw_noise, events_noise, 1, from_t,
            to_t, proj=True, picks=picks, baseline=None, preload=True)
        for event, event_id in events_id.items():
            # then make sure the number of epochs is the same
            epochs = mne.read_epochs(get_cond_fname(EPO, event))
            epochs_noise_dic[event] = epochs_noise[:len(epochs.events)]
            epochs_noise_dic[event].save(get_cond_fname(EPO_NOISE, event))
    return epochs_noise_dic 

def _handle_picks(info, picks):
    """Pick the data channls or return picks."""
    if picks is None:
        out = mne.pick_types(
            info, meg=True, eeg=True, ref_meg=False, exclude='bads')
    else:
        out = picks
    return out 

def test_global_autoreject():
    """Test global autoreject."""
    event_id = None
    tmin, tmax = -0.2, 0.5
    events = mne.find_events(raw)

    picks = mne.pick_types(raw.info, meg=True, eeg=True, stim=False,
                           eog=True, exclude=[])
    # raise error if preload is false
    epochs = mne.Epochs(raw, events, event_id, tmin, tmax,
                        picks=picks, baseline=(None, 0),
                        reject=None, preload=False)

    # Test get_rejection_thresholds.
    reject1 = get_rejection_threshold(epochs, decim=1, random_state=42)
    reject2 = get_rejection_threshold(epochs, decim=1, random_state=42)
    reject3 = get_rejection_threshold(epochs, decim=2, random_state=42)
    tols = dict(eeg=5e-6, eog=5e-6, grad=10e-12, mag=5e-15)
    assert reject1, isinstance(reject1, dict)
    for key, value in list(reject1.items()):
        assert reject1[key] == reject2[key]
        assert abs(reject1[key] - reject3[key]) < tols[key]

    reject = get_rejection_threshold(epochs, decim=4, ch_types='eeg')
    assert 'eog' not in reject
    assert 'eeg' in reject
    pytest.raises(ValueError, get_rejection_threshold, epochs,
                  decim=4, ch_types=5) 

def test_io():
    """Test IO functionality."""
    event_id = None
    tmin, tmax = -0.2, 0.5
    events = mne.find_events(raw)
    savedir = _TempDir()
    fname = op.join(savedir, 'autoreject.hdf5')

    include = [u'EEG %03d' % i for i in range(1, 45, 3)]
    picks = mne.pick_types(raw.info, meg=False, eeg=False, stim=False,
                           eog=True, include=include, exclude=[])

    # raise error if preload is false
    epochs = mne.Epochs(raw, events, event_id, tmin, tmax,
                        picks=picks, baseline=(None, 0), decim=4,
                        reject=None, preload=True)[:10]
    ar = AutoReject(cv=2, random_state=42, n_interpolate=[1],
                    consensus=[0.5], verbose=False)
    ar.save(fname)  # save without fitting

    # check that fit after saving is the same as fit
    # without saving
    ar2 = read_auto_reject(fname)
    ar.fit(epochs)
    ar2.fit(epochs)
    assert np.sum([ar.threshes_[k] - ar2.threshes_[k]
                   for k in ar.threshes_.keys()]) == 0.

    pytest.raises(ValueError, ar.save, fname)
    ar.save(fname, overwrite=True)
    ar3 = read_auto_reject(fname)
    epochs_clean1, reject_log1 = ar.transform(epochs, return_log=True)
    epochs_clean2, reject_log2 = ar3.transform(epochs, return_log=True)
    assert_array_equal(epochs_clean1.get_data(), epochs_clean2.get_data())
    assert_array_equal(reject_log1.labels, reject_log2.labels) 

def test_utils():
    """Test utils."""
    event_id = {'Visual/Left': 3}
    tmin, tmax = -0.2, 0.5
    events = mne.find_events(raw)
    picks = mne.pick_channels(raw.info['ch_names'],
                              ['MEG 2443', 'MEG 2442', 'MEG 2441'])
    epochs = mne.Epochs(raw, events, event_id, tmin, tmax,
                        picks=picks, baseline=(None, 0),
                        reject=None, preload=True)

    this_epoch = epochs.copy()
    epochs_clean = clean_by_interp(this_epoch)
    assert_array_equal(this_epoch.get_data(), epochs.get_data())
    pytest.raises(AssertionError, assert_array_equal, epochs_clean.get_data(),
                  this_epoch.get_data())

    picks_meg = mne.pick_types(evoked.info, meg='grad', eeg=False, exclude=[])
    picks_eeg = mne.pick_types(evoked.info, meg=False, eeg=True, exclude=[])
    picks_bad_meg = mne.pick_channels(evoked.ch_names, include=['MEG 2443'])
    picks_bad_eeg = mne.pick_channels(evoked.ch_names, include=['EEG 053'])
    evoked_orig = evoked.copy()
    for picks, picks_bad in zip([picks_meg, picks_eeg],
                                [picks_bad_meg, picks_bad_eeg]):
        evoked_autoreject = interpolate_bads(evoked, picks=picks,
                                             reset_bads=False)
        evoked.interpolate_bads(reset_bads=False)
        assert_array_equal(evoked.data[picks_bad],
                           evoked_autoreject.data[picks_bad])
        pytest.raises(AssertionError, assert_array_equal,
                      evoked_orig.data[picks_bad], evoked.data[picks_bad])

    # test that autoreject EEG interpolation code behaves the same as MNE
    evoked_ar = evoked_orig.copy()
    evoked_mne = evoked_orig.copy()

    origin = _check_origin('auto', evoked_ar.info)
    _interpolate_bads_eeg(evoked_ar, picks=None)
    mne.channels.interpolation._interpolate_bads_eeg(evoked_mne, origin=origin)
    assert_array_almost_equal(evoked_ar.data, evoked_mne.data) 

def test_interpolate_bads():
    """Test interpolate bads."""
    event_id = None
    events = mne.find_events(raw)
    tmin, tmax = -0.2, 0.5
    for ii, ch_name in enumerate(raw.info['ch_names'][:14]):
        raw.set_channel_types({ch_name: 'bio'})
        raw.rename_channels({ch_name: 'BIO%02d' % ii})

    picks = mne.pick_types(raw.info, meg='grad', eeg=False, eog=False)
    epochs = mne.Epochs(raw, events, event_id, tmin, tmax,
                        baseline=(None, 0), decim=10,
                        reject=None, preload=True)[:10]
    epochs.info['bads'] = ['MEG 2212']
    interpolate_bads(epochs, picks) 

def _calc_sensors_power_parallel(p):
    from mne.time_frequency import tfr_array_morlet

    subject, window_fname, modality, inverse_method, bad_channels, downsample, high_gamma_max, overwrite = p

    root_dir = utils.make_dir(op.join(EEG_DIR if modality == 'eeg' else MEG_DIR, subject))
    output_fname_template = op.join(root_dir, '{}-epilepsy-{}-{}-{}-sensors_power.npy'.format(
        subject, inverse_method, modality, '{window}'))

    window = utils.namebase(window_fname)
    output_fname = output_fname_template.format(window=window)
    if op.isfile(output_fname) and not overwrite:
        print('{} already exist'.format(output_fname))
        return
    evoked = mne.read_evokeds(window_fname)[0]
    # evoked = evoked.resample(1000)

    if modality == 'eeg':
        picks = mne.pick_types(evoked.info, meg=False, eeg=True, exclude=bad_channels)
    elif modality == 'meg':
        picks = mne.pick_types(evoked.info, meg=True, eeg=False, exclude=bad_channels)
    elif modality == 'meeg':
        picks = mne.pick_types(evoked.info, meg=True, eeg=True, exclude=bad_channels)
    else:
        raise Exception('Wrong modality!')

    evoked_data = evoked.data[np.newaxis, picks, :]
    freqs, n_cycles = calc_morlet_freqs(evoked, high_gamma_max)
    powers = tfr_array_morlet(
        evoked_data, sfreq=evoked.info['sfreq'], freqs=freqs, n_cycles=n_cycles, output='power')
    powers = np.squeeze(powers)
    if powers.shape[2] % 2 == 1:
        powers = powers[:, :, :-1]
    if downsample > 1:
        powers = utils.downsample_3d(powers, downsample)
    powers_db = 10 * np.log10(powers)  # dB/Hz should be baseline corrected!!!
    print('Saving {}'.format(output_fname))
    np.save(output_fname, powers_db.astype(np.float16)) 

def init_data():
    data_path = sample.data_path()
    raw_fname = data_path + '/MEG/sample/sample_audvis_raw.fif'
    fname_inv = data_path + '/MEG/sample/sample_audvis-meg-oct-6-meg-inv.fif'
    tmin, tmax, event_id = -0.2, 0.5, 1

    # Setup for reading the raw data
    raw = io.read_raw_fif(raw_fname)
    events = mne.find_events(raw, stim_channel='STI 014')
    inverse_operator = read_inverse_operator(fname_inv)

    # Setting the label
    label = mne.read_label(data_path + '/MEG/sample/labels/Aud-lh.label')

    include = []
    raw.info['bads'] += ['MEG 2443', 'EEG 053']  # bads + 2 more

    # picks MEG gradiometers
    picks = mne.pick_types(raw.info, meg=True, eeg=False, eog=True,
                           stim=False, include=include, exclude='bads')

    # Load condition 1
    event_id = 1
    # Use linear detrend to reduce any edge artifacts
    epochs = mne.Epochs(raw, events, event_id, tmin, tmax, picks=picks,
                        baseline=(None, 0), reject=dict(grad=4000e-13, eog=150e-6),
                        preload=True, detrend=1)
    return epochs, inverse_operator, label 

def _calc_epoches(params):
    subject, events_id, tmin, tmax = params
    out_file = op.join(LOCAL_ROOT_DIR, 'epo', '{}_ecr_nTSSS_conflict-epo.fif'.format(subject))
    if not op.isfile(out_file):
        events = mne.read_events(op.join(REMOTE_ROOT_DIR, 'events', '{}_ecr_nTSSS_conflict-eve.fif'.format(subject)))
        raw = mne.io.Raw(op.join(REMOTE_ROOT_DIR, 'raw', '{}_ecr_nTSSS_raw.fif'.format(subject)), preload=False)
        picks = mne.pick_types(raw.info, meg=True)
        epochs = find_epoches(raw, picks, events, events_id, tmin=tmin, tmax=tmax)
        epochs.save(out_file)
    else:
        epochs = mne.read_epochs(out_file)
    return epochs 

def test_set(_bids_validate):
    """Test write_raw_bids conversion for EEGLAB data."""
    # standalone .set file with associated .fdt
    bids_root = _TempDir()
    data_path = op.join(testing.data_path(), 'EEGLAB')
    raw_fname = op.join(data_path, 'test_raw.set')

    raw = mne.io.read_raw_eeglab(raw_fname)

    # embedded - test mne-version assertion
    tmp_version = mne.__version__
    mne.__version__ = '0.16'
    with pytest.raises(ValueError, match='Your version of MNE is too old.'):
        write_raw_bids(raw, bids_basename, bids_root)
    mne.__version__ = tmp_version

    # proceed with the actual test for EEGLAB data
    write_raw_bids(raw, bids_basename, bids_root, overwrite=False)
    read_raw_bids(bids_basename=bids_basename, bids_root=bids_root, kind='eeg')

    with pytest.raises(TypeError, match="unexpected keyword argument 'foo'"):
        read_raw_bids(bids_basename=bids_basename, bids_root=bids_root,
                      kind='eeg', extra_params=dict(foo='bar'))

    with pytest.raises(FileExistsError, match="already exists"):  # noqa: F821
        write_raw_bids(raw, bids_basename, bids_root=bids_root,
                       overwrite=False)
    _bids_validate(bids_root)

    # check events.tsv is written
    # XXX: only from 0.18 onwards because events_from_annotations
    # is broken for earlier versions
    events_tsv_fname = op.join(bids_root, 'sub-' + subject_id,
                               'ses-' + session_id, 'eeg',
                               str(bids_basename) + '_events.tsv')
    if check_version('mne', '0.18'):
        assert op.exists(events_tsv_fname)

    # Also cover iEEG
    # We use the same data and pretend that eeg channels are ecog
    raw.set_channel_types({raw.ch_names[i]: 'ecog'
                           for i in mne.pick_types(raw.info, eeg=True)})
    bids_root = _TempDir()
    write_raw_bids(raw, bids_basename, bids_root)
    _bids_validate(bids_root)

    # test anonymize and convert
    if check_version('mne', '0.20') and check_version('pybv', '0.2.0'):
        output_path = _test_anonymize(raw, bids_basename)
        _bids_validate(output_path) 

def _interpolate_bads_eeg(inst, picks=None, verbose=None):
    """ Interpolate bad EEG channels.

    Operates in place.

    Parameters
    ----------
    inst : mne.io.Raw, mne.Epochs or mne.Evoked
        The data to interpolate. Must be preloaded.
    picks: np.ndarray, shape(n_channels, ) | list | None
        The channel indices to be used for interpolation.
    """
    from mne.bem import _fit_sphere
    from mne.utils import logger, warn
    from mne.channels.interpolation import _do_interp_dots
    from mne.channels.interpolation import _make_interpolation_matrix
    import numpy as np

    if picks is None:
        picks = pick_types(inst.info, meg=False, eeg=True, exclude=[])

    bads_idx = np.zeros(len(inst.ch_names), dtype=np.bool)
    goods_idx = np.zeros(len(inst.ch_names), dtype=np.bool)

    inst.info._check_consistency()
    bads_idx[picks] = [inst.ch_names[ch] in inst.info['bads'] for ch in picks]

    if len(picks) == 0 or bads_idx.sum() == 0:
        return

    goods_idx[picks] = True
    goods_idx[bads_idx] = False

    pos = inst._get_channel_positions(picks)

    # Make sure only good EEG are used
    bads_idx_pos = bads_idx[picks]
    goods_idx_pos = goods_idx[picks]
    pos_good = pos[goods_idx_pos]
    pos_bad = pos[bads_idx_pos]

    # test spherical fit
    radius, center = _fit_sphere(pos_good)
    distance = np.sqrt(np.sum((pos_good - center) ** 2, 1))
    distance = np.mean(distance / radius)
    if np.abs(1. - distance) > 0.1:
        warn('Your spherical fit is poor, interpolation results are '
             'likely to be inaccurate.')

    logger.info('Computing interpolation matrix from {0} sensor '
                'positions'.format(len(pos_good)))

    interpolation = _make_interpolation_matrix(pos_good, pos_bad)

    logger.info('Interpolating {0} sensors'.format(len(pos_bad)))
    _do_interp_dots(inst, interpolation, goods_idx, bads_idx) 

def test_select_vertices_in_sensor_range(fwd, src):
    """Test selecting and restricting vertices in the sensor range"""
    fwd_r = restrict_forward_to_vertices(fwd, ([1170, 1609], [2159]))

    verts = select_vertices_in_sensor_range(fwd_r, 0.05)
    assert_array_equal(verts[0], np.array([1170]))
    assert_array_equal(verts[1], np.array([]))
    # Test indices
    verts = select_vertices_in_sensor_range(fwd_r, 0.07, indices=True)
    assert_array_equal(verts, np.array([0, 1, 2]))

    # Test restricting
    fwd_rs = restrict_forward_to_sensor_range(fwd_r, 0.05)
    assert_array_equal(fwd_rs['src'][0]['vertno'], np.array([1170]))
    assert_array_equal(fwd_rs['src'][1]['vertno'], np.array([]))

    verts = select_vertices_in_sensor_range(fwd_r, 0.07)
    assert_array_equal(verts[0], np.array([1170, 1609]))
    assert_array_equal(verts[1], np.array([2159]))

    src_r = restrict_src_to_vertices(src, ([1170, 1609], [2159]))

    with pytest.raises(ValueError):  # info missing
        select_vertices_in_sensor_range(src_r, 0.07)
    info = _info()
    with pytest.raises(ValueError):  # trans missing
        select_vertices_in_sensor_range(src_r, 0.07, info=info)

    # Correct input
    trans = _trans()
    verts2 = select_vertices_in_sensor_range(src_r, 0.05, info=info,
                                             trans=trans)
    assert_array_equal(verts2[0], np.array([1170]))
    verts2 = select_vertices_in_sensor_range(src_r, 0.07, info=info,
                                             trans=trans)
    assert_array_equal(verts[0], verts2[0])
    assert_array_equal(verts[1], verts2[1])
    # Indices
    verts2 = select_vertices_in_sensor_range(src_r, 0.07, info=info,
                                             trans=trans,  indices=True)
    assert_array_equal(verts2, np.array([0, 1, 2]))
    # Try with only EEG
    info = mne.pick_info(info, sel=mne.pick_types(info, meg=False, eeg=True))
    verts2 = select_vertices_in_sensor_range(src_r, 0.05, info=info,
                                             trans=trans)
    assert_array_equal(verts2[0], np.array([1170, 1609]))
    assert_array_equal(verts2[1], np.array([2159])) 

def __init__(
        self, instance, montage_kind="standard_1020", low_cut=0.01, high_cut=50.0
    ):
        """Initialize the class."""
        # Make sure that we got an MNE object
        assert isinstance(instance, mne.io.BaseRaw)

        # The data that we are targeting
        # and a modifiable copy
        self.raw_mne = instance
        self.raw_copy = self.raw_mne.copy()

        # Set montage, pick data type, get data and transform to uVolts
        # We also filter all data at `low_cut` Hz highpass and obtain some data
        # bandpassed between `low_cut` and `high_cut` Hz.
        if LooseVersion(mne.__version__) < LooseVersion("0.20"):
            montage = mne.channels.read_montage(
                kind=montage_kind, ch_names=self.raw_copy.ch_names
            )
        else:
            montage = mne.channels.make_standard_montage(montage_kind)

        self.raw_copy.set_montage(montage)
        self.raw_copy.pick_types(eeg=True, stim=False)
        self.raw_copy.filter(
            l_freq=low_cut,
            h_freq=None,
            method="fir",
            fir_design="firwin",
            verbose=False,
        )
        self.x = self.raw_copy.get_data() * 1e6
        self.raw_copy.filter(
            l_freq=None,
            h_freq=high_cut,
            method="fir",
            fir_design="firwin",
            verbose=False,
        )
        self.x_bp = self.raw_copy.get_data() * 1e6
        self.ch_names = np.asarray(self.raw_copy.ch_names)
        self.n_chans = len(self.ch_names)
        self.signal_len = len(self.raw_copy.times)
        self.sfreq = self.raw_copy.info["sfreq"]
        self.chn_pos = self.raw_copy._get_channel_positions()

        # The identified bad channels
        self.bad_by_flat = []
        self.bad_by_nan = []
        self.bad_by_deviation = []
        self.bad_by_hf_noise = []
        self.bad_by_correlation = []
        self.bad_by_ransac = []