Python mne.find_events() Examples

def test_viz():
    """Test viz."""
    import matplotlib.pyplot as plt

    set_matplotlib_defaults(plt)

    events = mne.find_events(raw)
    picks = mne.pick_channels(raw.info['ch_names'],
                              ['MEG 2443', 'MEG 2442', 'MEG 2441'])
    epochs = mne.Epochs(raw, events, picks=picks, baseline=(None, 0),
                        reject=None, preload=True,
                        event_id={'1': 1, '2': 2, '3': 3, '4': 4})
    bad_epochs_idx = [0, 1, 3]
    n_epochs, n_channels, _ = epochs.get_data().shape
    bad_epochs = np.zeros(n_epochs, dtype=bool)
    bad_epochs[bad_epochs_idx] = True

    labels = np.zeros((n_epochs, n_channels))
    reject_log = autoreject.RejectLog(bad_epochs, labels, epochs.ch_names)
    reject_log.plot_epochs(epochs)
    reject_log.plot()
    reject_log.plot(orientation='horizontal')
    pytest.raises(ValueError, reject_log.plot_epochs, epochs[:2])
    pytest.raises(ValueError, reject_log.plot, 'down')
    plt.close('all') 

def createStimChannel(events):
    """
    Create stim channel from events.

    Parameters
    ----------
    events : instance of pandas.core.DataFrame
        Dataframe containing list of events obtained with mne.find_events(raw)
       .

    Returns:

    stim : instance of pandas.core.series.Series
        Series containing the stimulus channel reconstructed from events.
    """

    newEvents = events.copy()
    newEvents[0] = newEvents[0].round(0).astype('int')
    chan = np.zeros(int(events.iloc[-1, 0]))
    chan[newEvents.iloc[:-2, 0]] = 8
    stim = pd.Series(chan, columns = ['STI 014'])
    return stim 

def downsampleEvents(events, oldFS, newFS):
    """
    Modify the timestamps of events to match a new sampling frequency.

    Parameters
    ----------
    events : instance of pandas.core.DataFrame
        Dataframe containing list of events obtained with mne.find_events(raw)
       .
    oldFS : float
        The sampling frequency of the input events.
    newFS : float
        The sampling frequency to the output events.

    Returns:

    newEvents : instance of pandas.DataFrame
        DataFrame containing the downsampled events.
    """

    newEvents = events.copy()
    newEvents[0] = (events[0]/oldFS)*newFS
    newEvents[0] = newEvents[0].round(0).astype('int')
    return newEvents 

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_ransac():
    """Some basic tests for ransac."""
    event_id = {'Visual/Left': 3}
    tmin, tmax = -0.2, 0.5

    events = mne.find_events(raw)
    epochs = mne.Epochs(raw, events, event_id, tmin, tmax,
                        baseline=(None, 0), decim=8,
                        reject=None, preload=True)
    # normal case
    picks = mne.pick_types(epochs.info, meg='mag', eeg=False, stim=False,
                           eog=False, exclude=[])

    ransac = Ransac(picks=picks)
    epochs_clean = ransac.fit_transform(epochs)
    assert len(epochs_clean) == len(epochs)
    # Pass numpy instead of epochs
    X = epochs.get_data()
    pytest.raises(AttributeError, ransac.fit, X)
    #
    # should not contain both channel types
    picks = mne.pick_types(epochs.info, meg=True, eeg=False, stim=False,
                           eog=False, exclude=[])
    ransac = Ransac(picks=picks)
    pytest.raises(ValueError, ransac.fit, epochs)
    #
    # should not contain other channel types.
    picks = mne.pick_types(raw.info, meg=False, eeg=True, stim=True,
                           eog=False, exclude=[])
    ransac = Ransac(picks=picks)
    pytest.raises(ValueError, ransac.fit, epochs) 

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 load_one_record(data_type, subject, run_index, sfreq=300, epoch=None,
                    filter_params=[5., None], n_jobs=1):
    # Load the record and correct the sensor space to get proper visualization
    print(f"subject={subject}, data_type={data_type}, run_index={run_index}, "
          f"hcp_path={HCP_DIR}")
    raw = hcp.read_raw(subject, data_type=data_type, run_index=run_index,
                       hcp_path=HCP_DIR, verbose=0)
    raw.load_data()
    hcp.preprocessing.map_ch_coords_to_mne(raw)
    raw.pick_types(meg='mag', eog=False, stim=True)

    # filter the electrical and low frequency components
    raw.notch_filter([60, 120], n_jobs=n_jobs)
    raw.filter(*filter_params, n_jobs=n_jobs)

    # Resample to the requested sfreq
    if sfreq is not None:
        raw.resample(sfreq=sfreq, n_jobs=n_jobs)

    events = mne.find_events(raw)
    raw.pick_types(meg='mag', stim=False)
    events[:, 0] -= raw.first_samp

    # Deep copy before modifying info to avoid issues when saving EvokedArray
    info = deepcopy(raw.info)
    info['events'] = events
    info['event_id'] = np.unique(events[:, 2])

    # Return the data
    return raw.get_data(), info 

