{ "nbformat_minor": 0, "nbformat": 4, "cells": [ { "execution_count": null, "cell_type": "code", "source": [ "%matplotlib inline" ], "outputs": [], "metadata": { "collapsed": false } }, { "source": [ "\n# Compute source space connectivity and visualize it using a circular graph\n\n\nThis example computes the all-to-all connectivity between 68 regions in\nsource space based on dSPM inverse solutions and a FreeSurfer cortical\nparcellation. The connectivity is visualized using a circular graph which\nis ordered based on the locations of the regions.\n\n" ], "cell_type": "markdown", "metadata": {} }, { "execution_count": null, "cell_type": "code", "source": [ "# Authors: Martin Luessi \n# Alexandre Gramfort \n# Nicolas P. Rougier (graph code borrowed from his matplotlib gallery)\n#\n# License: BSD (3-clause)\n\nimport numpy as np\nimport matplotlib.pyplot as plt\n\nimport mne\nfrom mne.datasets import sample\nfrom mne.minimum_norm import apply_inverse_epochs, read_inverse_operator\nfrom mne.connectivity import spectral_connectivity\nfrom mne.viz import circular_layout, plot_connectivity_circle\n\nprint(__doc__)\n\ndata_path = sample.data_path()\nsubjects_dir = data_path + '/subjects'\nfname_inv = data_path + '/MEG/sample/sample_audvis-meg-oct-6-meg-inv.fif'\nfname_raw = data_path + '/MEG/sample/sample_audvis_filt-0-40_raw.fif'\nfname_event = data_path + '/MEG/sample/sample_audvis_filt-0-40_raw-eve.fif'\n\n# Load data\ninverse_operator = read_inverse_operator(fname_inv)\nraw = mne.io.read_raw_fif(fname_raw)\nevents = mne.read_events(fname_event)\n\n# Add a bad channel\nraw.info['bads'] += ['MEG 2443']\n\n# Pick MEG channels\npicks = mne.pick_types(raw.info, meg=True, eeg=False, stim=False, eog=True,\n exclude='bads')\n\n# Define epochs for left-auditory condition\nevent_id, tmin, tmax = 1, -0.2, 0.5\nepochs = mne.Epochs(raw, events, event_id, tmin, tmax, picks=picks,\n baseline=(None, 0), reject=dict(mag=4e-12, grad=4000e-13,\n eog=150e-6))\n\n# Compute inverse solution and for each epoch. By using \"return_generator=True\"\n# stcs will be a generator object instead of a list.\nsnr = 1.0 # use lower SNR for single epochs\nlambda2 = 1.0 / snr ** 2\nmethod = \"dSPM\" # use dSPM method (could also be MNE or sLORETA)\nstcs = apply_inverse_epochs(epochs, inverse_operator, lambda2, method,\n pick_ori=\"normal\", return_generator=True)\n\n# Get labels for FreeSurfer 'aparc' cortical parcellation with 34 labels/hemi\nlabels = mne.read_labels_from_annot('sample', parc='aparc',\n subjects_dir=subjects_dir)\nlabel_colors = [label.color for label in labels]\n\n# Average the source estimates within each label using sign-flips to reduce\n# signal cancellations, also here we return a generator\nsrc = inverse_operator['src']\nlabel_ts = mne.extract_label_time_course(stcs, labels, src, mode='mean_flip',\n return_generator=True)\n\n# Now we are ready to compute the connectivity in the alpha band. Notice\n# from the status messages, how mne-python: 1) reads an epoch from the raw\n# file, 2) applies SSP and baseline correction, 3) computes the inverse to\n# obtain a source estimate, 4) averages the source estimate to obtain a\n# time series for each label, 5) includes the label time series in the\n# connectivity computation, and then moves to the next epoch. This\n# behaviour is because we are using generators and allows us to\n# compute connectivity in computationally efficient manner where the amount\n# of memory (RAM) needed is independent from the number of epochs.\nfmin = 8.\nfmax = 13.\nsfreq = raw.info['sfreq'] # the sampling frequency\ncon_methods = ['pli', 'wpli2_debiased']\ncon, freqs, times, n_epochs, n_tapers = spectral_connectivity(\n label_ts, method=con_methods, mode='multitaper', sfreq=sfreq, fmin=fmin,\n fmax=fmax, faverage=True, mt_adaptive=True, n_jobs=1)\n\n# con is a 3D array, get the connectivity for the first (and only) freq. band\n# for each method\ncon_res = dict()\nfor method, c in zip(con_methods, con):\n con_res[method] = c[:, :, 0]\n\n# Now, we visualize the connectivity using a circular graph layout\n\n# First, we reorder the labels based on their location in the left hemi\nlabel_names = [label.name for label in labels]\n\nlh_labels = [name for name in label_names if name.endswith('lh')]\n\n# Get the y-location of the label\nlabel_ypos = list()\nfor name in lh_labels:\n idx = label_names.index(name)\n ypos = np.mean(labels[idx].pos[:, 1])\n label_ypos.append(ypos)\n\n# Reorder the labels based on their location\nlh_labels = [label for (yp, label) in sorted(zip(label_ypos, lh_labels))]\n\n# For the right hemi\nrh_labels = [label[:-2] + 'rh' for label in lh_labels]\n\n# Save the plot order and create a circular layout\nnode_order = list()\nnode_order.extend(lh_labels[::-1]) # reverse the order\nnode_order.extend(rh_labels)\n\nnode_angles = circular_layout(label_names, node_order, start_pos=90,\n group_boundaries=[0, len(label_names) / 2])\n\n# Plot the graph using node colors from the FreeSurfer parcellation. We only\n# show the 300 strongest connections.\nplot_connectivity_circle(con_res['pli'], label_names, n_lines=300,\n node_angles=node_angles, node_colors=label_colors,\n title='All-to-All Connectivity left-Auditory '\n 'Condition (PLI)')\nplt.savefig('circle.png', facecolor='black')\n\n# Plot connectivity for both methods in the same plot\nfig = plt.figure(num=None, figsize=(8, 4), facecolor='black')\nno_names = [''] * len(label_names)\nfor ii, method in enumerate(con_methods):\n plot_connectivity_circle(con_res[method], no_names, n_lines=300,\n node_angles=node_angles, node_colors=label_colors,\n title=method, padding=0, fontsize_colorbar=6,\n fig=fig, subplot=(1, 2, ii + 1))\n\nplt.show()" ], "outputs": [], "metadata": { "collapsed": false } } ], "metadata": { "kernelspec": { "display_name": "Python 2", "name": "python2", "language": "python" }, "language_info": { "mimetype": "text/x-python", "nbconvert_exporter": "python", "name": "python", "file_extension": ".py", "version": "2.7.13", "pygments_lexer": "ipython2", "codemirror_mode": { "version": 2, "name": "ipython" } } } }