.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_tutorials/tutorial_07_anatomy_inference_accuracy.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_tutorials_tutorial_07_anatomy_inference_accuracy.py: Anatomical Predictors of Laminar Inference Accuracy =================================================== In [Szul et al., 2025, Beyond deep versus superficial: true laminar inference with MEG](https://doi.org/10.1101/2025.05.28.656642), we showed that even under ideal conditions (high SNR, perfect co-registration), the accuracy of laminar MEG source reconstruction varies regionally across the cortex. This tutorial demonstrates how to compute vertex-wise anatomical predictors that explain this spatial variability using `laMEG` utilities. .. GENERATED FROM PYTHON SOURCE LINES 9-12 Compute anatomical features --------------------------- .. GENERATED FROM PYTHON SOURCE LINES 12-45 .. code-block:: default import os import shutil import tempfile from IPython.display import Image import base64 import numpy as np from scipy.stats import zscore from lameg.invert import coregister, invert_ebb, get_lead_field_rms_diff from lameg.surf import LayerSurfaceSet from lameg.util import get_fiducial_coords from lameg.viz import show_surface, color_map import spm_standalone # Subject information for data to use subj_id = 'sub-104' ses_id = 'ses-01' # Fiducial coil coordinates fid_coords = get_fiducial_coords(subj_id, '../test_data/participants.tsv') # Data file to base simulations on data_file = os.path.join( '../test_data', subj_id, 'meg', ses_id, f'spm/spm-converted_autoreject-{subj_id}-{ses_id}-001-btn_trial-epo.mat' ) spm = spm_standalone.initialize() .. GENERATED FROM PYTHON SOURCE LINES 46-47 We'll use a forward model using the multilayer mesh .. GENERATED FROM PYTHON SOURCE LINES 47-50 .. code-block:: default surf_set = LayerSurfaceSet(subj_id, 11) .. GENERATED FROM PYTHON SOURCE LINES 51-52 We're going to copy the data file to a temporary directory and direct all output there. .. GENERATED FROM PYTHON SOURCE LINES 52-73 .. code-block:: default # Extract base name and path of data file data_path, data_file_name = os.path.split(data_file) data_base = os.path.splitext(data_file_name)[0] # Where to put simulated data tmp_dir = tempfile.mkdtemp() # Copy data files to tmp directory shutil.copy( os.path.join(data_path, f'{data_base}.mat'), os.path.join(tmp_dir, f'{data_base}.mat') ) shutil.copy( os.path.join(data_path, f'{data_base}.dat'), os.path.join(tmp_dir, f'{data_base}.dat') ) # Construct base file name for simulations base_fname = os.path.join(tmp_dir, f'{data_base}.mat') .. GENERATED FROM PYTHON SOURCE LINES 74-75 First we need to run source reconstruction in order to get the lead field information .. GENERATED FROM PYTHON SOURCE LINES 75-98 .. code-block:: default # Patch size to use for inversion (in this case it matches the simulated patch size) patch_size = 5 # Number of temporal modes to use for EBB inversion n_temp_modes = 4 # Coregister data to multilayer mesh coregister( fid_coords, base_fname, surf_set, spm_instance=spm ) # Run inversion [_,_] = invert_ebb( base_fname, surf_set, patch_size=patch_size, n_temp_modes=n_temp_modes, spm_instance=spm ) .. GENERATED FROM PYTHON SOURCE LINES 99-100 Compute vertex-wise anatomical predictors .. GENERATED FROM PYTHON SOURCE LINES 100-106 .. code-block:: default thickness = surf_set.get_cortical_thickness() # Cortical thickness (mm) lf_rms_diff = get_lead_field_rms_diff(base_fname, surf_set) # Lead-field variability across depth orientations = surf_set.get_radiality_to_scalp() # Column orientation relative to scalp distances = surf_set.get_distance_to_scalp() # Cortical distance to scalp (mm) .. GENERATED FROM PYTHON SOURCE LINES 107-108 Each predictor captures a distinct anatomical constraint on laminar inference: .. GENERATED FROM PYTHON SOURCE LINES 110-111 **Cortical thickness:** thicker patches of cortices have more separable laminar fields. .. GENERATED FROM PYTHON SOURCE LINES 111-134 .. code-block:: default c_range = [np.percentile(thickness,1), np.percentile(thickness,99)] # Plot change in power on each surface colors,_ = color_map( thickness, "Spectral_r", c_range[0], c_range[1], norm='N' ) cam_view = [335, 9.5, 51, 60, 37, 17, 0, 0, 1] plot = show_surface( surf_set, vertex_colors=colors, info=True, camera_view=cam_view ) .. GENERATED FROM PYTHON SOURCE LINES 135-138 .. image:: ../_static/tutorial_07_thickness.png :width: 800 :alt: .. GENERATED FROM PYTHON SOURCE LINES 140-141 **Lead-field RMS difference**: quantifies the difference in lead field patterns between the deepest and most superficial vertices; higher values indicate stronger layer separation. .. GENERATED FROM PYTHON SOURCE LINES 141-164 .. code-block:: default c_range = [np.percentile(lf_rms_diff,1), np.percentile(lf_rms_diff,99)] # Plot change in power on each surface colors,_ = color_map( lf_rms_diff, "Spectral_r", c_range[0], c_range[1], norm='N' ) cam_view = [335, 9.5, 51, 60, 37, 17, 0, 0, 1] plot = show_surface( surf_set, vertex_colors=colors, info=True, camera_view=cam_view ) .. GENERATED FROM PYTHON SOURCE LINES 165-168 .. image:: ../_static/tutorial_07_lf_rms_diff.png :width: 800 :alt: .. GENERATED FROM PYTHON SOURCE LINES 170-171 **Orientation**: measures dipole alignment with scalp normals (1 = radial, 0 = tangential). .. GENERATED FROM PYTHON SOURCE LINES 171-194 .. code-block:: default c_range = [np.percentile(orientations,1), np.percentile(orientations,99)] # Plot change in power on each surface colors,_ = color_map( orientations, "Spectral_r", c_range[0], c_range[1], norm='N' ) cam_view = [335, 9.5, 51, 60, 37, 17, 0, 0, 1] plot = show_surface( surf_set, vertex_colors=colors, info=True, camera_view=cam_view ) .. GENERATED FROM PYTHON SOURCE LINES 195-198 .. image:: ../_static/tutorial_07_orientation.png :width: 800 :alt: .. GENERATED FROM PYTHON SOURCE LINES 200-201 **Distance**: reflects the distance to the scalp .. GENERATED FROM PYTHON SOURCE LINES 201-224 .. code-block:: default c_range = [np.percentile(distances,1), np.percentile(distances,99)] # Plot change in power on each surface colors,_ = color_map( distances, "Spectral_r", c_range[0], c_range[1], norm='N' ) cam_view = [335, 9.5, 51, 60, 37, 17, 0, 0, 1] plot = show_surface( surf_set, vertex_colors=colors, info=True, camera_view=cam_view ) .. GENERATED FROM PYTHON SOURCE LINES 225-228 .. image:: ../_static/tutorial_07_distance.png :width: 800 :alt: .. GENERATED FROM PYTHON SOURCE LINES 230-233 Normalize and combine predictors -------------------------------- To assess joint contributions of these features, we standardize (z-score) each anatomical predictor and form a composite 'anatomical score': .. GENERATED FROM PYTHON SOURCE LINES 233-265 .. code-block:: default # Z-score normalization (invert distance since smaller = closer to sensors) z_thickness = zscore(thickness) z_lf_rms_diff = zscore(lf_rms_diff) z_orient = zscore(orientations) z_inv_dist = zscore(-distances) # Composite anatomical score anatomical_score = z_thickness + z_lf_rms_diff + z_orient + z_inv_dist c_range = [np.percentile(anatomical_score,1), np.percentile(anatomical_score,99)] # Plot change in power on each surface colors,_ = color_map( anatomical_score, "Spectral_r", c_range[0], c_range[1], norm='N' ) cam_view = [335, 9.5, 51, 60, 37, 17, 0, 0, 1] plot = show_surface( surf_set, vertex_colors=colors, info=True, camera_view=cam_view ) .. GENERATED FROM PYTHON SOURCE LINES 266-269 .. image:: ../_static/tutorial_07_anatomical_score.png :width: 800 :alt: .. GENERATED FROM PYTHON SOURCE LINES 271-272 We can find the best vertex over the whole brain for laminar inference. It happens to be the same one we've been using in the tutorial simulations. Isn't that a coincidence!? .. GENERATED FROM PYTHON SOURCE LINES 272-277 .. code-block:: default # Vertex with maximal predicted laminar discriminability best_vertex = np.argmax(anatomical_score) print(f"Highest anatomical score at vertex {best_vertex}") .. GENERATED FROM PYTHON SOURCE LINES 278-286 Interpretation -------------- Vertices with higher anatomical_score values are predicted to yield more accurate laminar inference. That is, lower reconstruction error when differentiating deep vs. superficial sources. In Szul et al., 2025, this multivariate anatomical score closely paralleled empirical reconstruction accuracy maps, highlighting regions where laminar separation is anatomically favorable. Selecting anatomical priors within a region of interest ------------------------------------------------------- In practical laminar inference, it is often desirable to constrain model comparison to a subset of vertices, e.g., within an anatomical ROI or around vertices showing maximal task-related power or signal amplitude. Within such a region, the anatomical score can guide the choice of the most anatomically suitable prior location for laminar inference: .. GENERATED FROM PYTHON SOURCE LINES 286-291 .. code-block:: default # Among the vertex in roi_idx, select the vertex with the best anatomical suitability candidate_anat_scores = anatomical_score[roi_idx] prior_idx = roi_idx[np.argmax(candidate_anat_scores)] .. GENERATED FROM PYTHON SOURCE LINES 292-297 .. code-block:: default spm.terminate() # Delete simulation files shutil.rmtree(tmp_dir) .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.000 seconds) .. _sphx_glr_download_auto_tutorials_tutorial_07_anatomy_inference_accuracy.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: tutorial_07_anatomy_inference_accuracy.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: tutorial_07_anatomy_inference_accuracy.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_