def discriminateEvents(events, threshold):
    """
    Discriminate triggers when different kind of events are on the same channel.
    A time threshold is used to determine if two events are from the same trial.

    Parameters
    ----------
    events : instance of pandas.core.DataFrame
        Dataframe containing the list of events obtained with
        mne.find_events(raw).
    threshold : float
        Time threshold in milliseconds. Keeps an event if the time difference
        with the next one is superior than threshold.

    Returns:

    newData : instance of pandas.series.Series
        List of trial number filling the requirements.
    """

    # calculate the rolling difference (between n and n+1)
    events['diff'] = events[0].diff()
    # replace the nan with the first value
    events['diff'].iloc[0] = events.iloc[0, 0]
    # select events with time distance superior to threshold
    events = events[events['diff']>threshold]
    events = events.reset_index(drop=True)
    del events['diff']
    return events 

def getEvents(raw, eventCode, shortest_event=None):
    """
    Get the events corresponding to `eventCode`.

    Parameters
    ----------
    raw : instance of mne.io.edf.edf.RawEDF
        RawEDF object from the MNE library containing data from the .bdf files.
    eventCode : int
        Code corresponding to a specific events. For instance, with a biosemi
        device, the triggers are coded 65284, 65288 and 65296 respectively on
        the first, second and third channel.

    Returns:

    startEvents : instance of pandas.core.DataFrame
        Dataframe containing the list of timing corresponding to the event code
        in the first column. The second column contains the code before the event
        and the third the code of the selected event.
    """
    if shortest_event:
      events = mne.find_events(raw, shortest_event=shortest_event)
    else:
      events = mne.find_events(raw)

    eventsDf = pd.DataFrame(events)

    # Keep only event corresponding to the event code
    startEvents = eventsDf.loc[eventsDf[2]==eventCode]

    startEvents = startEvents.set_index([np.arange(len(startEvents))])
    startEvents.index.name = 'start'
    return startEvents 

def find_events(self, stim_channel, consecutive=True, initial_event=True,
                    uint_cast=True, min_duration=0, shortest_event=0):
        """Find events in raw data."""
        events = mne.find_events(self.current["data"],
                                 stim_channel=stim_channel,
                                 consecutive=consecutive,
                                 initial_event=initial_event,
                                 uint_cast=uint_cast,
                                 min_duration=min_duration,
                                 shortest_event=shortest_event)
        if events.shape[0] > 0:  # if events were found
            self.current["events"] = events
            self.history.append("events = mne.find_events(data)") 

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 _read_events(events_data, event_id, raw, ext, verbose=None):
    """Read in events data.

    Parameters
    ----------
    events_data : str | array | None
        The events file. If a string, a path to the events file. If an array,
        the MNE events array (shape n_events, 3). If None, events will be
        inferred from the stim channel using `find_events`.
    event_id : dict
        The event id dict used to create a 'trial_type' column in events.tsv,
        mapping a description key to an integer valued event code.
    raw : instance of Raw
        The data as MNE-Python Raw object.
    ext : str
        The extension of the original data file.
    verbose : bool | str | int | None
        If not None, override default verbose level (see :func:`mne.verbose`).

    Returns
    -------
    events : array, shape = (n_events, 3)
        The first column contains the event time in samples and the third
        column contains the event id. The second column is ignored for now but
        typically contains the value of the trigger channel either immediately
        before the event or immediately after.

    """
    if isinstance(events_data, str):
        events = read_events(events_data, verbose=verbose).astype(int)
    elif isinstance(events_data, np.ndarray):
        if events_data.ndim != 2:
            raise ValueError('Events must have two dimensions, '
                             'found %s' % events_data.ndim)
        if events_data.shape[1] != 3:
            raise ValueError('Events must have second dimension of length 3, '
                             'found %s' % events_data.shape[1])
        events = events_data
    elif 'stim' in raw:
        events = find_events(raw, min_duration=0.001, initial_event=True,
                             verbose=verbose)
    elif ext in ['.vhdr', '.set'] and check_version('mne', '0.18'):
        events, event_id = events_from_annotations(raw, event_id,
                                                   verbose=verbose)
    else:
        warn('No events found or provided. Please make sure to'
             ' set channel type using raw.set_channel_types'
             ' or provide events_data.')
        events = None
    return events, event_id 

