mne_decoding.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# This part is from https://mne.tools/stable/auto_tutorials/machine-learning/plot_sensors_decoding.html \n",
    "\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "from sklearn.pipeline import make_pipeline\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "\n",
    "import mne\n",
    "from mne.datasets import sample\n",
    "from mne.decoding import (SlidingEstimator, GeneralizingEstimator, Scaler,\n",
    "                          cross_val_multiscore, LinearModel, get_coef,\n",
    "                          Vectorizer, CSP)\n",
    "\n",
    "data_path = sample.data_path()\n",
    "\n",
    "subjects_dir = data_path + '/subjects'\n",
    "raw_fname = data_path + '/MEG/sample/sample_audvis_raw.fif'\n",
    "tmin, tmax = -0.200, 0.500\n",
    "event_id = {'Auditory/Left': 1, 'Visual/Left': 3}  # just use two\n",
    "raw = mne.io.read_raw_fif(raw_fname, preload=True)\n",
    "\n",
    "# The subsequent decoding analyses only capture evoked responses, so we can\n",
    "# low-pass the MEG data. Usually a value more like 40 Hz would be used,\n",
    "# but here low-pass at 20 so we can more heavily decimate, and allow\n",
    "# the examlpe to run faster. The 2 Hz high-pass helps improve CSP.\n",
    "raw.filter(2, 20)\n",
    "events = mne.find_events(raw, 'STI 014')\n",
    "\n",
    "# Set up pick list: EEG + MEG - bad channels (modify to your needs)\n",
    "raw.info['bads'] += ['MEG 2443', 'EEG 053']  # bads + 2 more\n",
    "\n",
    "# Read epochs\n",
    "epochs = mne.Epochs(raw, events, event_id, tmin, tmax, proj=True,\n",
    "                    picks=('mag', 'grad', 'eog'), baseline=(None, 0.), preload=True,\n",
    "                    reject=dict(grad=4000e-13, eog=150e-6), decim=10)\n",
    "epochs.pick_types(meg=True, exclude='bads')  # remove stim and EOG\n",
    "del raw\n",
    "\n",
    "X = epochs.get_data()  # MEG signals: n_epochs, n_meg_channels, n_times\n",
    "y = epochs.events[:, 2]  # target: Audio left or right"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "epochs.info"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y[y==3] = 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.model_selection import StratifiedShuffleSplit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cv = StratifiedShuffleSplit(n_splits=5, test_size=0.2, random_state=124)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "clf = make_pipeline(StandardScaler(), LogisticRegression())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "time_decod = SlidingEstimator(clf, scoring='roc_auc')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scores = cross_val_multiscore(time_decod, X, y, cv=cv)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scores.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(epochs.times, scores.mean(axis=0));"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "clf = make_pipeline(StandardScaler(),\n",
    "                    LinearModel(LogisticRegression()))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "time_decod = SlidingEstimator(clf, scoring='roc_auc')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "time_decod.fit(X, y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "coef = get_coef(time_decod, 'filters_', inverse_transform=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "coef.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evoked_coef = mne.EvokedArray(coef,\n",
    "                              epochs.info,\n",
    "                              tmin=epochs.times[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "evoked_coef.plot_joint(); "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.naive_bayes import GaussianNB"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "clf_nb = make_pipeline(StandardScaler(), GaussianNB())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "time_decod_nb = SlidingEstimator(clf_nb, scoring='roc_auc')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scores_nb = cross_val_multiscore(time_decod_nb, X, y, cv=cv)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(epochs.times, scores.mean(axis=0), 'b', label='LR')\n",
    "plt.plot(epochs.times, scores_nb.mean(axis=0), 'r', label='GNB');\n",
    "plt.legend();\n",
    "plt.grid();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scores.mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scores_nb.mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(epochs.times ,scores.std(axis=0), 'b', label='LR')\n",
    "plt.plot(epochs.times ,scores_nb.std(axis=0), 'r', label='GNB')\n",
    "plt.legend();\n",
    "plt.grid();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}