<jupyter_start><jupyter_text>BIDS dataset first and second level analysis============================================by Nilearn & Nistats Teamsmodified/added by <NAME>Full step-by-step example of fitting a GLM to perform a first and second levelanalysis in a BIDS dataset and visualizing the results. Details about the BIDSstandard can be consulted at http://bids.neuroimaging.io/More specifically:1. Download an fMRI BIDS dataset with two language conditions to contrast.2. Extract automatically from the BIDS dataset first level model objects3. Fit a second level model on the fitted first level models. Notice that   in this case the preprocessed bold images were already normalized to the   same MNI space.<jupyter_code>import warnings
warnings.filterwarnings('ignore')

import matplotlib.pyplot as plt

%matplotlib inline<jupyter_output><empty_output><jupyter_text>Fetch example BIDS dataset--------------------------We download an simplified BIDS dataset made available for illustrativepurposes. It contains only the necessaryinformation to run a statistical analysis using Nistats. The raw datasubject folders only contain bold.json and events.tsv files, while thederivatives folder with preprocessed files contain preproc.nii andconfounds.tsv files.<jupyter_code>from nistats.datasets import fetch_language_localizer_demo_dataset
data_dir, _ = fetch_language_localizer_demo_dataset()
data_dir<jupyter_output><empty_output><jupyter_text>Obtain automatically FirstLevelModel objects and fit arguments--------------------------------------------------------------From the dataset directory we obtain automatically FirstLevelModel objectswith their subject_id filled from the BIDS dataset. Moreover we obtainfor each model a dictionary with run_imgs, events and confounder regressorssince in this case a confounds.tsv file is available in the BIDS dataset.To get the first level models we only have to specify the dataset directoryand the task_label as specified in the file names.<jupyter_code>from nistats.first_level_model import first_level_models_from_bids
task_label = 'languagelocalizer'
models, models_run_imgs, models_events, models_confounds = \
    first_level_models_from_bids(
        data_dir, task_label,
        img_filters=[('desc', 'preproc')])<jupyter_output><empty_output><jupyter_text>Quick sanity check on fit arguments-----------------------------------Additional checks or information extraction from pre-processed data canbe made here We just expect one run img per subject.<jupyter_code>import os
print([os.path.basename(run) for run in models_run_imgs[0]])<jupyter_output>['sub-01_task-languagelocalizer_desc-preproc_bold.nii.gz']<jupyter_text>The only confounds stored are regressors obtained from motion correction. Aswe can verify from the column headers of the confounds table correspondingto the only run_img present<jupyter_code>print(models_confounds[0][0].columns)
models_confounds[2][0].values.shape
%matplotlib inline
plt.plot(models_confounds[2][0].values[:,0:3])<jupyter_output><empty_output><jupyter_text>***Exercise:***What is this? What does it mean?Can you modify the code above to check what happened with this subject in the other confounds, and perhaps even in other subjects? If you succeed in doing so, you will learn about the datastructure. During this acquisition the subject read blocks of sentences andconsonant strings. So these are our only two conditions in events.We verify there are 12 blocks for each condition.<jupyter_code>print(models_events[0][0]['trial_type'].value_counts())<jupyter_output>language    12
string      12
Name: trial_type, dtype: int64<jupyter_text>***Exercise:***What again was a block? First level model estimation----------------------------Now we simply fit each first level model and plot for each subject thecontrast that reveals the language network (language - string). Notice thatwe can define a contrast using the names of the conditions specified in theevents dataframe. Sum, substraction and scalar multiplication are allowed. Set the threshold as the z-variate with an uncorrected p-value of 0.001<jupyter_code>from scipy.stats import norm
p001_unc = norm.isf(0.001)<jupyter_output><empty_output><jupyter_text>***Exercise:***What does norm.isf do? Prepare figure for concurrent plot of individual maps (might take a while)<jupyter_code>from nilearn import plotting
import matplotlib.pyplot as plt

fig, axes = plt.subplots(nrows=2, ncols=5, figsize=(8, 4.5))
model_and_args = zip(models, models_run_imgs, models_events, models_confounds)
for midx, (model, imgs, events, confounds) in enumerate(model_and_args):
    # fit the GLM
    model.fit(imgs, events, confounds)
    # compute the contrast of interest
    zmap = model.compute_contrast('language-string')
    plotting.plot_glass_brain(zmap, colorbar=False, threshold=p001_unc,
                              title=('sub-' + model.subject_label),
                              axes=axes[int(midx / 5), int(midx % 5)],
                              plot_abs=False, display_mode='x')
fig.suptitle('subjects z_map language network (unc p<0.001)')
plotting.show()<jupyter_output><empty_output><jupyter_text>Second level model estimation-----------------------------We just have to provide the list of fitted FirstLevelModel objects to the SecondLevelModel object for estimation. We can do this because all subjects share a similar design matrix (same variables reflected in column names)<jupyter_code>from nistats.second_level_model import SecondLevelModel
second_level_input = models
second_level_input<jupyter_output><empty_output><jupyter_text>***Exercise:***How many subjects again? Can you explain to your neighbor what will happen next?  Fitting the 2nd Level ModelNow we fit the model (drumroll)Note that we apply a smoothing of 8mm.<jupyter_code>second_level_model = SecondLevelModel(smoothing_fwhm=8.0)
second_level_model = second_level_model.fit(second_level_input)<jupyter_output><empty_output><jupyter_text>Computing contrasts at the second level is as simple as at the first levelSince we are not providing confounders we are performing an one-sample testat the second level with the images determined by the specified first levelcontrast.<jupyter_code>zmap = second_level_model.compute_contrast(
    first_level_contrast='language-string')<jupyter_output><empty_output><jupyter_text>The group level contrast reveals a left lateralized fronto-temporallanguage network<jupyter_code>p001_unc
plotting.view_img(zmap, colorbar=True, threshold=3,
                          title='Group language network (unc p<0.0001)',
                          plot_abs=False, display_mode='x')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.09740357196458052, 'fraction_numerical': 0.009305117814798139, 'mean_word_length': 4.91829484902309, 'pattern_counts': {'":': 0, '<': 50, '<?xml version=': 0, '>': 48, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '12', 'dir': 'jupyter-structured-clean-dedup', 'id': '97902', 'license': 'Apache-2.0', 'n_tokens_mistral': 1872, 'n_tokens_neox': 1798, 'n_words': 637, 'path': 'scripts/02_BIDS_1st_2nd_Level_Analysis.ipynb', 'repo_name': 'smgrady-git/medianeuroclass', 'size': '6672'}