def getTrialsAverage(data, events, trialDur=None, trialNumList=None,
    baselineDur=0.1, normalize=False, fs=2048., startOffset=0, noiseAve=None):
    """
    Get the average across trials (from `trialNumList`) based on time-locking
    provided by `events`.

    Parameters
    ----------
    data : instance of pandas.core.DataFrame
        Data containing values across time (not epoched).
    events : instance of pandas.core.DataFrame
        Dataframe containing the list of events obtained with
        mne.find_events(raw).
    trialDur : float | None
        Trial duration in seconds.
    trialNumList : array-like of int | None
        List of all trials to use. If None, all trials are taken.
    baselineDur : float, defaults to 0.1
        Duration of the baseline in seconds. If normalize is True, normalize
        each electrode with a baseline of duration `baselineDur`.
    normalize : bool, defaults to False
        If True data will be normalized.
    fs : float
        Sampling frequency of data in Hz.

    Returns:

    meanTrials : instance of pandas.series.Series
        Series containing the averaged values across trials.
    allTrials : instance of pandas.core.DataFrame
        Dataframe containing the values of each trial (1 column = 1 trial).
    """

    # Calculate average across trials for this Series
    if (trialNumList is None):
        trialNumList = [events.shape[0]]

    allTrials = pd.DataFrame()
    for trialNum in trialNumList:
        trialData = getTrialData(data, events=events, trialNum=trialNum,
            baselineDur=baselineDur, trialDur=trialDur, fs=fs,
            startOffset=startOffset)
        if normalize:
            trialData = normalizeFromBaseline(trialData,
                baselineDur=baselineDur, fs=fs)

        trialData = trialData.reset_index(drop=True)

        if noiseAve is not None:
            trialData = trialData - noiseAve

        allTrials['trial%d' % trialNum] = pd.Series(trialData)

    # convert baselineDur in frames
    start = int(np.round(baselineDur *fs))

    # Change index to have x-axis in seconds
    allTrials = allTrials.set_index(np.arange(-start, allTrials.shape[0]-start)/fs)

    meanTrials = allTrials.mean(axis=1)
    return meanTrials, allTrials 

def plotERPElectrodes(data, trialNumList, events, trialDur=None, fs=2048.,
    baselineDur=0.1, electrodes='Fp1', normalize=False, facet=False,
    startOffset=0):
    """
    Plot the ERP (average across trials time-locked to specific events) of
    each electrode as single lines on the same figure with or without facetting.

    Parameters
    ----------
    data : instance of pandas.core.DataFrame
        Data containing the time series to transform and plot. Each column is an
        electrode.
    trialNumList : array-like of int
        List of all trials to use to compute the FFT.
    events : instance of pandas.core.DataFrame
        Dataframe containing the list of events obtained with
        mne.find_events(raw).
    trialDur : float
        Trial duration in seconds.
    fs : float
        Sampling frequency of data in Hz.
    baselineDur : float, defaults to 0.1
        Duration of the baseline in seconds. If normalize is True, normalize
        each electrode with a baseline of duration `baselineDur`.
    electrodes : int | array-like of int, default to 'Fp1'
        List of electrodes to use to compute the FFT.
    normalize : bool, defaults to False
        If True data will be normalized.
    facet : bool, default to False
        If True, each electrode will be plotted on a different facet.

    Returns:

    fig : instance of matplotlib.figure.Figure
        The figure of the ERP.
    """

    print 'Average of %d trials' % len(trialNumList)
    meanTrials = pd.DataFrame()
    for electrode in electrodes:
        meanTrials[electrode], allTrials = getTrialsAverage(data=data[electrode],
            events=events, trialDur=trialDur, trialNumList=trialNumList,
            baselineDur=baselineDur, normalize=normalize, fs=fs, startOffset=startOffset)

    if (facet):
        print 'Faceting...'
        meanTrials.plot(subplots=True)
    else:
        plt.figure()
        plt.plot(meanTrials)
        plt.axvline(x=0, color='grey', linestyle='dotted')
        plt.axvspan(-baselineDur, 0, alpha=0.3, color='grey')
        plt.xlabel('Time (s)')
        # plt.legend(meanTrials.columns, bbox_to_anchor=(1, 1), ncol=4)
        plt.show()