Fixed effects analysis with Cifti surface data

rahul-brito · May 9, 2023, 7:31pm

Summary of what happened:

I have effects sizes from my first level analysis of preprocessed input cifti data. In other words, I now have a cifti output file with the effects size for a given contrast for each run in my dataset (which corresponds to 91K array of data) and I want to do a fixed-effects analysis to combine runs per subject/task/session before group analysis.

nilearn compute_fixed_effects (nilearn.glm.compute_fixed_effects - Nilearn) does not accept CIFTI data

I found the below example through this neurostars post (GLM on surface data in NiLearn - #6 by bthirion) but it uses Gifti so I was going to do modifications, but if someone already had code that would help a ton!

github.com

individual-brain-charting/public_analysis_code/blob/master/ibc_public/utils_pipeline.py#L556

      
        
                    write(ffx_variance, os.path.join(
                        write_dir, 'effect_variance_maps/%s_%s.gii' %
                        (contrast, side)))
                    write(ffx_stat, os.path.join(
                        write_dir, 'stat_maps/%s_%s.gii' % (contrast, side)))
            else:
                effect_size_maps, effect_variance_maps, data_available =\
                    _load_summary_stats(
                        subject_dic['output_dir'], session_paradigm, contrast,
                        data_available=True)
                shape = load(effect_size_maps[0]).shape
                if len(shape) > 3:
                    if shape[3] > 1:  # F contrast, skipping
                        continue
                ffx_effects, ffx_variance, ffx_stat = fixed_effects_img(
                    effect_size_maps, effect_variance_maps, mask_img)
                
                save(ffx_effects, os.path.join(
                    write_dir, 'effect_size_maps/%s.nii.gz' % contrast))
                save(ffx_variance, os.path.join(
                    write_dir, 'effect_variance_maps/%s.nii.gz' % contrast))

effigies · May 11, 2023, 12:48pm

Hey, @rahul-brito. Here’s our solution in FitLins for handling NIfTI and CIFTI-2:

github.com

poldracklab/fitlins/blob/ec655c7ec71d92cc81141467bdb44e2bcc253425/fitlins/interfaces/nistats.py#L409-L432

      
        
            if is_cifti:
                ffx_cont, ffx_var, ffx_t = _compute_fixed_effects_params(
                    np.squeeze(
                        [nb.load(fname).get_fdata(dtype='f4') for fname in contrast_imgs]
                    ),
                    np.squeeze(
                        [nb.load(fname).get_fdata(dtype='f4') for fname in variance_imgs]
                    ),
                    precision_weighted=False,
                )
                img = nb.load(filtered_effects[0])
                maps = {
                    'effect_size': dscalar_from_cifti(img, ffx_cont, "effect_size"),
                    'effect_variance': dscalar_from_cifti(img, ffx_var, "effect_variance"),
                    'stat': dscalar_from_cifti(img, ffx_t, "stat"),
                }
            
            
else:
                ffx_res = compute_fixed_effects(contrast_imgs, variance_imgs)
                maps = {

This file has been truncated. show original

rahul-brito · May 11, 2023, 4:51pm

Thank you!!! I will check that out:)

rahul-brito · May 16, 2023, 4:14pm

Thanks again @effigies. Is there also a fitlins implementation for random effects models with cifti? I couldn’t find that. Thank you!

effigies · May 16, 2023, 4:45pm

The random effects model is just the GLM.

github.com

poldracklab/fitlins/blob/ec655c7ec71d92cc81141467bdb44e2bcc253425/fitlins/interfaces/nistats.py#L374-L388

      
        
            if model_type != 'Meta':
                if len(filtered_effects) < 2:
                    raise RuntimeError(
                        "At least two inputs are required for a 't' for 'F' " "second level contrast"
                    )
                if is_cifti:
                    effect_data = np.squeeze(
                        [nb.load(effect).get_fdata(dtype='f4') for effect in filtered_effects]
                    )
                    labels, estimates = level1.run_glm(
                        effect_data, spec['X'].values, noise_model='ols'
                    )
                else:
                    model = level2.SecondLevelModel(smoothing_fwhm=smoothing_fwhm)
                    model.fit(filtered_effects, design_matrix=spec['X'])

rahul-brito · May 16, 2023, 10:34pm

Sorry, really stupid follow-up question. How do I use this? As far as I can tell, the main way to install fitlins is through a singularity/docker container. If I want to integrate this into my current pipeline, I am naively unsure how (without cloning the whole repo?)

effigies · May 16, 2023, 11:17pm

I would probably just copy these lines. It’s not really split into useful functions for importing.

If there is something from fitlins you do want to use, you can just pip install fitlins into the environment you’re working in and then import fitlins in your script. It doesn’t really require Docker. That’s just a way to bundle FSL and other non-Python dependencies. If you’re only using nilearn bits, a conda environment is as good as a Docker image and less cumbersome.

rahul-brito · May 17, 2023, 3:09pm

Awesome that’s what I will do. Usually that feels “wrong” to do so I try to go the importing route, but makes sense its not always the universal solution.

Thank you!