From 94f1356d962e6535cad8f6cc693c1dbcf52a8a47 Mon Sep 17 00:00:00 2001 From: Alexandre Gramfort Date: Sun, 4 Jul 2021 17:34:26 +0200 Subject: [PATCH 1/6] add plot_anat_landmarks function add plot_anat_landmarks function --- doc/Makefile | 13 ++++++ doc/api.rst | 17 +++++++- examples/convert_mri_and_trans.py | 38 ++++++++--------- mne_bids/__init__.py | 1 + mne_bids/tests/test_viz.py | 36 ++++++++++++++++ mne_bids/viz.py | 70 +++++++++++++++++++++++++++++++ 6 files changed, 155 insertions(+), 20 deletions(-) create mode 100644 mne_bids/tests/test_viz.py create mode 100644 mne_bids/viz.py diff --git a/doc/Makefile b/doc/Makefile index 2136bcd33..0e3384dee 100644 --- a/doc/Makefile +++ b/doc/Makefile @@ -9,6 +9,11 @@ PAPER = SOURCEDIR = . BUILDDIR = _build +# Internal variables. +PAPEROPT_a4 = -D latex_paper_size=a4 +PAPEROPT_letter = -D latex_paper_size=letter +ALLSPHINXOPTS = -d _build/doctrees $(PAPEROPT_$(PAPER)) $(SPHINXOPTS) . + # Put it first so that "make" without argument is like "make help". help: @$(SPHINXBUILD) -M help "$(SOURCEDIR)" "$(BUILDDIR)" $(SPHINXOPTS) $(O) @@ -23,3 +28,11 @@ help: # View the built documentation view: @python -c "import webbrowser; webbrowser.open_new_tab('file://$(PWD)/_build/html/index.html')" + +clean: + -rm -rf _build auto_examples generated + +html-noplot: + @$(SPHINXBUILD) -D plot_gallery=0 -b html $(ALLSPHINXOPTS) _build/html + @echo + @echo "Build finished. The HTML pages are in _build/html." diff --git a/doc/api.rst b/doc/api.rst index 762775a9a..ae2b4a8b2 100644 --- a/doc/api.rst +++ b/doc/api.rst @@ -53,7 +53,6 @@ mne_bids.stats count_events - mne_bids.copyfiles ------------------ @@ -74,3 +73,19 @@ mne_bids.copyfiles copyfile_ctf copyfile_bti copyfile_kit + +mne_bids.viz +-------------- + +:py:mod:`mne_bids.viz`: + +.. automodule:: mne_bids.viz + :no-members: + :no-inherited-members: + +.. currentmodule:: mne_bids.viz + +.. autosummary:: + :toctree: generated/ + + plot_anat_landmarks diff --git a/examples/convert_mri_and_trans.py b/examples/convert_mri_and_trans.py index bb1380b05..e53b17f75 100644 --- a/examples/convert_mri_and_trans.py +++ b/examples/convert_mri_and_trans.py @@ -51,6 +51,7 @@ from mne_bids import (write_raw_bids, BIDSPath, write_anat, get_head_mri_trans, print_dir_tree) +from mne_bids.viz import plot_anat_landmarks ############################################################################### # We will be using the `MNE sample data `_ and write a basic @@ -134,16 +135,20 @@ anat_dir = t1w_bids_path.directory ############################################################################### -# Let's have another look at our BIDS directory +# Let's have another look at our BIDS directory and plot the written landmarks print_dir_tree(output_path) +plot_anat_landmarks(t1w_bids_path, vmax=160) +plt.suptitle('T1 MRI') + ############################################################################### # Our BIDS dataset is now ready to be shared. We can easily estimate the -# transformation matrix using ``MNE-BIDS`` and the BIDS dataset. +# transformation matrix using ``MNE-BIDS`` and the BIDS dataset since we +# have now anatomical landmarks. estim_trans = get_head_mri_trans(bids_path=bids_path) ############################################################################### -# Finally, let's use the T1 weighted MRI image and plot the anatomical +# One can also directly use the T1 weighted MRI image and plot the anatomical # landmarks Nasion, LPA, and RPA onto the brain image. For that, we can # extract the location of Nasion, LPA, and RPA from the MEG file, apply our # transformation matrix :code:`trans`, and plot the results. @@ -169,7 +174,7 @@ t1_nii_fname = op.join(anat_dir, 'sub-01_ses-01_T1w.nii.gz') # Plot it -fig, axs = plt.subplots(3, 1, figsize=(7, 7), facecolor='k') +fig, axs = plt.subplots(3, 1, figsize=(7, 7)) for point_idx, label in enumerate(('LPA', 'NAS', 'RPA')): plot_anat(t1_nii_fname, axes=axs[point_idx], cut_coords=mri_pos[point_idx, :], @@ -214,9 +219,8 @@ t1_nii_fname = op.join(anat_dir, 'sub-01_ses-01_T1w.nii.gz') # Plot it -fig, ax = plt.subplots() -plot_anat(t1_nii_fname, axes=ax, title='Defaced', vmax=160) -plt.show() +plot_anat_landmarks(t1w_bids_path, vmax=160) +plt.suptitle('T1 Defaced') ############################################################################### # Writing defaced and anonymized FLASH MRI image @@ -250,9 +254,8 @@ flash_nii_fname = op.join(anat_dir, 'sub-01_ses-01_FLASH.nii.gz') # Plot it -fig, ax = plt.subplots() -plot_anat(flash_nii_fname, axes=ax, title='Defaced', vmax=700) -plt.show() +plot_anat_landmarks(flash_bids_path, vmax=700) +plt.suptitle('Flash Defaced') ############################################################################### # Option 2 : Use manual landmarks coordinates in scanner RAS for FLASH image @@ -282,9 +285,8 @@ ) # Plot it -fig, ax = plt.subplots() -plot_anat(flash_nii_fname, axes=ax, title='Defaced', vmax=700) -plt.show() +plot_anat_landmarks(flash_bids_path, vmax=700) +plt.suptitle('Flash Defaced') ############################################################################### # Option 3 : Compute the landmarks in scanner RAS or mri voxel space from trans @@ -324,9 +326,8 @@ ) # Plot it -fig, ax = plt.subplots() -plot_anat(flash_nii_fname, axes=ax, title='Defaced', vmax=700) -plt.show() +plot_anat_landmarks(flash_bids_path, vmax=700) +plt.suptitle('Flash Defaced') ############################################################################## # Let's now pass it in voxel coordinates @@ -352,9 +353,8 @@ ) # Plot it -fig, ax = plt.subplots() -plot_anat(flash_nii_fname, axes=ax, title='Defaced', vmax=700) -plt.show() +plot_anat_landmarks(flash_bids_path, vmax=700) +plt.suptitle('Flash Defaced') ############################################################################### # .. LINKS diff --git a/mne_bids/__init__.py b/mne_bids/__init__.py index 3897cb480..706ccc767 100644 --- a/mne_bids/__init__.py +++ b/mne_bids/__init__.py @@ -12,3 +12,4 @@ write_meg_calibration, write_meg_crosstalk) from mne_bids.sidecar_updates import update_sidecar_json from mne_bids.inspect import inspect_dataset +from mne_bids import viz diff --git a/mne_bids/tests/test_viz.py b/mne_bids/tests/test_viz.py new file mode 100644 index 000000000..546bf2c1e --- /dev/null +++ b/mne_bids/tests/test_viz.py @@ -0,0 +1,36 @@ +from pathlib import Path + +import pytest + +import mne +from mne.datasets import testing +from mne_bids.viz import plot_anat_landmarks +from mne_bids import BIDSPath, write_anat + +from mne.utils import requires_nibabel, requires_version + + +@requires_nibabel() +@requires_version('nilearn', '0.6') +def test_plot_anat_landmarks(tmpdir): + """Test writing anatomical data with pathlib.Paths.""" + data_path = Path(testing.data_path()) + raw_fname = data_path / 'MEG' / 'sample' / 'sample_audvis_trunc_raw.fif' + trans_fname = str(raw_fname).replace('_raw.fif', '-trans.fif') + raw = mne.io.read_raw_fif(raw_fname) + trans = mne.read_trans(trans_fname) + + bids_root = Path(tmpdir) + t1w_mgh_fname = Path(data_path) / 'subjects' / 'sample' / 'mri' / 'T1.mgz' + bids_path = BIDSPath(subject='01', session='mri', root=bids_root) + bids_path = write_anat(t1w_mgh_fname, bids_path=bids_path, overwrite=True) + + with pytest.warns(RuntimeWarning, match='No landmarks available'): + fig = plot_anat_landmarks(bids_path, show=False) + assert fig is None + + bids_path = write_anat(t1w_mgh_fname, bids_path=bids_path, raw=raw, + trans=trans, deface=True, verbose=True, + overwrite=True) + fig = plot_anat_landmarks(bids_path, show=False) + assert len(fig.axes) == 12 # 3 subplots + 3 x 3 MRI slices diff --git a/mne_bids/viz.py b/mne_bids/viz.py new file mode 100644 index 000000000..7a551bcbc --- /dev/null +++ b/mne_bids/viz.py @@ -0,0 +1,70 @@ +import json + +import numpy as np +from scipy import linalg + +import nibabel as nib +from mne.utils import warn +import mne + + +def plot_anat_landmarks(bids_path, vmax=None, show=True): + """Plot anatomical landmarks attached to an MRI image. + + Parameters + ---------- + bids_path : mne_bids.BIDSPath + Path of the MRI image. + vmax : float + Maximum colormap value. + show : bool + Call pyplot.show() at the end. Defaults to True. + + Returns + ------- + fig : matplotlib.figure.Figure | None + The figure object containing the plot. None if no landmarks + are avaiable. + """ + import matplotlib.pyplot as plt + from nilearn.plotting import plot_anat + + nii = nib.load(str(bids_path)) + + json_path = bids_path.copy().update(extension=".json") + + n_landmarks = 0 + if json_path.fpath.exists(): + json_content = json.load(open(json_path)) + coords_dict = json_content.get("AnatomicalLandmarkCoordinates", dict()) + n_landmarks = len(coords_dict) + + if not n_landmarks: + warn("No landmarks available with the image") + return + + # Move the coords_dict from MRI Voxel to RAS + mgh = nib.MGHImage(nii.dataobj, nii.affine) + ras2vox = mgh.header.get_ras2vox() + vox2ras = linalg.inv(ras2vox) + for label in coords_dict: + vox_pos = np.array(coords_dict[label]) + ras_pos = mne.transforms.apply_trans(vox2ras, vox_pos) + coords_dict[label] = ras_pos + + ######################################################################## + # Plot it with nilearn + fig, axs = plt.subplots( + n_landmarks, 1, figsize=(6, 2.3 * n_landmarks), + facecolor="w") + + for point_idx, (label, ras_pos) in enumerate(coords_dict.items()): + plot_anat( + str(bids_path), axes=axs[point_idx], cut_coords=ras_pos, + title=label, vmax=vmax, + ) + + if show: + plt.show() + + return fig From bc728d9b797d8972c58dfa89934bdd74c3207b36 Mon Sep 17 00:00:00 2001 From: Alexandre Gramfort Date: Thu, 8 Jul 2021 14:24:13 +0200 Subject: [PATCH 2/6] Apply suggestions from code review MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Richard Höchenberger --- examples/convert_mri_and_trans.py | 2 +- mne_bids/tests/test_viz.py | 3 +-- mne_bids/viz.py | 6 +++--- 3 files changed, 5 insertions(+), 6 deletions(-) diff --git a/examples/convert_mri_and_trans.py b/examples/convert_mri_and_trans.py index e53b17f75..a3563287f 100644 --- a/examples/convert_mri_and_trans.py +++ b/examples/convert_mri_and_trans.py @@ -148,7 +148,7 @@ estim_trans = get_head_mri_trans(bids_path=bids_path) ############################################################################### -# One can also directly use the T1 weighted MRI image and plot the anatomical +# One can also directly use the T1-weighted MR image and plot the anatomical # landmarks Nasion, LPA, and RPA onto the brain image. For that, we can # extract the location of Nasion, LPA, and RPA from the MEG file, apply our # transformation matrix :code:`trans`, and plot the results. diff --git a/mne_bids/tests/test_viz.py b/mne_bids/tests/test_viz.py index 546bf2c1e..c0e8430e4 100644 --- a/mne_bids/tests/test_viz.py +++ b/mne_bids/tests/test_viz.py @@ -4,11 +4,10 @@ import mne from mne.datasets import testing +from mne.utils import requires_nibabel, requires_version from mne_bids.viz import plot_anat_landmarks from mne_bids import BIDSPath, write_anat -from mne.utils import requires_nibabel, requires_version - @requires_nibabel() @requires_version('nilearn', '0.6') diff --git a/mne_bids/viz.py b/mne_bids/viz.py index 7a551bcbc..ef500cfc3 100644 --- a/mne_bids/viz.py +++ b/mne_bids/viz.py @@ -2,10 +2,10 @@ import numpy as np from scipy import linalg - import nibabel as nib -from mne.utils import warn + import mne +from mne.utils import warn def plot_anat_landmarks(bids_path, vmax=None, show=True): @@ -18,7 +18,7 @@ def plot_anat_landmarks(bids_path, vmax=None, show=True): vmax : float Maximum colormap value. show : bool - Call pyplot.show() at the end. Defaults to True. + Whether to show the figure after plotting. Defaults to ``True``. Returns ------- From 8ebeadbf846607f0955a95408eaf1838a70af811 Mon Sep 17 00:00:00 2001 From: Alexandre Gramfort Date: Thu, 8 Jul 2021 14:24:39 +0200 Subject: [PATCH 3/6] address some comments --- examples/convert_mri_and_trans.py | 7 ++----- mne_bids/viz.py | 12 ++++++------ 2 files changed, 8 insertions(+), 11 deletions(-) diff --git a/examples/convert_mri_and_trans.py b/examples/convert_mri_and_trans.py index 7589d986c..484b4e7e0 100644 --- a/examples/convert_mri_and_trans.py +++ b/examples/convert_mri_and_trans.py @@ -147,14 +147,14 @@ # Our BIDS dataset is now ready to be shared. We can easily estimate the # transformation matrix using ``MNE-BIDS`` and the BIDS dataset since we # have now anatomical landmarks. -estim_trans = get_head_mri_trans(bids_path=bids_path) +estim_trans = get_head_mri_trans( + bids_path=bids_path, fs_subject='sample', fs_subjects_dir=fs_subjects_dir) ############################################################################### # Let's have another look at our BIDS directory and plot the written landmarks print_dir_tree(output_path) plot_anat_landmarks(t1w_bids_path, vmax=160) -plt.suptitle('T1 MRI') # %% # Our BIDS dataset is now ready to be shared. We can easily estimate the @@ -242,7 +242,6 @@ # Plot it plot_anat_landmarks(t1w_bids_path, vmax=160) -plt.suptitle('T1 Defaced') # %% # Writing defaced and anonymized FLASH MRI image @@ -274,7 +273,6 @@ # Plot it plot_anat_landmarks(flash_bids_path, vmax=700) -plt.suptitle('Flash Defaced') # %% # Using manual landmark coordinates in scanner RAS @@ -310,7 +308,6 @@ # Plot it plot_anat_landmarks(flash_bids_path, vmax=700) -plt.suptitle('Flash Defaced') # %% # .. LINKS diff --git a/mne_bids/viz.py b/mne_bids/viz.py index 7a551bcbc..95ab2a043 100644 --- a/mne_bids/viz.py +++ b/mne_bids/viz.py @@ -22,9 +22,8 @@ def plot_anat_landmarks(bids_path, vmax=None, show=True): Returns ------- - fig : matplotlib.figure.Figure | None - The figure object containing the plot. None if no landmarks - are avaiable. + fig : matplotlib.figure.Figure + The figure object containing the plot. """ import matplotlib.pyplot as plt from nilearn.plotting import plot_anat @@ -40,8 +39,7 @@ def plot_anat_landmarks(bids_path, vmax=None, show=True): n_landmarks = len(coords_dict) if not n_landmarks: - warn("No landmarks available with the image") - return + raise ValueError("No landmarks available with the image") # Move the coords_dict from MRI Voxel to RAS mgh = nib.MGHImage(nii.dataobj, nii.affine) @@ -64,7 +62,9 @@ def plot_anat_landmarks(bids_path, vmax=None, show=True): title=label, vmax=vmax, ) + plt.suptitle(bids_path.fpath.name) + if show: - plt.show() + fig.show() return fig From ca17020354630c7b13174c9e9f6d5e26f92eeb68 Mon Sep 17 00:00:00 2001 From: Alexandre Gramfort Date: Thu, 8 Jul 2021 14:28:33 +0200 Subject: [PATCH 4/6] nitpicks --- examples/convert_mri_and_trans.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/examples/convert_mri_and_trans.py b/examples/convert_mri_and_trans.py index 484b4e7e0..fd14dbf90 100644 --- a/examples/convert_mri_and_trans.py +++ b/examples/convert_mri_and_trans.py @@ -143,14 +143,14 @@ ) anat_dir = t1w_bids_path.directory -############################################################################### +# %% # Our BIDS dataset is now ready to be shared. We can easily estimate the # transformation matrix using ``MNE-BIDS`` and the BIDS dataset since we # have now anatomical landmarks. estim_trans = get_head_mri_trans( bids_path=bids_path, fs_subject='sample', fs_subjects_dir=fs_subjects_dir) -############################################################################### +# %% # Let's have another look at our BIDS directory and plot the written landmarks print_dir_tree(output_path) From 3e0deb035f2cf556ceb31565681f60fee00da5f1 Mon Sep 17 00:00:00 2001 From: Alexandre Gramfort Date: Thu, 8 Jul 2021 14:40:36 +0200 Subject: [PATCH 5/6] simplify + fix test --- mne_bids/tests/test_viz.py | 19 +++++++++++-------- mne_bids/viz.py | 7 +------ 2 files changed, 12 insertions(+), 14 deletions(-) diff --git a/mne_bids/tests/test_viz.py b/mne_bids/tests/test_viz.py index c0e8430e4..93529168c 100644 --- a/mne_bids/tests/test_viz.py +++ b/mne_bids/tests/test_viz.py @@ -6,7 +6,7 @@ from mne.datasets import testing from mne.utils import requires_nibabel, requires_version from mne_bids.viz import plot_anat_landmarks -from mne_bids import BIDSPath, write_anat +from mne_bids import BIDSPath, write_anat, get_anat_landmarks @requires_nibabel() @@ -18,18 +18,21 @@ def test_plot_anat_landmarks(tmpdir): trans_fname = str(raw_fname).replace('_raw.fif', '-trans.fif') raw = mne.io.read_raw_fif(raw_fname) trans = mne.read_trans(trans_fname) + fs_subjects_dir = Path(data_path) / 'subjects' bids_root = Path(tmpdir) - t1w_mgh_fname = Path(data_path) / 'subjects' / 'sample' / 'mri' / 'T1.mgz' + t1w_mgh_fname = fs_subjects_dir / 'sample' / 'mri' / 'T1.mgz' bids_path = BIDSPath(subject='01', session='mri', root=bids_root) bids_path = write_anat(t1w_mgh_fname, bids_path=bids_path, overwrite=True) - with pytest.warns(RuntimeWarning, match='No landmarks available'): - fig = plot_anat_landmarks(bids_path, show=False) - assert fig is None + with pytest.raises(ValueError, match='No landmarks available'): + plot_anat_landmarks(bids_path, show=False) - bids_path = write_anat(t1w_mgh_fname, bids_path=bids_path, raw=raw, - trans=trans, deface=True, verbose=True, - overwrite=True) + landmarks = get_anat_landmarks( + t1w_mgh_fname, raw.info, trans, fs_subject='sample', + fs_subjects_dir=fs_subjects_dir) + + bids_path = write_anat(t1w_mgh_fname, bids_path=bids_path, + landmarks=landmarks, overwrite=True) fig = plot_anat_landmarks(bids_path, show=False) assert len(fig.axes) == 12 # 3 subplots + 3 x 3 MRI slices diff --git a/mne_bids/viz.py b/mne_bids/viz.py index b5d6f536a..cb5a82cd5 100644 --- a/mne_bids/viz.py +++ b/mne_bids/viz.py @@ -5,7 +5,6 @@ import nibabel as nib import mne -from mne.utils import warn def plot_anat_landmarks(bids_path, vmax=None, show=True): @@ -41,13 +40,9 @@ def plot_anat_landmarks(bids_path, vmax=None, show=True): if not n_landmarks: raise ValueError("No landmarks available with the image") - # Move the coords_dict from MRI Voxel to RAS - mgh = nib.MGHImage(nii.dataobj, nii.affine) - ras2vox = mgh.header.get_ras2vox() - vox2ras = linalg.inv(ras2vox) for label in coords_dict: vox_pos = np.array(coords_dict[label]) - ras_pos = mne.transforms.apply_trans(vox2ras, vox_pos) + ras_pos = mne.transforms.apply_trans(nii.affine, vox_pos) coords_dict[label] = ras_pos ######################################################################## From 2e7ee25f90763cbff764d495215edfb245ebc9fc Mon Sep 17 00:00:00 2001 From: Alexandre Gramfort Date: Thu, 8 Jul 2021 15:51:58 +0200 Subject: [PATCH 6/6] flake8 --- mne_bids/viz.py | 1 - 1 file changed, 1 deletion(-) diff --git a/mne_bids/viz.py b/mne_bids/viz.py index cb5a82cd5..fb1fb4855 100644 --- a/mne_bids/viz.py +++ b/mne_bids/viz.py @@ -1,7 +1,6 @@ import json import numpy as np -from scipy import linalg import nibabel as nib import mne