replacing data in the specific time range with NaN

[oconnell_E]

Hi everyone,

• MNE version: 1.0.3
• operating system: macOS 11.1

I’m looking for a way to extract EEG data in the certain time range and replace them with NaN.

The reason why I want to do it is that I want to get rid of TMS-induced artifacts (marked by red line) and interpolate the gap.

The data is first epoched -

epoched = mne.Epochs(raw,
                    events=events_TMS,
                    event_id=eventsid_TMS,
                    tmin=-0.2,
                    tmax=0.5, baseline=(-0.1, 0.05), preload = True)

• based on this event information. 10001 (named DC Correction) indicates the timing of TMS pulse.

eventsid_TMS
Out[9]: 
{'DC Correction/': 10001,
 'New Segment/': 99999,
 'Stimulus/S  1': 1,
 'Stimulus/S  3': 3}

(This is the array of events_TMS)

And the time of onset/offset of the TMS-induced artefacts are also specified. I want to first extract the data between the onset and offset of the artefacts to replace them with NAN.

tms_artefact_pre =0.005  # for rejecting artefact, before TMS pulse
tms_artefact_post=0.015   # for rejecting artefact, after TMS pulse

Could someone provide me a clue how to perform this on a “epoched data”?
Thank you so much.

[richard (Richard Höchenberger)]

Hello @oconnell_E and welcome to the forum!

Why don’t you just mark the respective data segments as “bad” on the raw data and throw them away? I don’t think there’s any use in interpolating them, as I assume all channels are affected by the TMS pulses, no?

Best wishes,
Richard

[oconnell_E]

Hi Richard,

thanks for your reply!
I understand what you suggested - but we’d like to avoid the discontinuity in the data as much as possible, that’s why I posted the question above.
I’d appreciate if you could provide any idea!
Thanks a lot.

[richard (Richard Höchenberger)]

So you don’t intend to epoch the data? I’m asking because in your example, you’re creating epochs around the TMS events. Could you briefly explain what you intend to so, so we can first assess whether it’s feasible?

For example, I don’t think MNE has any machinery for interpolation across time (only across space), unless I’m mistaken…

Best wishes,
Richard

[mscheltienne (Mathieu Scheltienne)]

I don’t think MNE has any functionalities to do temporal interpolation. You will have to do this manually:

• get the underlying numpy data array with raw.get_data().
• determine the corrupted segment’s start/stop indices.
• interpolate.

I do know researchers in Geneva doing EEG + TMS, and doing something very similar. They try to use high-sampling rate recordings (5 kHz) and interpolate the TMS-corrupted segments. Maybe this is a required feature for TMS studies? But I’m a bit concerned about the validity of this approach.

[agramfort (Alexandre Gramfort)]

you have this https://mne.tools/stable/generated/mne.preprocessing.fix_stim_artifact.html

but it’s not very generic

Alex

[oconnell_E]

Hi,

The data set we have (it’s on here: Dealing with TMS-EEG datasets - FieldTrip toolbox) consists of 2 conditions (marked as 1 and 3 in the data set). At the start of the each condition, there is a TMS pulse (marked as 10001) which causes huge artifacts in the data and we want to remove them -

There are several approaches we want to try to deal with that artifacts:

1. Replace with baseline (with EEG before the TMS pulse)
2. Interpolation
3. Replace with zero or 1 (with constant values)

After applying one of these approaches, we want to epoch the continuous EEG data based on the TMS events, therefore, epoching can also be done after removing and replacing the TMS artifacts.

We want to perform any kinds of preprocessing (resampling, filtering…etc.) only after dealing with the TMS artifacts, since it could result in additional noise in the data. For that reason we want to deal with the TMS artifacts in the first place.

To do that, we want to know first whether it’s feasible with MNE to remove the data segment of the TMS artifacts (by detecting the segment with the help of the TMS onset marker ‘10001’ and the time range before TMS and after TMS pulse we want to remove - depicted as tms_artefact_pre =0.005 and tms_artefact_post=0.015 in the first question above) and replace the values in this data segment with NaN.

Epoching and Interpolation are therefore not the scope of the current question!

Thanks a lot.

[richard (Richard Höchenberger)]

@oconnell_E

we want to know first whether it’s feasible with MNE to remove the data segment of the TMS artifacts (by detecting the segment with the help of the TMS onset marker ‘10001’ and the time range before TMS and after TMS pulse we want to remove - depicted as tms_artefact_pre =0.005 and tms_artefact_post=0.015 in the first question above) and replace the values in this data segment with NaN.

Thanks for the clarification, @oconnell_E!

Yes, this is possible. I would do it on the raw data level by first creating annotations for the affected segments, and then filling them with np.nan.

Let us know if you need any guidance there!

Best wishes,
Richard

[oconnell_E]

Thank you for your suggestion!
Actually, I’m having a problem with indexing on the raw data (mne.io.array.array.RawArray).
The information of the raw data is this (the data contains EEG data of 63 channels) :
Out[114]: <RawArray | 63 x 479060 (958.1 s), ~230.3 MB, data loaded>

What I’m trying to do is to replace values in the corresponding columns (affected segments) of the raw data with 0 by using the array below indicating the timing of TMS onset, and then fill them with np.nan.

TMS_pulse_index
Out[118]: 
array([   770,   4690,  52270,  98717, 145114, 191817, 239887, 288074,
       333588, 382353, 427734])

I tried with this following code, but it did’t work.

for column in TMS_pulse_index:
    for row in range(len(raw)):
        raw[row][column] = 0

I’d appreciate any feedback. Thanks!

[richard (Richard Höchenberger)]

Here’s a MWE, this should help get you started.

Best wishes,
Richard

# %%
import numpy as np
import mne

sample_dir = mne.datasets.sample.data_path()
sample_fname = sample_dir / 'MEG' / 'sample' / 'sample_audvis_raw.fif'

raw = mne.io.read_raw_fif(sample_fname)
raw.crop(tmax=60)

events = mne.find_events(raw, stim_channel='STI 014')

data = raw.get_data()

# %%
event_onset_samples = events[:, 0] - raw.first_samp  # don't forget first_samp!
for event_onset_sample in event_onset_samples:
    samp_start = event_onset_sample - 100
    samp_end = event_onset_sample + 100
    data[:, samp_start:samp_end] = np.nan

raw_with_nans = mne.io.RawArray(
    data=data,
    info=raw.info,
    first_samp=raw.first_samp,  # don't forget first_samp!
)
raw_with_nans.plot(events=events)

[oconnell_E]

Thank you! It also worked well with my data too.
I want to confirm one thing…:

samp_start = event_onset_sample - 100
samp_end = event_onset_sample + 100

In this part, ‘-100’ and ‘+100’ should be replaced with the time range which need to be rejected, am I right?
Assuming we want to reject the data segment between TMS onset - 1 ms and TMS onset + 6 ms, and the sampling frequency is 600Hz, then the time range won’t be integers anymore.

samp_start = event_onset_sample - 600 * 0.001
samp_end = event_onset_sample + 600 * 0.006

What should I do in this case to reject the data segment? Please let me know if I’m mistaken.
Thanks a lot.

[richard (Richard Höchenberger)]

Hello,

@oconnell_E

samp_start = event_onset_sample - 600 * 0.001
samp_end = event_onset_sample + 600 * 0.006

this is correct in principle, but samp_start and samp_end must be integers as they refer to sample numbers. So you should probably do

samp_start  = int(round(event_onset_sample - 600 * 0.001))

and the equivalent for samp_end.

Best wishes,
Richard