.. only:: html
.. note::
:class: sphx-glr-download-link-note
Click :ref:`here ` to download the full example code
.. rst-class:: sphx-glr-example-title
.. _sphx_glr_auto_examples_plot_mne_epoch_masking.py:
Interface with MNE-python
-------------------------
This example shows how to use this package with MNE-python.
It relies on the function ``raw_to_mask``, which takes as input a MNE.Raw
instance and an events array, and returns the corresponding input signals
and masks for the ``Comodulogram.fit`` method.
.. code-block:: default
import mne
import numpy as np
import os.path as op
import matplotlib.pyplot as plt
from pactools import simulate_pac, raw_to_mask, Comodulogram, MaskIterator
.. rst-class:: sphx-glr-script-out
Out:
.. code-block:: none
/home/tom/.conda/envs/py38/lib/python3.8/site-packages/mne/viz/backends/renderer.py:32: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__
_mod = importlib.__import__(
Simulate a dataset and save it out
.. code-block:: default
n_events = 100
mu = 1. # mean onset of PAC in seconds
sigma = 0.25 # standard deviation of onset of PAC in seconds
trial_len = 2. # len of the simulated trial in seconds
first_samp = 5 # seconds before the first sample and after the last
fs = 200. # Hz
high_fq = 50.0 # Hz
low_fq = 3.0 # Hz
low_fq_width = 2.0 # Hz
n_points = int(trial_len * fs)
noise_level = 0.4
def gaussian1d(array, mu, sigma):
return np.exp(-0.5 * ((array - mu) / sigma) ** 2)
# leave one channel for events and make signal as long as events
# with a bit of room on either side so events don't get cut off
signal = np.zeros((2, int(n_points * n_events + 2 * first_samp * fs)))
events = np.zeros((n_events, 3), dtype=int)
events[:, 0] = np.arange((first_samp + mu) * fs,
n_points * n_events + (first_samp + mu) * fs,
trial_len * fs, dtype=int)
events[:, 2] = np.ones((n_events))
mod_fun = gaussian1d(np.arange(n_points), mu * fs, sigma * fs)
for i in range(n_events):
signal_no_pac = simulate_pac(n_points=n_points, fs=fs,
high_fq=high_fq, low_fq=low_fq,
low_fq_width=low_fq_width,
noise_level=1.0, random_state=i)
signal_pac = simulate_pac(n_points=n_points, fs=fs,
high_fq=high_fq, low_fq=low_fq,
low_fq_width=low_fq_width,
noise_level=noise_level, random_state=i)
signal[0, i * n_points:(i + 1) * n_points] = \
signal_pac * mod_fun + signal_no_pac * (1 - mod_fun)
info = mne.create_info(['Ch1', 'STI 014'], fs, ['eeg', 'stim'])
raw = mne.io.RawArray(signal, info)
raw.add_events(events, stim_channel='STI 014')
.. rst-class:: sphx-glr-script-out
Out:
.. code-block:: none
Creating RawArray with float64 data, n_channels=2, n_times=42000
Range : 0 ... 41999 = 0.000 ... 209.995 secs
Ready.
Let's plot the signal and its power spectral density to visualize the data.
As shown in the plots below, there is a peak for the driver frequency at
3 Hz and a peak for the carrier frequency at 50 Hz but phase-amplitude
coupling cannot be seen in the evoked plot by eye because the signal is
averaged over different phases for each epoch.
.. code-block:: default
raw.plot_psd(fmax=60)
epochs = mne.Epochs(raw, events, tmin=-3, tmax=3)
epochs.average().plot()
.. rst-class:: sphx-glr-horizontal
*
.. image:: /auto_examples/images/sphx_glr_plot_mne_epoch_masking_001.png
:class: sphx-glr-multi-img
*
.. image:: /auto_examples/images/sphx_glr_plot_mne_epoch_masking_002.png
:class: sphx-glr-multi-img
.. rst-class:: sphx-glr-script-out
Out:
.. code-block:: none
Effective window size : 10.240 (s)
Need more than one channel to make topography for eeg. Disabling interactivity.
/home/tom/work/github/pactools/examples/plot_mne_epoch_masking.py:71: RuntimeWarning: Channel locations not available. Disabling spatial colors.
raw.plot_psd(fmax=60)
100 matching events found
Applying baseline correction (mode: mean)
Not setting metadata
0 projection items activated
Need more than one channel to make topography for eeg. Disabling interactivity.