Resting state with dataset and example (#400)

* In some places, the virtual reality dataset code was wrong. * fix: PC data not downloading. fix: inversion 12 blocks of 5 repetitions * push example from Pedro * fix error with datframe initialization * [pre-commit.ci] auto fixes from pre-commit.com hooks * add whats new * add test * [pre-commit.ci] auto fixes from pre-commit.com hooks * fix pytest/unittest * [pre-commit.ci] auto fixes from pre-commit.com hooks * replace logging by warnings library * move docstring to the top * [pre-commit.ci] auto fixes from pre-commit.com hooks * test completed * [pre-commit.ci] auto fixes from pre-commit.com hooks * leftover * typo >< * Update examples/vr_pc_p300_different_epoch_size.py Co-authored-by: Sylvain Chevallier <sylvain.chevallier@universite-paris-saclay.fr> * rename into plot_vr_pc_p300_different_epoch_size.py * - Add figure plot - add comments * [pre-commit.ci] auto fixes from pre-commit.com hooks * Update plot_vr_pc_p300_different_epoch_size.py * [pre-commit.ci] auto fixes from pre-commit.com hooks * Create resting_state.py * push resting state * add dataset * push example * couple of bug fixes * add a condition to p300 to ignore Target/NonTarget check Fix loading of the mat file * working example * improve doc * [pre-commit.ci] auto fixes from pre-commit.com hooks * Update whats_new.rst * Update phmd_ml.py * Update plot_phmd_ml_spectrum.py flake8 * complete documentation * improve lisibility * push test * [pre-commit.ci] auto fixes from pre-commit.com hooks * fix tests * event_list missing in initialization. Correct typo. * [pre-commit.ci] auto fixes from pre-commit.com hooks * fix typo * [pre-commit.ci] auto fixes from pre-commit.com hooks * Applying and improving small details inside the tutorial * [pre-commit.ci] auto fixes from pre-commit.com hooks --------- Co-authored-by: Gregoire Cattan <gregoire.cattan@ibm.com> Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com> Co-authored-by: Bru <a.bruno@aluno.ufabc.edu.br> Co-authored-by: Sylvain Chevallier <sylvain.chevallier@universite-paris-saclay.fr>
NeuroTechX · Jun 22, 2023 · 2821954 · 2821954
1 parent 7559bca
commit 2821954
Show file tree

Hide file tree

Showing 10 changed files with 371 additions and 5 deletions.
diff --git a/docs/source/dataset_summary.rst b/docs/source/dataset_summary.rst
@@ -78,6 +78,19 @@ SSVEP
    Wang2016,34,62,40,6,5s,250Hz,1
 
 
+Resting States
+======================
+
+Include neuro experiments where the participant is not actively doing something.
+For example, recoding the EEG of a subject while s/he is having the eye closed or opened
+is a resting state experiment.
+
+.. csv-table::
+   :header: Dataset, #Subj, #Chan, #Classes, #Blocks / class, Trials length, Sampling rate, #Sessions
+   :class: sortable
+
+   HeadMountedDisplay,12,16,2,10,60s,512Hz,1
+
 
 Submit a new dataset
 ~~~~~~~~~~~~~~~~~~~~

diff --git a/docs/source/datasets.rst b/docs/source/datasets.rst
@@ -75,6 +75,17 @@ SSVEP Datasets
     Lee2019_SSVEP
 
 
+----------------------
+Resting State Datasets
+----------------------
+
+.. autosummary::
+    :toctree: generated/
+    :template: class.rst
+
+    HeadMountedDisplay
+
+
 ------------
 Base & Utils
 ------------

diff --git a/docs/source/whats_new.rst b/docs/source/whats_new.rst
@@ -25,6 +25,7 @@ Enhancements
 - Adding second deployment of the documentation (:gh:`374` by `Bruno Aristimunha`_)
 - Adding Parallel evaluation for :func:`moabb.evaluations.WithinSessionEvaluation` , :func:`moabb.evaluations.CrossSessionEvaluation` (:gh:`364` by `Bruno Aristimunha`_)
 - Add example with VirtualReality BrainInvaders dataset (:gh:`393` by `Gregoire Cattan`_ and `Pedro L. C. Rodrigues`_)
+- Add resting state paradigm with dataset and example (:gh:`400` by `Gregoire Cattan`_ and `Pedro L. C. Rodrigues`_)
 
 Bugs
 ~~~~

diff --git a/examples/plot_phmd_ml_spectrum.py b/examples/plot_phmd_ml_spectrum.py
@@ -0,0 +1,74 @@
+"""
+================================
+Spectral analysis of the trials
+================================
+
+This example demonstrates how to perform spectral
+analysis on epochs extracted from a specific subject
+within the :class:`moabb.datasets.HeadMountedDisplay`  dataset.
+
+"""
+
+# Authors: Pedro Rodrigues <pedro.rodrigues01@gmail.com>
+# Modified by: Gregoire Cattan <gcattan@hotmail.fr>
+# License: BSD (3-clause)
+
+import warnings
+
+import matplotlib.pyplot as plt
+import numpy as np
+from scipy.signal import welch
+
+from moabb.datasets import HeadMountedDisplay
+from moabb.paradigms import RestingStateToP300Adapter
+
+
+warnings.filterwarnings("ignore")
+
+###############################################################################
+# Initialization
+# ---------------
+#
+# 1) Specify the channel and subject to compute the power spectrum.
+# 2) Create an instance of the :class:`moabb.datasets.HeadMountedDisplay` dataset.
+# 3) Create an instance of the :class:`moabb.paradigms.RestingStateToP300Adapter`  paradigm.
+#    By default, the data is filtered between 1-35 Hz,
+#    and epochs are extracted from 10 to 50 seconds after event tagging.
+
+# Select channel and subject for the remaining of the example.
+channel = "Cz"
+subject = 1
+
+dataset = HeadMountedDisplay()
+events = ["on", "off"]
+paradigm = RestingStateToP300Adapter(events=events, channels=[channel])
+
+
+###############################################################################
+# Estimate Power Spectral Density
+# ---------------
+# 1) Obtain the epochs for the specified subject.
+# 2) Use Welch's method to estimate the power spectral density.
+
+X, y, _ = paradigm.get_data(dataset, [subject])
+f, S = welch(X, axis=-1, nperseg=1024, fs=paradigm.resample)
+
+###############################################################################
+# Display of the data
+# ---------------
+#
+# Plot the averaged Power Spectral Density (PSD) for each label condition,
+# using the selected channel specified at the beginning of the script.
+
+fig, ax = plt.subplots(facecolor="white", figsize=(8.2, 5.1))
+for condition in events:
+    mean_power = np.mean(S[y == condition], axis=0).flatten()
+    ax.plot(f, 10 * np.log10(mean_power), label=condition)
+
+ax.set_xlim(paradigm.fmin, paradigm.fmax)
+ax.set_ylim(100, 135)
+ax.set_ylabel("Spectrum Magnitude (dB)", fontsize=14)
+ax.set_xlabel("Frequency (Hz)", fontsize=14)
+ax.set_title("PSD for Channel " + channel, fontsize=16)
+ax.legend()
+fig.show()
diff --git a/moabb/datasets/__init__.py b/moabb/datasets/__init__.py
@@ -35,6 +35,7 @@
 from .Lee2019 import Lee2019_ERP, Lee2019_MI, Lee2019_SSVEP
 from .mpi_mi import MunichMI
 from .neiry import DemonsP300
+from .phmd_ml import HeadMountedDisplay
 from .physionet_mi import PhysionetMI
 from .schirrmeister2017 import Schirrmeister2017
 from .sosulski2019 import Sosulski2019

diff --git a/moabb/datasets/phmd_ml.py b/moabb/datasets/phmd_ml.py
@@ -0,0 +1,124 @@
+import os
+
+import mne
+import numpy as np
+from scipy.io import loadmat
+
+from . import download as dl
+from .base import BaseDataset
+
+
+HEADMOUNTED_URL = "https://zenodo.org/record/2617085/files/"
+
+
+class HeadMountedDisplay(BaseDataset):
+    """
+    Passive Head Mounted Display with Music Listening dataset.
+
+    .. admonition:: Dataset summary
+
+
+        =================  =======  =======  ==========  =================  ============  ===============  ===========
+        Name                 #Subj    #Chan    #Classes    #Blocks/class     Trials len    Sampling rate    #Sessions
+        ================== =======  =======  ==========  =================  ============  ===============  ===========
+        HeadMountedDisplay      12       16           2                 10   60s            512Hz                   1
+        ================== =======  =======  ==========  =================  ============  ===============  ===========
+
+    We describe the experimental procedures for a dataset that we have made publicly available
+    at https://doi.org/10.5281/zenodo.2617084 in mat (Mathworks, Natick, USA) and csv formats.
+    This dataset contains electroencephalographic recordings of 12 subjects listening to music
+    with and without a passive head-mounted display, that is, a head-mounted display which does
+    not include any electronics at the exception of a smartphone. The electroencephalographic
+    headset consisted of 16 electrodes. Data were recorded during a pilot experiment taking
+    place in the GIPSA-lab, Grenoble, France, in 2017 (Cattan and al, 2018).
+    The ID of this dataset is PHMDML.EEG.2017-GIPSA.
+
+    **full description of the experiment**
+    https://hal.archives-ouvertes.fr/hal-02085118
+
+    **Link to the data**
+    https://doi.org/10.5281/zenodo.2617084
+
+    **Authors**
+    Principal Investigator: Eng. Grégoire Cattan
+    Technical Supervisors: Eng. Pedro L. C. Rodrigues
+    Scientific Supervisor: Dr. Marco Congedo
+
+    **ID of the dataset**
+    PHMDML.EEG.2017-GIPSA
+
+    Notes
+    -----
+
+    .. versionadded:: 0.6.0
+
+    References
+    ----------
+
+    .. [1] G. Cattan, P. L. Coelho Rodrigues, and M. Congedo,
+        ‘Passive Head-Mounted Display Music-Listening EEG dataset’,
+        Gipsa-Lab ; IHMTEK, Research Report 2, Mar. 2019. doi: 10.5281/zenodo.2617084.
+
+
+    """
+
+    def __init__(self):
+        super().__init__(
+            subjects=list(range(1, 12 + 1)),
+            sessions_per_subject=1,
+            events=dict(on=1, off=2),
+            code="PHMD-ML",
+            interval=[0, 1],
+            paradigm="rstate",
+            doi="https://doi.org/10.5281/zenodo.2617084 ",
+        )
+        self._chnames = [
+            "Fp1",
+            "Fp2",
+            "Fc5",
+            "Fz",
+            "Fc6",
+            "T7",
+            "Cz",
+            "T8",
+            "P7",
+            "P3",
+            "Pz",
+            "P4",
+            "P8",
+            "O1",
+            "Oz",
+            "O2",
+            "stim",
+        ]
+        self._chtypes = ["eeg"] * 16 + ["stim"]
+
+    def _get_single_subject_data(self, subject):
+        """return data for a single subject"""
+
+        filepath = self.data_path(subject)[0]
+        data = loadmat(os.path.join(filepath, os.listdir(filepath)[0]))
+
+        first_channel = 1
+        last_channel = 17
+        S = data["data"][:, first_channel:last_channel]
+        stim = data["data"][:, -1]
+
+        X = np.concatenate([S, stim[:, None]], axis=1).T
+
+        info = mne.create_info(
+            ch_names=self._chnames, sfreq=512, ch_types=self._chtypes, verbose=False
+        )
+        raw = mne.io.RawArray(data=X, info=info, verbose=False)
+        return {"session_0": {"run_0": raw}}
+
+    def data_path(
+        self, subject, path=None, force_update=False, update_path=None, verbose=None
+    ):
+        if subject not in self.subject_list:
+            raise (ValueError("Invalid subject number"))
+
+        url = "{:s}subject_{:02d}.mat".format(HEADMOUNTED_URL, subject)
+        file_path = dl.data_path(url, "HEADMOUNTED")
+
+        return [file_path]
diff --git a/moabb/paradigms/__init__.py b/moabb/paradigms/__init__.py
@@ -9,4 +9,5 @@
 
 # flake8: noqa
 from moabb.paradigms.p300 import *
+from moabb.paradigms.resting_state import *
 from moabb.paradigms.ssvep import *
diff --git a/moabb/paradigms/p300.py b/moabb/paradigms/p300.py
@@ -168,11 +168,15 @@ def process_raw(  # noqa: C901
 
         # pick events, based on event_id
         try:
-            if type(event_id["Target"]) is list and type(event_id["NonTarget"]) == list:
-                event_id_new = dict(Target=1, NonTarget=0)
-                events = mne.merge_events(events, event_id["Target"], 1)
-                events = mne.merge_events(events, event_id["NonTarget"], 0)
-                event_id = event_id_new
+            if "Target" in event_id and "NonTarget" in event_id:
+                if (
+                    type(event_id["Target"]) is list
+                    and type(event_id["NonTarget"]) == list
+                ):
+                    event_id_new = dict(Target=1, NonTarget=0)
+                    events = mne.merge_events(events, event_id["Target"], 1)
+                    events = mne.merge_events(events, event_id["NonTarget"], 0)
+                    event_id = event_id_new
             events = mne.pick_events(events, include=list(event_id.values()))
         except RuntimeError:
             # skip raw if no event found
@@ -317,6 +321,14 @@ def __init__(self, fmin=1, fmax=24, **kwargs):
             raise (ValueError("P300 does not take argument filters"))
         super().__init__(filters=[[fmin, fmax]], **kwargs)
 
+    @property
+    def fmax(self):
+        return self.filters[0][1]
+
+    @property
+    def fmin(self):
+        return self.filters[0][0]
+
 
 class P300(SinglePass):
     """P300 for Target/NonTarget classification

diff --git a/moabb/paradigms/resting_state.py b/moabb/paradigms/resting_state.py
@@ -0,0 +1,81 @@
+"""Resting state Paradigms
+
+Regroups paradigms for experience where we record the EEG
+and the participant is not doing an active task, such
+as focusing, counting or speaking.
+
+Typically, a open/close eye experiment, where we
+record the EEG of a subject while he is having the eye open or close
+is a resting state experiment.
+"""
+
+from moabb.paradigms.p300 import SinglePass
+
+
+class RestingStateToP300Adapter(SinglePass):
+    """Adapter to the P300 paradigm for resting state experiments.
+    It implements a SinglePass processing as for P300, except that:
+    - the name of the event is free (it is not enforced to Target/NonTarget as for P300)
+    - the default values are different. In particular, the length of the epochs is larger.
+
+    Parameters
+    ----------
+    fmin: float (default 1)
+        cutoff frequency (Hz) for the high pass filter
+
+    fmax: float (default 35)
+        cutoff frequency (Hz) for the low pass filter
+
+    events: List of str | None (default None)
+        event to use for epoching. If None, default to all events defined in
+        the dataset.
+
+    tmin: float (default 10s)
+        Start time (in second) of the epoch, relative to the dataset specific
+        task interval e.g. tmin = 1 would mean the epoch will start 1 second
+        after the beginning of the task as defined by the dataset.
+
+    tmax: float | None, (default 50s)
+        End time (in second) of the epoch, relative to the beginning of the
+        dataset specific task interval. tmax = 5 would mean the epoch will end
+        5 second after the beginning of the task as defined in the dataset. If
+        None, use the dataset value.
+
+    resample: float | None (default 128)
+        If not None, resample the eeg data with the sampling rate provided.
+
+    baseline: None | tuple of length 2
+            The time interval to consider as “baseline” when applying baseline
+            correction. If None, do not apply baseline correction.
+            If a tuple (a, b), the interval is between a and b (in seconds),
+            including the endpoints.
+            Correction is applied by computing the mean of the baseline period
+            and subtracting it from the data (see mne.Epochs)
+
+    channels: list of str | None (default None)
+        list of channel to select. If None, use all EEG channels available in
+        the dataset.
+    """
+
+    def __init__(self, fmin=1, fmax=35, tmin=10, tmax=50, resample=128, **kwargs):
+        super().__init__(
+            fmin=fmin, fmax=fmax, tmin=tmin, tmax=tmax, resample=resample, **kwargs
+        )
+
+    def used_events(self, dataset):
+        return {ev: dataset.event_id[ev] for ev in self.events}
+
+    def is_valid(self, dataset):
+        ret = True
+        if not (dataset.paradigm == "rstate"):
+            ret = False
+
+        if self.events:
+            if not set(self.events) <= set(dataset.event_id.keys()):
+                ret = False
+
+        return ret
+
+    @property
+    def scoring(self):
+        return "roc_auc"