.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/05_glm_second_level/plot_oasis.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 or to run this example in your browser via Binder .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_05_glm_second_level_plot_oasis.py: Voxel-Based Morphometry on OASIS dataset ======================================== This example uses voxel-based morphometry (:term:`VBM`) to study the relationship between aging, sex, and gray matter density. The data come from the `OASIS `_ project. If you use it, you need to agree with the data usage agreement available on the website. It has been run through a standard :term:`VBM` pipeline (using SPM8 and NewSegment) to create :term:`VBM` maps, which we study here. VBM analysis of aging --------------------- We run a standard :term:`GLM` analysis to study the association between age and gray matter density from the :term:`VBM` data. We use only 100 subjects from the OASIS dataset to limit the memory usage. Note that more power would be obtained from using a larger sample of subjects. .. Original authors: - Bertrand Thirion, , July 2018 - Elvis Dhomatob, Apr. 2014 - Virgile Fritsch, Apr 2014 - Gael Varoquaux, Apr 2014 .. GENERATED FROM PYTHON SOURCE LINES 37-39 Load Oasis dataset ------------------ .. GENERATED FROM PYTHON SOURCE LINES 39-50 .. code-block:: Python from nilearn import datasets, plotting n_subjects = 100 # more subjects requires more memory oasis_dataset = datasets.fetch_oasis_vbm( n_subjects=n_subjects, legacy_format=False, ) gray_matter_map_filenames = oasis_dataset.gray_matter_maps age = oasis_dataset.ext_vars["age"].astype(float) .. GENERATED FROM PYTHON SOURCE LINES 51-52 Sex is encoded as 'M' or 'F'. Hence, we make it a binary variable. .. GENERATED FROM PYTHON SOURCE LINES 52-54 .. code-block:: Python sex = oasis_dataset.ext_vars["mf"] == "F" .. GENERATED FROM PYTHON SOURCE LINES 55-56 Print basic information on the dataset. .. GENERATED FROM PYTHON SOURCE LINES 56-65 .. code-block:: Python print( "First gray-matter anatomy image (3D) is located at: " f"{oasis_dataset.gray_matter_maps[0]}" ) print( "First white-matter anatomy image (3D) is located at: " f"{oasis_dataset.white_matter_maps[0]}" ) .. rst-class:: sphx-glr-script-out .. code-block:: none First gray-matter anatomy image (3D) is located at: /home/remi/nilearn_data/oasis1/OAS1_0001_MR1/mwrc1OAS1_0001_MR1_mpr_anon_fslswapdim_bet.nii.gz First white-matter anatomy image (3D) is located at: /home/remi/nilearn_data/oasis1/OAS1_0001_MR1/mwrc2OAS1_0001_MR1_mpr_anon_fslswapdim_bet.nii.gz .. GENERATED FROM PYTHON SOURCE LINES 66-67 Get a mask image: A mask of the cortex of the ICBM template. .. GENERATED FROM PYTHON SOURCE LINES 67-69 .. code-block:: Python gm_mask = datasets.fetch_icbm152_brain_gm_mask() .. GENERATED FROM PYTHON SOURCE LINES 70-71 Resample the mask, since this mask has a different resolution. .. GENERATED FROM PYTHON SOURCE LINES 71-79 .. code-block:: Python from nilearn.image import resample_to_img mask_img = resample_to_img( gm_mask, gray_matter_map_filenames[0], interpolation="nearest", ) .. GENERATED FROM PYTHON SOURCE LINES 80-84 Analyse data ------------ First, we create an adequate design matrix with three columns: 'age', 'sex', and 'intercept'. .. GENERATED FROM PYTHON SOURCE LINES 84-95 .. code-block:: Python import numpy as np import pandas as pd intercept = np.ones(n_subjects) design_matrix = pd.DataFrame( np.vstack((age, sex, intercept)).T, columns=["age", "sex", "intercept"], ) from matplotlib import pyplot as plt .. GENERATED FROM PYTHON SOURCE LINES 96-97 Let's plot the design matrix. .. GENERATED FROM PYTHON SOURCE LINES 97-102 .. code-block:: Python fig, ax1 = plt.subplots(1, 1, figsize=(4, 8)) ax = plotting.plot_design_matrix(design_matrix, ax=ax1) ax.set_ylabel("maps") fig.suptitle("Second level design matrix") .. image-sg:: /auto_examples/05_glm_second_level/images/sphx_glr_plot_oasis_001.png :alt: Second level design matrix :srcset: /auto_examples/05_glm_second_level/images/sphx_glr_plot_oasis_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Text(0.5, 0.98, 'Second level design matrix') .. GENERATED FROM PYTHON SOURCE LINES 103-105 Next, we specify and fit the second-level model when loading the data and also smooth a little bit to improve statistical behavior. .. GENERATED FROM PYTHON SOURCE LINES 105-115 .. code-block:: Python from nilearn.glm.second_level import SecondLevelModel second_level_model = SecondLevelModel( smoothing_fwhm=2.0, mask_img=mask_img, n_jobs=2 ) second_level_model.fit( gray_matter_map_filenames, design_matrix=design_matrix, ) .. raw:: html
SecondLevelModel(mask_img=<nibabel.nifti1.Nifti1Image object at 0x7fb3c6b23110>,
                     n_jobs=2, smoothing_fwhm=2.0)
