.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/plot_02_sphara_basis_eeg.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_plot_02_sphara_basis_eeg.py: .. _sphara_basis_eeg: Determination of the SPHARA basis functions for an EEG sensor setup =================================================================== .. topic:: Section contents This tutorial introduces the steps necessary to determine a generalized spatial Fourier basis for an :term:`EEG` sensor setup using SpharaPy. The special properties of the different discretization approaches of the Laplace–Beltrami operator will be discussed. Introduction ------------ A Fourier basis can be obtained as a solution to Laplace's eigenvalue problem .. math:: \mathbf{L} \boldsymbol{\varphi}_i = \lambda_i \boldsymbol{\varphi}_i\,,\qquad\qquad(1) where :math:`\mathbf{L}` is a discrete :term:`Laplace–Beltrami operator` in matrix notation, the eigenvectors :math:`\boldsymbol{\varphi}_i` contain the spatial harmonic functions and the eigenvalues :math:`\lambda_i \ge 0` are real-valued. When :math:`\mathbf{L}` is obtained from a finite element (FEM) discretization of the Laplace–Beltrami operator, the quantities :math:`\sqrt{\lambda_i}` can be interpreted as spatial angular frequencies; see :ref:`introduction` for details. By solving a Laplace eigenvalue problem, it is also possible to determine a basis for a spatial Fourier analysis. Often in practical applications, a measured quantity to be subjected to a Fourier analysis is only known at spatially discrete sampling points (the sensor positions). An arbitrary arrangement of sample points on a surface in three-dimensional space can be described by means of a :term:`triangular mesh`. In the case of an :term:`EEG` system, the sample positions (the vertices of the triangular mesh) are the locations of the sensors arranged on the head surface. A SPHARA basis is a solution of a Laplace eigenvalue problem for the given triangle mesh, that can be obtained by discretizing a Laplace–Beltrami operator for the mesh and solving the Laplace eigenvalue problem in equation (1). The SpharaPy package provides three methods for the discretization of the Laplace–Beltrami operator; unit weighting of the edges and weighting with the inverse of the Euclidean distance of the edges, and a FEM approach. For more detailed information please refer to :ref:`introduction` and :ref:`eigensystems_of_lb_operators`. .. GENERATED FROM PYTHON SOURCE LINES 56-58 At the beginning we import three modules of the SpharaPy package as well as several other packages and single functions of packages. .. GENERATED FROM PYTHON SOURCE LINES 58-72 .. code-block:: Python # Code source: Uwe Graichen # License: BSD 3 clause # import modules from spharapy package # import additional modules used in this tutorial import matplotlib.pyplot as plt import numpy as np from mpl_toolkits.mplot3d import Axes3D # noqa: F401 (registers 3D) import spharapy.datasets as sd import spharapy.spharabasis as sb import spharapy.trimesh as tm .. GENERATED FROM PYTHON SOURCE LINES 73-84 Specification of the spatial configuration of the EEG sensors ------------------------------------------------------------- Import information about EEG sensor setup of the sample data set ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ In this tutorial we will determine a SPHARA basis for a 256 channel :term:`EEG` system with equidistant layout. The data set is one of the example data sets contained in the SpharaPy toolbox, see :mod:`spharapy.datasets` and :func:`spharapy.datasets.load_eeg_256_channel_study`. .. GENERATED FROM PYTHON SOURCE LINES 84-88 .. code-block:: Python # loading the 256 channel EEG dataset from spharapy sample datasets mesh_in = sd.load_eeg_256_channel_study() .. GENERATED FROM PYTHON SOURCE LINES 89-93 The dataset includes lists of vertices, triangles, and sensor labels, as well as :term:`EEG` data from previously performed experiment addressing the cortical activation related to somatosensory-evoked potentials (SEP). .. GENERATED FROM PYTHON SOURCE LINES 93-96 .. code-block:: Python print(mesh_in.keys()) .. rst-class:: sphx-glr-script-out .. code-block:: none dict_keys(['vertlist', 'trilist', 'labellist', 'eegdata']) .. GENERATED FROM PYTHON SOURCE LINES 97-99 The triangulation of the :term:`EEG` sensor setup consists of 256 vertices and 480 triangles. .. GENERATED FROM PYTHON SOURCE LINES 99-105 .. code-block:: Python vertlist = np.array(mesh_in["vertlist"]) trilist = np.array(mesh_in["trilist"]) print("vertices = ", vertlist.shape) print("triangles = ", trilist.shape) .. rst-class:: sphx-glr-script-out .. code-block:: none vertices = (256, 3) triangles = (482, 3) .. GENERATED FROM PYTHON SOURCE LINES 106-128 .. code-block:: Python fig = plt.figure() fig.subplots_adjust(left=0.02, right=0.98, top=0.98, bottom=0.02) ax = fig.add_subplot(111, projection="3d") ax.set_xlabel("x") ax.set_ylabel("y") ax.set_zlabel("z") ax.set_title("The triangulated EEG sensor setup") ax.view_init(elev=20.0, azim=80.0) ax.set_aspect("auto") ax.plot_trisurf( vertlist[:, 0], vertlist[:, 1], vertlist[:, 2], triangles=trilist, color="lightblue", edgecolor="black", linewidth=0.5, shade=True, ) plt.show() .. image-sg:: /auto_examples/images/sphx_glr_plot_02_sphara_basis_eeg_001.png :alt: The triangulated EEG sensor setup :srcset: /auto_examples/images/sphx_glr_plot_02_sphara_basis_eeg_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 129-134 Create a SpharaPy TriMesh instance ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ In the next step we create an instance of the class :class:`spharapy.trimesh.TriMesh` from the list of vertices and triangles. .. GENERATED FROM PYTHON SOURCE LINES 134-138 .. code-block:: Python # create an instance of the TriMesh class mesh_eeg = tm.TriMesh(trilist, vertlist) .. GENERATED FROM PYTHON SOURCE LINES 139-142 The class :class:`spharapy.trimesh.TriMesh` provides a number of methods to determine certain properties of the triangle mesh required to generate the SPHARA basis, listed below: .. GENERATED FROM PYTHON SOURCE LINES 142-146 .. code-block:: Python # print all implemented methods of the TriMesh class print([func for func in dir(tm.TriMesh) if not func.startswith("__")]) .. rst-class:: sphx-glr-script-out .. code-block:: none ['adjacent_tri', 'is_edge', 'laplacianmatrix', 'massmatrix', 'one_ring_neighborhood', 'remove_vertices', 'stiffnessmatrix', 'trilist', 'vertlist', 'weightmatrix'] .. GENERATED FROM PYTHON SOURCE LINES 147-159 Determining SPHARA bases using different discretisation approaches ------------------------------------------------------------------ Computing the basis functions ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ In the final step of the tutorial we will calculate SPHARA bases for the given :term:`EEG` sensor setup. For this we create three instances of the class :class:`spharapy.spharabasis.SpharaBasis`. We use the three discretization approaches implemented in this class for the Laplace–Beltrami operator: unit weighting ('unit') and inverse Euclidean weigthing ('inv_euclidean') of the edges of the triangular mesh as well as the FEM discretization ('fem') .. GENERATED FROM PYTHON SOURCE LINES 159-172 .. code-block:: Python # 'unit' discretization sphara_basis_unit = sb.SpharaBasis(mesh_eeg, "unit") basis_functions_unit, natural_frequencies_unit = sphara_basis_unit.basis() # 'inv_euclidean' discretization sphara_basis_ie = sb.SpharaBasis(mesh_eeg, "inv_euclidean") basis_functions_ie, natural_frequencies_ie = sphara_basis_ie.basis() # 'fem' discretization sphara_basis_fem = sb.SpharaBasis(mesh_eeg, "fem") basis_functions_fem, natural_frequencies_fem = sphara_basis_fem.basis() .. GENERATED FROM PYTHON SOURCE LINES 173-180 Visualization the basis functions ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ The first 15 spatially low-frequency SPHARA basis functions are shown below, starting with DC at the top left. SPHARA basis using the discretization approache 'unit' """""""""""""""""""""""""""""""""""""""""""""""""""""" .. GENERATED FROM PYTHON SOURCE LINES 180-215 .. code-block:: Python figsb1, axes1 = plt.subplots(nrows=5, ncols=3, figsize=(8, 12), subplot_kw={"projection": "3d"}) for i in range(np.size(axes1)): colors = np.mean(basis_functions_unit[trilist, i + 0], axis=1) ax = axes1.flat[i] ax.set_xlabel("x") ax.set_ylabel("y") ax.set_zlabel("z") ax.view_init(elev=60.0, azim=80.0) ax.set_aspect("auto") trisurfplot = ax.plot_trisurf( vertlist[:, 0], vertlist[:, 1], vertlist[:, 2], triangles=trilist, cmap=plt.cm.bwr, edgecolor="white", linewidth=0.0, ) trisurfplot.set_array(colors) trisurfplot.set_clim(-0.15, 0.15) cbar = figsb1.colorbar( trisurfplot, ax=axes1.ravel().tolist(), shrink=0.85, orientation="horizontal", fraction=0.05, pad=0.05, anchor=(0.5, -4.5), ) plt.subplots_adjust(left=0.0, right=1.0, bottom=0.08, top=1.0) plt.show() .. image-sg:: /auto_examples/images/sphx_glr_plot_02_sphara_basis_eeg_002.png :alt: plot 02 sphara basis eeg :srcset: /auto_examples/images/sphx_glr_plot_02_sphara_basis_eeg_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 217-219 SPHARA basis using the discretization approache 'inv_euclidean' """"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" .. GENERATED FROM PYTHON SOURCE LINES 219-254 .. code-block:: Python figsb1, axes1 = plt.subplots(nrows=5, ncols=3, figsize=(8, 12), subplot_kw={"projection": "3d"}) for i in range(np.size(axes1)): colors = np.mean(basis_functions_ie[trilist, i + 0], axis=1) ax = axes1.flat[i] ax.set_xlabel("x") ax.set_ylabel("y") ax.set_zlabel("z") ax.view_init(elev=60.0, azim=80.0) ax.set_aspect("auto") trisurfplot = ax.plot_trisurf( vertlist[:, 0], vertlist[:, 1], vertlist[:, 2], triangles=trilist, cmap=plt.cm.bwr, edgecolor="white", linewidth=0.0, ) trisurfplot.set_array(colors) trisurfplot.set_clim(-0.18, 0.18) cbar = figsb1.colorbar( trisurfplot, ax=axes1.ravel().tolist(), shrink=0.85, orientation="horizontal", fraction=0.05, pad=0.05, anchor=(0.5, -4.5), ) plt.subplots_adjust(left=0.0, right=1.0, bottom=0.08, top=1.0) plt.show() .. image-sg:: /auto_examples/images/sphx_glr_plot_02_sphara_basis_eeg_003.png :alt: plot 02 sphara basis eeg :srcset: /auto_examples/images/sphx_glr_plot_02_sphara_basis_eeg_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 255-257 SPHARA basis using the discretization approache 'fem' """"""""""""""""""""""""""""""""""""""""""""""""""""" .. GENERATED FROM PYTHON SOURCE LINES 257-291 .. code-block:: Python figsb1, axes1 = plt.subplots(nrows=5, ncols=3, figsize=(8, 12), subplot_kw={"projection": "3d"}) for i in range(np.size(axes1)): colors = np.mean(basis_functions_fem[trilist, i + 0], axis=1) ax = axes1.flat[i] ax.set_xlabel("x") ax.set_ylabel("y") ax.set_zlabel("z") ax.view_init(elev=60.0, azim=80.0) ax.set_aspect("auto") trisurfplot = ax.plot_trisurf( vertlist[:, 0], vertlist[:, 1], vertlist[:, 2], triangles=trilist, cmap=plt.cm.bwr, edgecolor="white", linewidth=0.0, ) trisurfplot.set_array(colors) trisurfplot.set_clim(-0.01, 0.01) cbar = figsb1.colorbar( trisurfplot, ax=axes1.ravel().tolist(), shrink=0.85, orientation="horizontal", fraction=0.05, pad=0.05, anchor=(0.5, -4.5), ) plt.subplots_adjust(left=0.0, right=1.0, bottom=0.08, top=1.0) plt.show() .. image-sg:: /auto_examples/images/sphx_glr_plot_02_sphara_basis_eeg_004.png :alt: plot 02 sphara basis eeg :srcset: /auto_examples/images/sphx_glr_plot_02_sphara_basis_eeg_004.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 2.811 seconds) .. _sphx_glr_download_auto_examples_plot_02_sphara_basis_eeg.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_02_sphara_basis_eeg.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_02_sphara_basis_eeg.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_02_sphara_basis_eeg.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_