Note
Go to the end to download the full example code.
Decoding (MVPA)#
Design philosophy#
Decoding (a.k.a. MVPA) in MNE largely follows the machine learning API of the
scikit-learn package.
Each estimator implements fit, transform, fit_transform, and
(optionally) inverse_transform methods. For more details on this design,
visit scikit-learn. For additional theoretical insights into the decoding
framework in MNE [1].
For ease of comprehension, we will denote instantiations of the class using the same name as the class but in small caps instead of camel cases.
Let’s start by loading data for a simple two-class problem:
importmatplotlib.pyplotasplt importnumpyasnp fromsklearn.linear_modelimport LogisticRegression fromsklearn.pipelineimport make_pipeline fromsklearn.preprocessingimport StandardScaler importmne frommne.datasetsimport sample frommne.decodingimport ( CSP , GeneralizingEstimator , LinearModel , Scaler , SlidingEstimator , Vectorizer , cross_val_multiscore , get_coef , get_spatial_filter_from_estimator , ) data_path = sample.data_path () subjects_dir = data_path / "subjects" meg_path = data_path / "MEG" / "sample" raw_fname = meg_path / "sample_audvis_filt-0-40_raw.fif" tmin , tmax = -0.200, 0.500 event_id = {"Auditory/Left": 1, "Visual/Left": 3} # just use two raw = mne.io.read_raw_fif (raw_fname ) raw.pick(picks=["grad", "stim", "eog"]) # The subsequent decoding analyses only capture evoked responses, so we can # low-pass the MEG data. Usually a value more like 40 Hz would be used, # but here low-pass at 20 so we can more heavily decimate, and allow # the example to run faster. The 2 Hz high-pass helps improve CSP. raw.load_data().filter(2, 20) events = mne.find_events (raw, "STI 014") # Set up bad channels (modify to your needs) raw.info["bads"] += ["MEG 2443"] # bads + 2 more # Read epochs epochs = mne.Epochs ( raw, events , event_id , tmin , tmax , proj=True, picks=("grad", "eog"), baseline=(None, 0.0), preload=True, reject=dict(grad=4000e-13, eog=150e-6), decim=3, verbose="error", ) epochs.pick(picks="meg", exclude="bads") # remove stim and EOG del raw X = epochs.get_data(copy=False) # MEG signals: n_epochs, n_meg_channels, n_times y = epochs.events [:, 2] # target: auditory left vs visual left
Opening raw data file /home/circleci/mne_data/MNE-sample-data/MEG/sample/sample_audvis_filt-0-40_raw.fif... Read a total of 4 projection items: PCA-v1 (1 x 102) idle PCA-v2 (1 x 102) idle PCA-v3 (1 x 102) idle Average EEG reference (1 x 60) idle Range : 6450 ... 48149 = 42.956 ... 320.665 secs Ready. Reading 0 ... 41699 = 0.000 ... 277.709 secs... Filtering raw data in 1 contiguous segment Setting up band-pass filter from 2 - 20 Hz FIR filter parameters --------------------- Designing a one-pass, zero-phase, non-causal bandpass filter: - Windowed time-domain design (firwin) method - Hamming window with 0.0194 passband ripple and 53 dB stopband attenuation - Lower passband edge: 2.00 - Lower transition bandwidth: 2.00 Hz (-6 dB cutoff frequency: 1.00 Hz) - Upper passband edge: 20.00 Hz - Upper transition bandwidth: 5.00 Hz (-6 dB cutoff frequency: 22.50 Hz) - Filter length: 249 samples (1.658 s) Finding events on: STI 014 319 events found on stim channel STI 014 Event IDs: [ 1 2 3 4 5 32]
Transformation classes#
Scaler#
The mne.decoding.Scaler will standardize the data based on channel
scales. In the simplest modes scalings=None or scalings=dict(...),
each data channel type (e.g., mag, grad, eeg) is treated separately and
scaled by a constant. This is the approach used by e.g.,
mne.compute_covariance() to standardize channel scales.
If scalings='mean' or scalings='median', each channel is scaled using
empirical measures. Each channel is scaled independently by the mean and
standand deviation, or median and interquartile range, respectively, across
all epochs and time points during fit
(during training). The transform() method is
called to transform data (training or test set) by scaling all time points
and epochs on a channel-by-channel basis. To perform both the fit and
transform operations in a single call, the
fit_transform() method may be used. To invert the
transform, inverse_transform() can be used. For
scalings='median', scikit-learn version 0.17+ is required.
Note
Using this class is different from directly applying
sklearn.preprocessing.StandardScaler or
sklearn.preprocessing.RobustScaler offered by
scikit-learn. These scale each classification feature, e.g.
each time point for each channel, with mean and standard
deviation computed across epochs, whereas
mne.decoding.Scaler scales each channel using mean and
standard deviation computed across all of its time points
and epochs.
Vectorizer#
Scikit-learn API provides functionality to chain transformers and estimators
by using sklearn.pipeline.Pipeline. We can construct decoding
pipelines and perform cross-validation and grid-search. However scikit-learn
transformers and estimators generally expect 2D data
(n_samples * n_features), whereas MNE transformers typically output data
with a higher dimensionality
(e.g. n_samples * n_channels * n_frequencies * n_times). A Vectorizer
therefore needs to be applied between the MNE and the scikit-learn steps
like:
# Uses all MEG sensors and time points as separate classification # features, so the resulting filters used are spatio-temporal clf = make_pipeline ( Scaler (epochs.info), Vectorizer (), LogisticRegression (solver="liblinear"), # liblinear is faster than lbfgs ) scores = cross_val_multiscore (clf , X , y , cv=5, n_jobs=None) # Mean scores across cross-validation splits score = np.mean (scores , axis=0) print(f"Spatio-temporal: {100*score :0.1f}%")
Spatio-temporal: 99.2%
PSDEstimator#
The mne.decoding.PSDEstimator
computes the power spectral density (PSD) using the multitaper
method. It takes a 3D array as input, converts it into 2D and computes the
PSD.
FilterEstimator#
The mne.decoding.FilterEstimator filters the 3D epochs data.
Spatial filters#
Just like temporal filters, spatial filters provide weights to modify the data along the sensor dimension. They are popular in the BCI community because of their simplicity and ability to distinguish spatially-separated neural activity.
Common spatial pattern#
mne.decoding.CSP is a technique to analyze multichannel data based
on recordings from two classes [2] (see also
https://en.wikipedia.org/wiki/Common_spatial_pattern).
Let \(X \in R^{C\times T}\) be a segment of data with \(C\) channels and \(T\) time points. The data at a single time point is denoted by \(x(t)\) such that \(X=[x(t), x(t+1), ..., x(t+T-1)]\). Common spatial pattern (CSP) finds a decomposition that projects the signal in the original sensor space to CSP space using the following transformation:
where each column of \(W \in R^{C\times C}\) is a spatial filter and each row of \(x_{CSP}\) is a CSP component. The matrix \(W\) is also called the de-mixing matrix in other contexts. Let \(\Sigma^{+} \in R^{C\times C}\) and \(\Sigma^{-} \in R^{C\times C}\) be the estimates of the covariance matrices of the two conditions. CSP analysis is given by the simultaneous diagonalization of the two covariance matrices
where \(\lambda^{C}\) is a diagonal matrix whose entries are the eigenvalues of the following generalized eigenvalue problem
Large entries in the diagonal matrix corresponds to a spatial filter which gives high variance in one class but low variance in the other. Thus, the filter facilitates discrimination between the two classes.
Note
The winning entry of the Grasp-and-lift EEG competition in Kaggle used
the CSP implementation in MNE and was featured as
a script of the week.
We can use CSP with these data with:
csp = CSP (n_components=3, norm_trace=False) clf_csp = make_pipeline (csp , LinearModel (LogisticRegression (solver="liblinear"))) scores = cross_val_multiscore (clf_csp , X , y , cv=5, n_jobs=None) print(f"CSP: {100*scores .mean():0.1f}%")
Computing rank from data with rank=None Using tolerance 6e-11 (2.2e-16 eps * 203 dim * 1.3e+03 max singular value) Estimated rank (data): 203 data: rank 203 computed from 203 data channels with 0 projectors Reducing data rank from 203 -> 203 Estimating class=1 covariance using EMPIRICAL Done. Estimating class=3 covariance using EMPIRICAL Done. Computing rank from data with rank=None Using tolerance 6e-11 (2.2e-16 eps * 203 dim * 1.3e+03 max singular value) Estimated rank (data): 203 data: rank 203 computed from 203 data channels with 0 projectors Reducing data rank from 203 -> 203 Estimating class=1 covariance using EMPIRICAL Done. Estimating class=3 covariance using EMPIRICAL Done. Computing rank from data with rank=None Using tolerance 6e-11 (2.2e-16 eps * 203 dim * 1.3e+03 max singular value) Estimated rank (data): 203 data: rank 203 computed from 203 data channels with 0 projectors Reducing data rank from 203 -> 203 Estimating class=1 covariance using EMPIRICAL Done. Estimating class=3 covariance using EMPIRICAL Done. Computing rank from data with rank=None Using tolerance 5.8e-11 (2.2e-16 eps * 203 dim * 1.3e+03 max singular value) Estimated rank (data): 203 data: rank 203 computed from 203 data channels with 0 projectors Reducing data rank from 203 -> 203 Estimating class=1 covariance using EMPIRICAL Done. Estimating class=3 covariance using EMPIRICAL Done. Computing rank from data with rank=None Using tolerance 5.9e-11 (2.2e-16 eps * 203 dim * 1.3e+03 max singular value) Estimated rank (data): 203 data: rank 203 computed from 203 data channels with 0 projectors Reducing data rank from 203 -> 203 Estimating class=1 covariance using EMPIRICAL Done. Estimating class=3 covariance using EMPIRICAL Done. CSP: 89.4%
Source power comodulation (SPoC)#
Source Power Comodulation (mne.decoding.SPoC)
[3] identifies the composition of
orthogonal spatial filters that maximally correlate with a continuous target.
SPoC can be seen as an extension of the CSP where the target is driven by a continuous variable rather than a discrete variable. Typical applications include extraction of motor patterns using EMG power or audio patterns using sound envelope.
xDAWN#
mne.preprocessing.Xdawn is a spatial filtering method designed to
improve the signal to signal + noise ratio (SSNR) of the ERP responses
[4]. Xdawn was originally
designed for P300 evoked potential by enhancing the target response with
respect to the non-target response. The implementation in MNE-Python is a
generalization to any type of ERP.
Effect-matched spatial filtering#
The result of mne.decoding.EMS is a spatial filter at each time
point and a corresponding time course [5].
Intuitively, the result gives the similarity between the filter at
each time point and the data vector (sensors) at that time point.
Patterns vs. filters#
When interpreting the components of the CSP (or spatial filters in general), it is often more intuitive to think about how \(x(t)\) is composed of the different CSP components \(x_{CSP}(t)\). In other words, we can rewrite Equation (1) as follows:
The columns of the matrix \((W^{-1})^T\) are called spatial patterns. This is also called the mixing matrix. The example Linear classifier on sensor data with plot patterns and filters discusses the difference between patterns and filters.
These can be plotted for every spatial filter including CSP, XdawnTransformer, SSD and SPoC:
# Fit CSP on full data, plot eigenvalues sorted based on mutual information, # and plot patterns and filters for the three components largest components. csp.fit (X , y ) spf = get_spatial_filter_from_estimator (csp , info=epochs.info) spf.plot_scree () spf.plot_patterns (components=[0, 1, 2]) spf.plot_filters (components=[0, 1, 2], scalings=1e-9)
- Scree plot
- Pattern0, Pattern1, Pattern2, AU
- Filter0, Filter1, Filter2, AU
Computing rank from data with rank=None Using tolerance 6.6e-11 (2.2e-16 eps * 203 dim * 1.5e+03 max singular value) Estimated rank (data): 203 data: rank 203 computed from 203 data channels with 0 projectors Reducing data rank from 203 -> 203 Estimating class=1 covariance using EMPIRICAL Done. Estimating class=3 covariance using EMPIRICAL Done.
Decoding over time#
This strategy consists in fitting a multivariate predictive model on each
time instant and evaluating its performance at the same instant on new
epochs. The mne.decoding.SlidingEstimator will take as input a
pair of features \(X\) and targets \(y\), where \(X\) has
more than 2 dimensions. For decoding over time the data \(X\)
is the epochs data of shape n_epochs ×ばつ n_channels ×ばつ n_times. As the
last dimension of \(X\) is the time, an estimator will be fit
on every time instant.
This approach is analogous to SlidingEstimator-based approaches in fMRI, where here we are interested in when one can discriminate experimental conditions and therefore figure out when the effect of interest happens.
When working with linear models as estimators, this approach boils down to estimating a discriminative spatial filter for each time instant.
Temporal decoding#
We’ll use a Logistic Regression for a binary classification as machine learning model.
# We will train the classifier on all left visual vs auditory trials on MEG clf = make_pipeline (StandardScaler (), LogisticRegression (solver="liblinear")) time_decod = SlidingEstimator (clf , n_jobs=None, scoring="roc_auc", verbose=True) # here we use cv=3 just for speed scores = cross_val_multiscore (time_decod , X , y , cv=3, n_jobs=None) # Mean scores across cross-validation splits scores = np.mean (scores , axis=0) # Plot fig , ax = plt.subplots () ax.plot (epochs.times, scores , label="score") ax.axhline (0.5, color="k", linestyle="--", label="chance") ax.set_xlabel ("Times") ax.set_ylabel ("AUC") # Area Under the Curve ax.legend () ax.axvline (0.0, color="k", linestyle="-") ax.set_title ("Sensor space decoding")
0%| | Fitting SlidingEstimator : 0/36 [00:00<?, ?it/s] 22%|██▏ | Fitting SlidingEstimator : 8/36 [00:00<00:00, 236.04it/s] 47%|████▋ | Fitting SlidingEstimator : 17/36 [00:00<00:00, 251.92it/s] 72%|███████▏ | Fitting SlidingEstimator : 26/36 [00:00<00:00, 257.43it/s] 94%|█████████▍| Fitting SlidingEstimator : 34/36 [00:00<00:00, 252.24it/s] 100%|██████████| Fitting SlidingEstimator : 36/36 [00:00<00:00, 257.06it/s] 0%| | Fitting SlidingEstimator : 0/36 [00:00<?, ?it/s] 22%|██▏ | Fitting SlidingEstimator : 8/36 [00:00<00:00, 235.91it/s] 47%|████▋ | Fitting SlidingEstimator : 17/36 [00:00<00:00, 249.75it/s] 72%|███████▏ | Fitting SlidingEstimator : 26/36 [00:00<00:00, 255.97it/s] 97%|█████████▋| Fitting SlidingEstimator : 35/36 [00:00<00:00, 259.06it/s] 100%|██████████| Fitting SlidingEstimator : 36/36 [00:00<00:00, 262.76it/s] 0%| | Fitting SlidingEstimator : 0/36 [00:00<?, ?it/s] 22%|██▏ | Fitting SlidingEstimator : 8/36 [00:00<00:00, 235.44it/s] 47%|████▋ | Fitting SlidingEstimator : 17/36 [00:00<00:00, 251.87it/s] 75%|███████▌ | Fitting SlidingEstimator : 27/36 [00:00<00:00, 266.47it/s] 97%|█████████▋| Fitting SlidingEstimator : 35/36 [00:00<00:00, 258.88it/s] 100%|██████████| Fitting SlidingEstimator : 36/36 [00:00<00:00, 262.45it/s]
You can retrieve the spatial filters and spatial patterns if you explicitly use a LinearModel
clf = make_pipeline ( StandardScaler (), LinearModel (LogisticRegression (solver="liblinear")) ) time_decod = SlidingEstimator (clf , n_jobs=None, scoring="roc_auc", verbose=True) time_decod.fit (X , y ) coef = get_coef (time_decod , "patterns_", inverse_transform=True) evoked_time_gen = mne.EvokedArray (coef , epochs.info, tmin =epochs.times[0]) joint_kwargs = dict(ts_args=dict(time_unit="s"), topomap_args=dict(time_unit="s")) evoked_time_gen.plot_joint ( times=np.arange (0.0, 0.500, 0.100), title="patterns", **joint_kwargs )
0%| | Fitting SlidingEstimator : 0/36 [00:00<?, ?it/s] 11%|█ | Fitting SlidingEstimator : 4/36 [00:00<00:00, 117.88it/s] 22%|██▏ | Fitting SlidingEstimator : 8/36 [00:00<00:00, 118.42it/s] 36%|███▌ | Fitting SlidingEstimator : 13/36 [00:00<00:00, 128.99it/s] 53%|█████▎ | Fitting SlidingEstimator : 19/36 [00:00<00:00, 142.28it/s] 69%|██████▉ | Fitting SlidingEstimator : 25/36 [00:00<00:00, 150.08it/s] 83%|████████▎ | Fitting SlidingEstimator : 30/36 [00:00<00:00, 149.76it/s] 94%|█████████▍| Fitting SlidingEstimator : 34/36 [00:00<00:00, 144.64it/s] 100%|██████████| Fitting SlidingEstimator : 36/36 [00:00<00:00, 147.67it/s] No projector specified for this dataset. Please consider the method self.add_proj.
Temporal generalization#
Temporal generalization is an extension of the decoding over time approach. It consists in evaluating whether the model estimated at a particular time instant accurately predicts any other time instant. It is analogous to transferring a trained model to a distinct learning problem, where the problems correspond to decoding the patterns of brain activity recorded at distinct time instants.
The object to for Temporal generalization is
mne.decoding.GeneralizingEstimator. It expects as input \(X\)
and \(y\) (similarly to SlidingEstimator) but
generates predictions from each model for all time instants. The class
GeneralizingEstimator is generic and will treat the
last dimension as the one to be used for generalization testing. For
convenience, here, we refer to it as different tasks. If \(X\)
corresponds to epochs data then the last dimension is time.
This runs the analysis used in [6] and further detailed in [7]:
# define the Temporal generalization object time_gen = GeneralizingEstimator (clf , n_jobs=None, scoring="roc_auc", verbose=True) # again, cv=3 just for speed scores = cross_val_multiscore (time_gen , X , y , cv=3, n_jobs=None) # Mean scores across cross-validation splits scores = np.mean (scores , axis=0) # Plot the diagonal (it's exactly the same as the time-by-time decoding above) fig , ax = plt.subplots () ax.plot (epochs.times, np.diag (scores ), label="score") ax.axhline (0.5, color="k", linestyle="--", label="chance") ax.set_xlabel ("Times") ax.set_ylabel ("AUC") ax.legend () ax.axvline (0.0, color="k", linestyle="-") ax.set_title ("Decoding MEG sensors over time")
0%| | Fitting GeneralizingEstimator : 0/36 [00:00<?, ?it/s] 17%|█▋ | Fitting GeneralizingEstimator : 6/36 [00:00<00:00, 177.76it/s] 39%|███▉ | Fitting GeneralizingEstimator : 14/36 [00:00<00:00, 208.01it/s] 61%|██████ | Fitting GeneralizingEstimator : 22/36 [00:00<00:00, 218.13it/s] 81%|████████ | Fitting GeneralizingEstimator : 29/36 [00:00<00:00, 215.46it/s] 100%|██████████| Fitting GeneralizingEstimator : 36/36 [00:00<00:00, 220.89it/s] 100%|██████████| Fitting GeneralizingEstimator : 36/36 [00:00<00:00, 219.64it/s] 0%| | Scoring GeneralizingEstimator : 0/1296 [00:00<?, ?it/s] 2%|▏ | Scoring GeneralizingEstimator : 22/1296 [00:00<00:01, 649.27it/s] 3%|▎ | Scoring GeneralizingEstimator : 45/1296 [00:00<00:01, 663.43it/s] 5%|▌ | Scoring GeneralizingEstimator : 67/1296 [00:00<00:01, 657.13it/s] 7%|▋ | Scoring GeneralizingEstimator : 89/1296 [00:00<00:01, 655.48it/s] 9%|▊ | Scoring GeneralizingEstimator : 113/1296 [00:00<00:01, 661.16it/s] 10%|█ | Scoring GeneralizingEstimator : 136/1296 [00:00<00:01, 660.72it/s] 12%|█▏ | Scoring GeneralizingEstimator : 159/1296 [00:00<00:01, 661.75it/s] 14%|█▍ | Scoring GeneralizingEstimator : 182/1296 [00:00<00:01, 663.14it/s] 16%|█▌ | Scoring GeneralizingEstimator : 204/1296 [00:00<00:01, 661.99it/s] 18%|█▊ | Scoring GeneralizingEstimator : 227/1296 [00:00<00:01, 664.59it/s] 19%|█▉ | Scoring GeneralizingEstimator : 250/1296 [00:00<00:01, 665.16it/s] 21%|██ | Scoring GeneralizingEstimator : 272/1296 [00:00<00:01, 662.70it/s] 23%|██▎ | Scoring GeneralizingEstimator : 295/1296 [00:00<00:01, 661.78it/s] 25%|██▍ | Scoring GeneralizingEstimator : 319/1296 [00:00<00:01, 664.86it/s] 26%|██▋ | Scoring GeneralizingEstimator : 342/1296 [00:00<00:01, 665.18it/s] 28%|██▊ | Scoring GeneralizingEstimator : 365/1296 [00:00<00:01, 665.16it/s] 30%|██▉ | Scoring GeneralizingEstimator : 387/1296 [00:00<00:01, 663.90it/s] 32%|███▏ | Scoring GeneralizingEstimator : 410/1296 [00:00<00:01, 664.71it/s] 33%|███▎ | Scoring GeneralizingEstimator : 433/1296 [00:00<00:01, 665.16it/s] 35%|███▌ | Scoring GeneralizingEstimator : 456/1296 [00:00<00:01, 665.02it/s] 37%|███▋ | Scoring GeneralizingEstimator : 478/1296 [00:00<00:01, 663.98it/s] 39%|███▊ | Scoring GeneralizingEstimator : 502/1296 [00:00<00:01, 665.39it/s] 41%|████ | Scoring GeneralizingEstimator : 525/1296 [00:00<00:01, 664.97it/s] 42%|████▏ | Scoring GeneralizingEstimator : 547/1296 [00:00<00:01, 664.27it/s] 44%|████▍ | Scoring GeneralizingEstimator : 570/1296 [00:00<00:01, 665.10it/s] 46%|████▌ | Scoring GeneralizingEstimator : 593/1296 [00:00<00:01, 664.40it/s] 48%|████▊ | Scoring GeneralizingEstimator : 616/1296 [00:00<00:01, 664.19it/s] 49%|████▉ | Scoring GeneralizingEstimator : 639/1296 [00:00<00:00, 664.59it/s] 51%|█████ | Scoring GeneralizingEstimator : 662/1296 [00:00<00:00, 664.23it/s] 53%|█████▎ | Scoring GeneralizingEstimator : 685/1296 [00:01<00:00, 664.37it/s] 55%|█████▍ | Scoring GeneralizingEstimator : 707/1296 [00:01<00:00, 663.08it/s] 56%|█████▋ | Scoring GeneralizingEstimator : 730/1296 [00:01<00:00, 664.07it/s] 58%|█████▊ | Scoring GeneralizingEstimator : 753/1296 [00:01<00:00, 663.66it/s] 60%|█████▉ | Scoring GeneralizingEstimator : 775/1296 [00:01<00:00, 661.48it/s] 62%|██████▏ | Scoring GeneralizingEstimator : 799/1296 [00:01<00:00, 663.07it/s] 63%|██████▎ | Scoring GeneralizingEstimator : 821/1296 [00:01<00:00, 662.57it/s] 65%|██████▌ | Scoring GeneralizingEstimator : 844/1296 [00:01<00:00, 663.00it/s] 67%|██████▋ | Scoring GeneralizingEstimator : 867/1296 [00:01<00:00, 662.56it/s] 68%|██████▊ | Scoring GeneralizingEstimator : 886/1296 [00:01<00:00, 656.03it/s] 70%|███████ | Scoring GeneralizingEstimator : 909/1296 [00:01<00:00, 656.41it/s] 72%|███████▏ | Scoring GeneralizingEstimator : 931/1296 [00:01<00:00, 655.70it/s] 74%|███████▎ | Scoring GeneralizingEstimator : 954/1296 [00:01<00:00, 656.80it/s] 75%|███████▌ | Scoring GeneralizingEstimator : 977/1296 [00:01<00:00, 657.10it/s] 77%|███████▋ | Scoring GeneralizingEstimator : 999/1296 [00:01<00:00, 656.50it/s] 79%|███████▉ | Scoring GeneralizingEstimator : 1022/1296 [00:01<00:00, 657.47it/s] 81%|████████ | Scoring GeneralizingEstimator : 1044/1296 [00:01<00:00, 656.67it/s] 82%|████████▏ | Scoring GeneralizingEstimator : 1067/1296 [00:01<00:00, 657.81it/s] 84%|████████▍ | Scoring GeneralizingEstimator : 1090/1296 [00:01<00:00, 657.98it/s] 86%|████████▌ | Scoring GeneralizingEstimator : 1113/1296 [00:01<00:00, 657.90it/s] 88%|████████▊ | Scoring GeneralizingEstimator : 1136/1296 [00:01<00:00, 657.84it/s] 89%|████████▉ | Scoring GeneralizingEstimator : 1158/1296 [00:01<00:00, 657.24it/s] 91%|█████████ | Scoring GeneralizingEstimator : 1181/1296 [00:01<00:00, 657.74it/s] 93%|█████████▎| Scoring GeneralizingEstimator : 1204/1296 [00:01<00:00, 657.98it/s] 95%|█████████▍| Scoring GeneralizingEstimator : 1226/1296 [00:01<00:00, 656.99it/s] 96%|█████████▋| Scoring GeneralizingEstimator : 1249/1296 [00:01<00:00, 657.88it/s] 98%|█████████▊| Scoring GeneralizingEstimator : 1271/1296 [00:01<00:00, 657.48it/s] 100%|█████████▉| Scoring GeneralizingEstimator : 1294/1296 [00:01<00:00, 658.15it/s] 100%|██████████| Scoring GeneralizingEstimator : 1296/1296 [00:01<00:00, 659.97it/s] 0%| | Fitting GeneralizingEstimator : 0/36 [00:00<?, ?it/s] 17%|█▋ | Fitting GeneralizingEstimator : 6/36 [00:00<00:00, 177.74it/s] 36%|███▌ | Fitting GeneralizingEstimator : 13/36 [00:00<00:00, 193.01it/s] 58%|█████▊ | Fitting GeneralizingEstimator : 21/36 [00:00<00:00, 208.69it/s] 78%|███████▊ | Fitting GeneralizingEstimator : 28/36 [00:00<00:00, 207.68it/s] 94%|█████████▍| Fitting GeneralizingEstimator : 34/36 [00:00<00:00, 201.23it/s] 100%|██████████| Fitting GeneralizingEstimator : 36/36 [00:00<00:00, 205.91it/s] 0%| | Scoring GeneralizingEstimator : 0/1296 [00:00<?, ?it/s] 2%|▏ | Scoring GeneralizingEstimator : 22/1296 [00:00<00:01, 638.03it/s] 3%|▎ | Scoring GeneralizingEstimator : 44/1296 [00:00<00:01, 645.59it/s] 5%|▌ | Scoring GeneralizingEstimator : 67/1296 [00:00<00:01, 651.49it/s] 7%|▋ | Scoring GeneralizingEstimator : 90/1296 [00:00<00:01, 652.80it/s] 9%|▊ | Scoring GeneralizingEstimator : 113/1296 [00:00<00:01, 654.19it/s] 10%|█ | Scoring GeneralizingEstimator : 136/1296 [00:00<00:01, 656.22it/s] 12%|█▏ | Scoring GeneralizingEstimator : 158/1296 [00:00<00:01, 650.65it/s] 14%|█▍ | Scoring GeneralizingEstimator : 182/1296 [00:00<00:01, 656.40it/s] 16%|█▌ | Scoring GeneralizingEstimator : 205/1296 [00:00<00:01, 656.99it/s] 18%|█▊ | Scoring GeneralizingEstimator : 227/1296 [00:00<00:01, 656.67it/s] 19%|█▉ | Scoring GeneralizingEstimator : 250/1296 [00:00<00:01, 658.70it/s] 21%|██ | Scoring GeneralizingEstimator : 273/1296 [00:00<00:01, 658.66it/s] 23%|██▎ | Scoring GeneralizingEstimator : 295/1296 [00:00<00:01, 658.17it/s] 25%|██▍ | Scoring GeneralizingEstimator : 318/1296 [00:00<00:01, 658.28it/s] 26%|██▋ | Scoring GeneralizingEstimator : 341/1296 [00:00<00:01, 658.16it/s] 28%|██▊ | Scoring GeneralizingEstimator : 363/1296 [00:00<00:01, 655.39it/s] 30%|██▉ | Scoring GeneralizingEstimator : 387/1296 [00:00<00:01, 658.35it/s] 32%|███▏ | Scoring GeneralizingEstimator : 410/1296 [00:00<00:01, 658.84it/s] 33%|███▎ | Scoring GeneralizingEstimator : 432/1296 [00:00<00:01, 658.49it/s] 35%|███▌ | Scoring GeneralizingEstimator : 455/1296 [00:00<00:01, 658.96it/s] 37%|███▋ | Scoring GeneralizingEstimator : 478/1296 [00:00<00:01, 659.12it/s] 39%|███▊ | Scoring GeneralizingEstimator : 501/1296 [00:00<00:01, 659.52it/s] 40%|████ | Scoring GeneralizingEstimator : 524/1296 [00:00<00:01, 659.60it/s] 42%|████▏ | Scoring GeneralizingEstimator : 547/1296 [00:00<00:01, 659.14it/s] 44%|████▍ | Scoring GeneralizingEstimator : 570/1296 [00:00<00:01, 658.80it/s] 46%|████▌ | Scoring GeneralizingEstimator : 593/1296 [00:00<00:01, 658.97it/s] 48%|████▊ | Scoring GeneralizingEstimator : 616/1296 [00:00<00:01, 658.84it/s] 49%|████▉ | Scoring GeneralizingEstimator : 639/1296 [00:00<00:00, 658.71it/s] 51%|█████ | Scoring GeneralizingEstimator : 662/1296 [00:01<00:00, 658.73it/s] 53%|█████▎ | Scoring GeneralizingEstimator : 685/1296 [00:01<00:00, 658.53it/s] 55%|█████▍ | Scoring GeneralizingEstimator : 708/1296 [00:01<00:00, 658.53it/s] 56%|█████▋ | Scoring GeneralizingEstimator : 731/1296 [00:01<00:00, 658.80it/s] 58%|█████▊ | Scoring GeneralizingEstimator : 754/1296 [00:01<00:00, 658.88it/s] 60%|█████▉ | Scoring GeneralizingEstimator : 776/1296 [00:01<00:00, 658.33it/s] 62%|██████▏ | Scoring GeneralizingEstimator : 799/1296 [00:01<00:00, 659.35it/s] 63%|██████▎ | Scoring GeneralizingEstimator : 822/1296 [00:01<00:00, 659.05it/s] 65%|██████▌ | Scoring GeneralizingEstimator : 845/1296 [00:01<00:00, 659.09it/s] 67%|██████▋ | Scoring GeneralizingEstimator : 868/1296 [00:01<00:00, 659.04it/s] 69%|██████▉ | Scoring GeneralizingEstimator : 891/1296 [00:01<00:00, 659.37it/s] 70%|███████ | Scoring GeneralizingEstimator : 913/1296 [00:01<00:00, 657.75it/s] 72%|███████▏ | Scoring GeneralizingEstimator : 937/1296 [00:01<00:00, 659.22it/s] 74%|███████▍ | Scoring GeneralizingEstimator : 960/1296 [00:01<00:00, 659.34it/s] 76%|███████▌ | Scoring GeneralizingEstimator : 983/1296 [00:01<00:00, 659.68it/s] 78%|███████▊ | Scoring GeneralizingEstimator : 1006/1296 [00:01<00:00, 659.46it/s] 79%|███████▉ | Scoring GeneralizingEstimator : 1029/1296 [00:01<00:00, 659.60it/s] 81%|████████ | Scoring GeneralizingEstimator : 1052/1296 [00:01<00:00, 659.56it/s] 83%|████████▎ | Scoring GeneralizingEstimator : 1075/1296 [00:01<00:00, 659.27it/s] 85%|████████▍ | Scoring GeneralizingEstimator : 1098/1296 [00:01<00:00, 659.38it/s] 86%|████████▋ | Scoring GeneralizingEstimator : 1121/1296 [00:01<00:00, 659.68it/s] 88%|████████▊ | Scoring GeneralizingEstimator : 1144/1296 [00:01<00:00, 659.54it/s] 90%|████████▉ | Scoring GeneralizingEstimator : 1166/1296 [00:01<00:00, 659.02it/s] 92%|█████████▏| Scoring GeneralizingEstimator : 1189/1296 [00:01<00:00, 660.32it/s] 93%|█████████▎| Scoring GeneralizingEstimator : 1211/1296 [00:01<00:00, 659.89it/s] 95%|█████████▌| Scoring GeneralizingEstimator : 1234/1296 [00:01<00:00, 660.02it/s] 97%|█████████▋| Scoring GeneralizingEstimator : 1257/1296 [00:01<00:00, 659.85it/s] 99%|█████████▉| Scoring GeneralizingEstimator : 1280/1296 [00:01<00:00, 659.75it/s] 100%|██████████| Scoring GeneralizingEstimator : 1296/1296 [00:01<00:00, 660.01it/s] 100%|██████████| Scoring GeneralizingEstimator : 1296/1296 [00:01<00:00, 659.10it/s] 0%| | Fitting GeneralizingEstimator : 0/36 [00:00<?, ?it/s] 17%|█▋ | Fitting GeneralizingEstimator : 6/36 [00:00<00:00, 177.82it/s] 36%|███▌ | Fitting GeneralizingEstimator : 13/36 [00:00<00:00, 193.29it/s] 58%|█████▊ | Fitting GeneralizingEstimator : 21/36 [00:00<00:00, 207.70it/s] 78%|███████▊ | Fitting GeneralizingEstimator : 28/36 [00:00<00:00, 207.72it/s] 97%|█████████▋| Fitting GeneralizingEstimator : 35/36 [00:00<00:00, 206.78it/s] 100%|██████████| Fitting GeneralizingEstimator : 36/36 [00:00<00:00, 209.02it/s] 0%| | Scoring GeneralizingEstimator : 0/1296 [00:00<?, ?it/s] 2%|▏ | Scoring GeneralizingEstimator : 21/1296 [00:00<00:02, 620.84it/s] 3%|▎ | Scoring GeneralizingEstimator : 44/1296 [00:00<00:01, 645.76it/s] 5%|▌ | Scoring GeneralizingEstimator : 67/1296 [00:00<00:01, 649.51it/s] 7%|▋ | Scoring GeneralizingEstimator : 90/1296 [00:00<00:01, 650.85it/s] 9%|▊ | Scoring GeneralizingEstimator : 113/1296 [00:00<00:01, 652.13it/s] 10%|█ | Scoring GeneralizingEstimator : 136/1296 [00:00<00:01, 653.59it/s] 12%|█▏ | Scoring GeneralizingEstimator : 159/1296 [00:00<00:01, 654.31it/s] 14%|█▍ | Scoring GeneralizingEstimator : 182/1296 [00:00<00:01, 654.76it/s] 16%|█▌ | Scoring GeneralizingEstimator : 205/1296 [00:00<00:01, 654.46it/s] 18%|█▊ | Scoring GeneralizingEstimator : 227/1296 [00:00<00:01, 654.41it/s] 19%|█▉ | Scoring GeneralizingEstimator : 249/1296 [00:00<00:01, 654.39it/s] 21%|██ | Scoring GeneralizingEstimator : 271/1296 [00:00<00:01, 654.41it/s] 23%|██▎ | Scoring GeneralizingEstimator : 294/1296 [00:00<00:01, 654.44it/s] 24%|██▍ | Scoring GeneralizingEstimator : 317/1296 [00:00<00:01, 654.65it/s] 26%|██▌ | Scoring GeneralizingEstimator : 340/1296 [00:00<00:01, 654.74it/s] 28%|██▊ | Scoring GeneralizingEstimator : 363/1296 [00:00<00:01, 654.58it/s] 30%|██▉ | Scoring GeneralizingEstimator : 385/1296 [00:00<00:01, 654.57it/s] 31%|███▏ | Scoring GeneralizingEstimator : 408/1296 [00:00<00:01, 655.03it/s] 33%|███▎ | Scoring GeneralizingEstimator : 431/1296 [00:00<00:01, 654.90it/s] 35%|███▌ | Scoring GeneralizingEstimator : 454/1296 [00:00<00:01, 655.20it/s] 37%|███▋ | Scoring GeneralizingEstimator : 477/1296 [00:00<00:01, 655.74it/s] 39%|███▊ | Scoring GeneralizingEstimator : 500/1296 [00:00<00:01, 656.09it/s] 40%|████ | Scoring GeneralizingEstimator : 522/1296 [00:00<00:01, 655.17it/s] 42%|████▏ | Scoring GeneralizingEstimator : 544/1296 [00:00<00:01, 653.75it/s] 44%|████▎ | Scoring GeneralizingEstimator : 566/1296 [00:00<00:01, 652.04it/s] 45%|████▌ | Scoring GeneralizingEstimator : 588/1296 [00:00<00:01, 650.10it/s] 47%|████▋ | Scoring GeneralizingEstimator : 611/1296 [00:00<00:01, 650.37it/s] 49%|████▉ | Scoring GeneralizingEstimator : 634/1296 [00:00<00:01, 651.21it/s] 51%|█████ | Scoring GeneralizingEstimator : 657/1296 [00:01<00:00, 651.29it/s] 52%|█████▏ | Scoring GeneralizingEstimator : 679/1296 [00:01<00:00, 651.18it/s] 54%|█████▍ | Scoring GeneralizingEstimator : 702/1296 [00:01<00:00, 651.56it/s] 56%|█████▌ | Scoring GeneralizingEstimator : 725/1296 [00:01<00:00, 652.16it/s] 58%|█████▊ | Scoring GeneralizingEstimator : 748/1296 [00:01<00:00, 652.19it/s] 59%|█████▉ | Scoring GeneralizingEstimator : 771/1296 [00:01<00:00, 652.60it/s] 61%|██████▏ | Scoring GeneralizingEstimator : 794/1296 [00:01<00:00, 652.81it/s] 63%|██████▎ | Scoring GeneralizingEstimator : 817/1296 [00:01<00:00, 652.81it/s] 65%|██████▍ | Scoring GeneralizingEstimator : 840/1296 [00:01<00:00, 653.05it/s] 67%|██████▋ | Scoring GeneralizingEstimator : 863/1296 [00:01<00:00, 653.06it/s] 68%|██████▊ | Scoring GeneralizingEstimator : 886/1296 [00:01<00:00, 653.03it/s] 70%|███████ | Scoring GeneralizingEstimator : 909/1296 [00:01<00:00, 653.36it/s] 72%|███████▏ | Scoring GeneralizingEstimator : 931/1296 [00:01<00:00, 652.70it/s] 74%|███████▎ | Scoring GeneralizingEstimator : 954/1296 [00:01<00:00, 652.72it/s] 75%|███████▌ | Scoring GeneralizingEstimator : 977/1296 [00:01<00:00, 653.02it/s] 77%|███████▋ | Scoring GeneralizingEstimator : 1000/1296 [00:01<00:00, 653.13it/s] 79%|███████▉ | Scoring GeneralizingEstimator : 1023/1296 [00:01<00:00, 653.61it/s] 81%|████████ | Scoring GeneralizingEstimator : 1046/1296 [00:01<00:00, 653.62it/s] 82%|████████▏ | Scoring GeneralizingEstimator : 1068/1296 [00:01<00:00, 653.47it/s] 84%|████████▍ | Scoring GeneralizingEstimator : 1091/1296 [00:01<00:00, 653.71it/s] 86%|████████▌ | Scoring GeneralizingEstimator : 1114/1296 [00:01<00:00, 653.98it/s] 88%|████████▊ | Scoring GeneralizingEstimator : 1137/1296 [00:01<00:00, 653.94it/s] 90%|████████▉ | Scoring GeneralizingEstimator : 1160/1296 [00:01<00:00, 654.56it/s] 91%|█████████▏| Scoring GeneralizingEstimator : 1183/1296 [00:01<00:00, 655.01it/s] 93%|█████████▎| Scoring GeneralizingEstimator : 1206/1296 [00:01<00:00, 655.19it/s] 95%|█████████▍| Scoring GeneralizingEstimator : 1229/1296 [00:01<00:00, 655.56it/s] 97%|█████████▋| Scoring GeneralizingEstimator : 1252/1296 [00:01<00:00, 655.76it/s] 98%|█████████▊| Scoring GeneralizingEstimator : 1275/1296 [00:01<00:00, 655.73it/s] 100%|██████████| Scoring GeneralizingEstimator : 1296/1296 [00:01<00:00, 656.61it/s] 100%|██████████| Scoring GeneralizingEstimator : 1296/1296 [00:01<00:00, 654.44it/s]
Plot the full (generalization) matrix:
fig , ax = plt.subplots (1, 1) im = ax.imshow ( scores , interpolation="lanczos", origin="lower", cmap="RdBu_r", extent=epochs.times[[0, -1, 0, -1]], vmin=0.0, vmax=1.0, ) ax.set_xlabel ("Testing Time (s)") ax.set_ylabel ("Training Time (s)") ax.set_title ("Temporal generalization") ax.axvline (0, color="k") ax.axhline (0, color="k") cbar = plt.colorbar (im , ax =ax ) cbar.set_label ("AUC")
Projecting sensor-space patterns to source space#
If you use a linear classifier (or regressor) for your data, you can also
project these to source space. For example, using our evoked_time_gen
from before:
cov = mne.compute_covariance (epochs, tmax =0.0) del epochs fwd = mne.read_forward_solution (meg_path / "sample_audvis-meg-eeg-oct-6-fwd.fif") inv = mne.minimum_norm.make_inverse_operator (evoked_time_gen.info , fwd, cov , loose=0.0) stc = mne.minimum_norm.apply_inverse (evoked_time_gen , inv, 1.0 / 9.0, "dSPM") del fwd, inv
Reducing data rank from 203 -> 203 Estimating covariance using EMPIRICAL Done. Number of samples used : 1353 [done] Reading forward solution from /home/circleci/mne_data/MNE-sample-data/MEG/sample/sample_audvis-meg-eeg-oct-6-fwd.fif... Reading a source space... Computing patch statistics... Patch information added... Distance information added... [done] Reading a source space... Computing patch statistics... Patch information added... Distance information added... [done] 2 source spaces read Desired named matrix (kind = 3523 (FIFF_MNE_FORWARD_SOLUTION_GRAD)) not available Read MEG forward solution (7498 sources, 306 channels, free orientations) Desired named matrix (kind = 3523 (FIFF_MNE_FORWARD_SOLUTION_GRAD)) not available Read EEG forward solution (7498 sources, 60 channels, free orientations) Forward solutions combined: MEG, EEG Source spaces transformed to the forward solution coordinate frame Computing inverse operator with 203 channels. 203 out of 366 channels remain after picking Selected 203 channels Creating the depth weighting matrix... 203 planar channels limit = 7262/7498 = 10.020866 scale = 2.58122e-08 exp = 0.8 Picked elements from a free-orientation depth-weighting prior into the fixed-orientation one Average patch normals will be employed in the rotation to the local surface coordinates.... Converting to surface-based source orientations... [done] Whitening the forward solution. Computing rank from covariance with rank=None Using tolerance 1.6e-13 (2.2e-16 eps * 203 dim * 3.6 max singular value) Estimated rank (grad): 203 GRAD: rank 203 computed from 203 data channels with 0 projectors Setting small GRAD eigenvalues to zero (without PCA) Creating the source covariance matrix Adjusting source covariance matrix. Computing SVD of whitened and weighted lead field matrix. largest singular value = 3.91709 scaling factor to adjust the trace = 6.26373e+18 (nchan = 203 nzero = 0) Preparing the inverse operator for use... Scaled noise and source covariance from nave = 1 to nave = 1 Created the regularized inverter The projection vectors do not apply to these channels. Created the whitener using a noise covariance matrix with rank 203 (0 small eigenvalues omitted) Computing noise-normalization factors (dSPM)... [done] Applying inverse operator to ""... Picked 203 channels from the data Computing inverse... Eigenleads need to be weighted ... Computing residual... Explained 76.4% variance dSPM... [done]
And this can be visualized using stc.plot:
brain = stc.plot ( hemi="split", views=("lat", "med"), initial_time=0.1, subjects_dir =subjects_dir )
Using control points [1.98776221 2.41838256 8.06628583]
Source-space decoding#
Source space decoding is also possible, but because the number of features can be much larger than in the sensor space, univariate feature selection using ANOVA f-test (or some other metric) can be done to reduce the feature dimension. Interpreting decoding results might be easier in source space as compared to sensor space.
Exercise#
Explore other datasets from MNE (e.g. Face dataset from SPM to predict Face vs. Scrambled)
References#
Total running time of the script: (0 minutes 17.950 seconds)