In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.


.. GENERATED FROM PYTHON SOURCE LINES 116-118 Estimating the :term:`contrast` is very simple. We can just provide the column name of the design matrix. .. GENERATED FROM PYTHON SOURCE LINES 118-123 .. code-block:: Python z_map = second_level_model.compute_contrast( second_level_contrast=[1, 0, 0], output_type="z_score", ) .. GENERATED FROM PYTHON SOURCE LINES 124-126 We threshold the second level :term:`contrast` at FDR-corrected p < 0.05 and plot it. .. GENERATED FROM PYTHON SOURCE LINES 126-143 .. code-block:: Python from nilearn.glm import threshold_stats_img _, threshold = threshold_stats_img(z_map, alpha=0.05, height_control="fdr") print(f"The FDR=.05-corrected threshold is: {threshold:03g}") fig = plt.figure(figsize=(5, 3)) display = plotting.plot_stat_map( z_map, threshold=threshold, colorbar=True, display_mode="z", cut_coords=[-4, 26], figure=fig, ) fig.suptitle("age effect on gray matter density (FDR = .05)") plotting.show() .. image-sg:: /auto_examples/05_glm_second_level/images/sphx_glr_plot_oasis_002.png :alt: age effect on gray matter density (FDR = .05) :srcset: /auto_examples/05_glm_second_level/images/sphx_glr_plot_oasis_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none The FDR=.05-corrected threshold is: 2.40175 .. GENERATED FROM PYTHON SOURCE LINES 144-146 We can also study the effect of sex by computing the contrast, thresholding it and plot the resulting map. .. GENERATED FROM PYTHON SOURCE LINES 146-158 .. code-block:: Python z_map = second_level_model.compute_contrast( second_level_contrast="sex", output_type="z_score", ) _, threshold = threshold_stats_img(z_map, alpha=0.05, height_control="fdr") plotting.plot_stat_map( z_map, threshold=threshold, colorbar=True, title="sex effect on gray matter density (FDR = .05)", ) .. image-sg:: /auto_examples/05_glm_second_level/images/sphx_glr_plot_oasis_003.png :alt: plot oasis :srcset: /auto_examples/05_glm_second_level/images/sphx_glr_plot_oasis_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 159-161 Note that there does not seem to be any significant effect of sex on gray matter density on that dataset. .. GENERATED FROM PYTHON SOURCE LINES 163-169 Generating a report ------------------- It can be useful to quickly generate a portable, ready-to-view report with most of the pertinent information. This is easy to do if you have a fitted model and the list of contrasts, which we do here. .. GENERATED FROM PYTHON SOURCE LINES 169-178 .. code-block:: Python from nilearn.reporting import make_glm_report icbm152_2009 = datasets.fetch_icbm152_2009() report = make_glm_report( model=second_level_model, contrasts=["age", "sex"], bg_img=icbm152_2009["t1"], ) .. GENERATED FROM PYTHON SOURCE LINES 179-180 We have several ways to access the report: .. GENERATED FROM PYTHON SOURCE LINES 180-189 .. code-block:: Python # report # This report can be viewed in a notebook # report.open_in_browser() # or we can save as an html file # from pathlib import Path # output_dir = Path.cwd() / "results" / "plot_oasis" # output_dir.mkdir(exist_ok=True, parents=True) # report.save_as_html(output_dir / 'report.html') .. rst-class:: sphx-glr-timing **Total running time of the script:** (2 minutes 38.932 seconds) **Estimated memory usage:** 1517 MB .. _sphx_glr_download_auto_examples_05_glm_second_level_plot_oasis.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: binder-badge .. image:: images/binder_badge_logo.svg :target: https://mybinder.org/v2/gh/nilearn/nilearn/0.10.3?urlpath=lab/tree/notebooks/auto_examples/05_glm_second_level/plot_oasis.ipynb :alt: Launch binder :width: 150 px .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_oasis.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_oasis.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_