1"""Example of explicit fixed effects fMRI model fitting 2==================================================== 3 4This example illustrates how to run a fixed effects model based on 5pre-computed statistics. This is helpful when the initial models 6have to be fit separately. 7 8For details on the data, please see: 9 10Dehaene-Lambertz G, Dehaene S, Anton JL, Campagne A, Ciuciu P, Dehaene 11G, Denghien I, Jobert A, LeBihan D, Sigman M, Pallier C, Poline 12JB. Functional segregation of cortical language areas by sentence 13repetition. Hum Brain Mapp. 2006: 27:360--371. 14http://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=2653076#R11 15 16Please see `Simple example of two-session :term:`fMRI` model fitting 17<https://nistats.github.io/auto_examples/02_first_levels/plot_fiac_analysis.html>`_ 18example for details. The main difference is that 19the fixed-effects model is run explicitly here, 20after GLM fitting on two sessions. 21 22""" 23 24######################################################################### 25# Prepare data and analysis parameters 26# -------------------------------------- 27# 28# Inspecting 'data', we note that there are two sessions 29 30from nilearn.datasets import func 31data = func.fetch_fiac_first_level() 32fmri_img = [data['func1'], data['func2']] 33 34######################################################################### 35# Create a mean image for plotting purpose 36from nilearn.image import mean_img 37mean_img_ = mean_img(fmri_img[0]) 38 39######################################################################### 40# The design matrices were pre-computed, we simply put them in a list of 41# DataFrames 42design_files = [data['design_matrix1'], data['design_matrix2']] 43import pandas as pd 44import numpy as np 45design_matrices = [pd.DataFrame(np.load(df)['X']) for df in design_files] 46 47######################################################################### 48# GLM estimation 49# ---------------------------------- 50# GLM specification. Note that the mask was provided in the dataset. 51# So we use it. 52 53from nilearn.glm.first_level import FirstLevelModel 54fmri_glm = FirstLevelModel(mask_img=data['mask'], smoothing_fwhm=5, 55 minimize_memory=True) 56 57######################################################################### 58# Compute fixed effects of the two runs and compute related images 59# For this, we first define the contrasts as we would do for a single session 60n_columns = design_matrices[0].shape[1] 61contrast_val = np.hstack(([-1, -1, 1, 1], np.zeros(n_columns - 4))) 62 63######################################################################### 64# Statistics for the first session 65from nilearn import plotting 66cut_coords = [-129, -126, 49] 67contrast_id = 'DSt_minus_SSt' 68 69fmri_glm = fmri_glm.fit(fmri_img[0], design_matrices=design_matrices[0]) 70summary_statistics_session1 = fmri_glm.compute_contrast( 71 contrast_val, output_type='all') 72plotting.plot_stat_map( 73 summary_statistics_session1['z_score'], 74 bg_img=mean_img_, threshold=3.0, cut_coords=cut_coords, 75 title='{0}, first session'.format(contrast_id)) 76 77######################################################################### 78# Statistics for the second session 79 80fmri_glm = fmri_glm.fit(fmri_img[1], design_matrices=design_matrices[1]) 81summary_statistics_session2 = fmri_glm.compute_contrast( 82 contrast_val, output_type='all') 83plotting.plot_stat_map( 84 summary_statistics_session2['z_score'], 85 bg_img=mean_img_, threshold=3.0, cut_coords=cut_coords, 86 title='{0}, second session'.format(contrast_id)) 87 88######################################################################### 89# Fixed effects statistics 90from nilearn.glm.contrasts import compute_fixed_effects 91 92contrast_imgs = [summary_statistics_session1['effect_size'], 93 summary_statistics_session2['effect_size']] 94variance_imgs = [summary_statistics_session1['effect_variance'], 95 summary_statistics_session2['effect_variance']] 96 97fixed_fx_contrast, fixed_fx_variance, fixed_fx_stat = compute_fixed_effects( 98 contrast_imgs, variance_imgs, data['mask']) 99plotting.plot_stat_map( 100 fixed_fx_stat, bg_img=mean_img_, threshold=3.0, cut_coords=cut_coords, 101 title='{0}, fixed effects'.format(contrast_id)) 102 103######################################################################### 104# Not unexpectedly, the fixed effects version displays higher peaks than the 105# input sessions. Computing fixed effects enhances the signal-to-noise ratio of 106# the resulting brain maps 107# Note however that, technically, the output maps of the fixed effects map is a 108# t statistic (not a z statistic) 109 110plotting.show() 111