.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/plot_power.py" .. LINE NUMBERS ARE GIVEN BELOW. .. 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_power.py: .. _power: =========================== Compute PSD on sensor space =========================== The power pipeline computes the power spectral density (PSD) on epochs or raw data on **sensor space** or **source space**. The **mean PSD** for each selected frequency band is also computed and saved in a numpy file. The input data shoud be in **fif** or **numpy** format. .. GENERATED FROM PYTHON SOURCE LINES 14-30 .. code-block:: default # Authors: Annalisa Pascarella # Mainak Jas # License: BSD (3-clause) # sphinx_gallery_thumbnail_number = 2 import os.path as op import numpy as np import nipype.pipeline.engine as pe import ephypype from ephypype.nodes import create_iterator, create_datagrabber from ephypype.datasets import fetch_omega_dataset .. GENERATED FROM PYTHON SOURCE LINES 31-32 Let us fetch the data first. It is around 675 MB download. .. GENERATED FROM PYTHON SOURCE LINES 32-36 .. code-block:: default base_path = op.join(op.dirname(ephypype.__file__), '..', 'examples') data_path = fetch_omega_dataset(base_path) .. GENERATED FROM PYTHON SOURCE LINES 37-40 then read the parameters for experiment and power analysis from a :download:`json ` file and print it .. GENERATED FROM PYTHON SOURCE LINES 40-58 .. code-block:: default import json # noqa import pprint # noqa params = json.load(open("params.json")) pprint.pprint({'experiment parameters': params["general"]}) subject_ids = params["general"]["subject_ids"] # sub-003 session_ids = params["general"]["session_ids"] # ses-0001 NJOBS = params["general"]["NJOBS"] pprint.pprint({'power parameters': params["power"]}) freq_band_names = params["power"]['freq_band_names'] freq_bands = params["power"]['freq_bands'] is_epoched = params["power"]['is_epoched'] fmin = params["power"]['fmin'] fmax = params["power"]['fmax'] power_method = params["power"]['method'] .. rst-class:: sphx-glr-script-out .. code-block:: none {'experiment parameters': {'NJOBS': 1, 'data_type': 'fif', 'session_ids': ['ses-0001'], 'subject_ids': ['sub-0003'], 'subjects_dir': 'FSF'}} {'power parameters': {'fmax': 150, 'fmin': 0.1, 'freq_band_names': ['theta', 'alpha', 'beta'], 'freq_bands': [[3, 6], [8, 13], [13, 30]], 'is_epoched': False, 'method': 'welch'}} .. GENERATED FROM PYTHON SOURCE LINES 59-61 Then, we create our workflow and specify the `base_dir` which tells nipype the directory in which to store the outputs. .. GENERATED FROM PYTHON SOURCE LINES 61-68 .. code-block:: default # workflow directory within the `base_dir` power_analysis_name = 'power_workflow' main_workflow = pe.Workflow(name=power_analysis_name) main_workflow.base_dir = data_path .. GENERATED FROM PYTHON SOURCE LINES 69-70 Then we create a node to pass input filenames to DataGrabber from nipype .. GENERATED FROM PYTHON SOURCE LINES 70-74 .. code-block:: default infosource = create_iterator(['subject_id', 'session_id'], [subject_ids, session_ids]) .. GENERATED FROM PYTHON SOURCE LINES 75-77 and a node to grab data. The template_args in this node iterate upon the values in the infosource node .. GENERATED FROM PYTHON SOURCE LINES 77-82 .. code-block:: default template_path = '*%s/%s/meg/%s*rest*0_60*ica.fif' template_args = [['subject_id', 'session_id', 'subject_id']] datasource = create_datagrabber(data_path, template_path, template_args) .. GENERATED FROM PYTHON SOURCE LINES 83-93 Ephypype creates for us a pipeline which can be connected to these nodes we created. The power pipeline in the **sensor space** is implemented by the function :func:`ephypype.pipelines.power.create_pipeline_power`, thus to instantiate this pipeline node, we import it and pass our parameters to it. The power pipeline contains only one node :class:`ephypype.interfaces.mne.power.Power` that wraps the MNE-Python functions :func:`mne.time_frequency.psd_welch` and :func:`mne.time_frequency.psd_multitaper` for computing the PSD using Welch's method and multitapers respectively. .. GENERATED FROM PYTHON SOURCE LINES 93-100 .. code-block:: default from ephypype.pipelines import create_pipeline_power # noqa power_workflow = create_pipeline_power(data_path, freq_bands, fmin=fmin, fmax=fmax, method=power_method, is_epoched=is_epoched) .. rst-class:: sphx-glr-script-out .. code-block:: none *** main_path -> /home/pasca/Tools/python/packages/neuropycon/ephypype/examples/sample_BIDS_omega *** .. GENERATED FROM PYTHON SOURCE LINES 101-104 We then connect the nodes two at a time. First, we connect the two outputs (subject_id and session_id) of the infosource node to the datasource node. So, these two nodes taken together can grab data. .. GENERATED FROM PYTHON SOURCE LINES 104-108 .. code-block:: default main_workflow.connect(infosource, 'subject_id', datasource, 'subject_id') main_workflow.connect(infosource, 'session_id', datasource, 'session_id') .. GENERATED FROM PYTHON SOURCE LINES 109-111 Similarly, for the inputnode of the power_workflow. Things will become clearer in a moment when we plot the graph of the workflow. .. GENERATED FROM PYTHON SOURCE LINES 111-116 .. code-block:: default main_workflow.connect(datasource, 'raw_file', power_workflow, 'inputnode.fif_file') .. GENERATED FROM PYTHON SOURCE LINES 117-118 To do so, we first write the workflow graph (optional) .. GENERATED FROM PYTHON SOURCE LINES 118-121 .. code-block:: default main_workflow.write_graph(graph2use='colored') # colored .. rst-class:: sphx-glr-script-out .. code-block:: none '/home/pasca/Tools/python/packages/neuropycon/ephypype/examples/sample_BIDS_omega/power_workflow/graph.png' .. GENERATED FROM PYTHON SOURCE LINES 122-124 and visualize it. Take a moment to pause and notice how the connections here correspond to how we connected the nodes. .. GENERATED FROM PYTHON SOURCE LINES 124-131 .. code-block:: default import matplotlib.pyplot as plt # noqa img = plt.imread(op.join(data_path, power_analysis_name, 'graph.png')) plt.figure(figsize=(6, 6)) plt.imshow(img) plt.axis('off') .. image-sg:: /auto_examples/images/sphx_glr_plot_power_001.png :alt: plot power :srcset: /auto_examples/images/sphx_glr_plot_power_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none (-0.5, 404.5, 498.5, -0.5) .. GENERATED FROM PYTHON SOURCE LINES 132-133 Finally, we are now ready to execute our workflow. .. GENERATED FROM PYTHON SOURCE LINES 133-139 .. code-block:: default main_workflow.config['execution'] = {'remove_unnecessary_outputs': 'false'} # Run workflow locally on 1 CPU main_workflow.run(plugin='MultiProc', plugin_args={'n_procs': NJOBS}) .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 140-148 The outputs are the **psd tensor and frequencies in .npz format** and the **mean PSD in .npy format** stored in the workflow directory defined by `base_dir` .. note:: The power pipeline in the **source space** is implemented by the function :func:`ephypype.pipelines.power.create_pipeline_power_src_space` and its Node :class:`ephypype.interfaces.mne.power.Power` compute the PSD by the welch function of the scipy package. .. GENERATED FROM PYTHON SOURCE LINES 150-189 .. code-block:: default from ephypype.gather import get_results # noqa from visbrain.objects import SourceObj, SceneObj, ColorbarObj # noqa from visbrain.utils import normalize # noqa from nipype.utils.filemanip import split_filename # noqa psd_files, channel_coo_files = get_results(main_workflow.base_dir, main_workflow.name, pipeline='power') sc = SceneObj(size=(1800, 500), bgcolor=(.1, .1, .1)) for psd_file, channel_coo_file in zip(psd_files, channel_coo_files): path_xyz, basename, ext = split_filename(psd_file) arch = np.load(psd_file) psds, freqs = arch['psds'], arch['freqs'] xyz = np.genfromtxt(channel_coo_file, dtype=float) freq_bands = np.asarray(freq_bands) clim = (psds.min(), psds.max()) cmap = 'cool' txtcolor = 'white' # Find indices of frequencies : idx_fplt = np.abs((freqs.reshape(1, 1, -1) - freq_bands[..., np.newaxis])).argmin(2) psdf = np.array([psds[:, k[0]:k[1]].mean(1) for k in idx_fplt]) radius = normalize(np.c_[psdf.min(1), psdf.max(1)], 5, 25).astype(float) for num, (fb, fbn, psd, rx) in enumerate(zip(freq_bands, freq_band_names, psdf, radius)): s_obj = SourceObj('s', xyz, data=psd, radius_min=rx[0], radius_max=rx[1]) # noqa s_obj.color_sources(data=psd, cmap=cmap, clim=clim) sc.add_to_subplot(s_obj, col=num, title=str(fb) + ' - ' + fbn, title_color=txtcolor, rotate='top', zoom=.6) cbar = ColorbarObj(s_obj, txtcolor=txtcolor, cblabel='PSD', txtsz=15, cbtxtsz=20) sc.add_to_subplot(cbar, col=len(freq_bands), width_max=200) sc.preview() .. image-sg:: /auto_examples/images/sphx_glr_plot_power_002.png :alt: plot power :srcset: /auto_examples/images/sphx_glr_plot_power_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 0 minutes 7.470 seconds) .. _sphx_glr_download_auto_examples_plot_power.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_power.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_power.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_