nipreps / fmriprep

fMRIPrep is a robust and easy-to-use pipeline for preprocessing of diverse fMRI data. The transparent workflow dispenses of manual intervention, thereby ensuring the reproducibility of the results.
https://fmriprep.org
Apache License 2.0
638 stars 294 forks source link

sub-<subjID>_run-<runID>_desc-aseg_dseg.nii.gz looks different from aseg.mgz #2642

Closed zhengchencai closed 2 years ago

zhengchencai commented 2 years ago

What happened?

Hi there,

Could you please help me to figure out the following problem?

Short story: sub-<subjID>_run-<runID>_desc-aseg_dseg.nii.gz under fmriprep/sub-<>/anat folder looks different to aseg.mgz under fmriprep/source/freesurfer/sub-<>/mri folder (see snapshot). I understand they are not in the same space, but the first one seems to get some extra "tissues" outside the brain mask. I really hope I am wrong.

Long story: I was trying to resample an individual segmentation (e.g. Schaefer2018_200Parcels_17Networks.mgz) in freesurfer space to BOLD resolution. Searching from this repo and NeuroStars, it seems that fmriprep used antsApplyTransforms. To make sure I don't make a mistake, I thought if I could resample aseg.mgz under the source/freesurfer folder to BOLD ref and end up with an identical result to _space-T1w_desc-aseg_dseg.nii.gz under the func folder, then I can proceed with confidence. However, I could not reproduce what fmriprep did and posted questions https://neurostars.org/t/resampling-volume-atlas-to-bold-space/20696 and https://github.com/nipreps/fmriprep/issues/2632. In fmriprep source code https://github.com/nipreps/fmriprep/blob/master/fmriprep/workflows/bold/registration.py line242 t1w_aseg is already projected on T1w from freesurfer space before applying antsApplyTransforms with --transform identity, tracking this down I couldn't find exactly how fmriprep did the registration to get t1w_aseg, because I am not familiar with Nipype. So I suppose (hopefully I am wrong) sub-<subjID>_run-<runID>_desc-aseg_dseg.nii.gz is the t1w_aseg and compared it with aseg.mgz, but they don't look similar from the snapshot. I also tried mri_label2vol --seg aseg.mgz --temp rawavg.mgz --o aseg-in-rawavg.mgz --regheader aseg.mgz suggested by freesurfer to convert aseg back to T1w and applied antsApplyTransforms with --input $SUBJECTS_DIR/$subject_name/mri/aseg-in-rawavg.mgz --transform identity --reference-image $fmriprep_DIR/$subject_name/func/${subject_name}_${task_name}_run-1_space-T1w_boldref.nii.gz, the result still looks different from _space-T1w_desc-aseg_dseg.nii.gz

Thanks lot

What command did you use?

antsApplyTransforms --default-value 0 --float 1 --interpolation MultiLabel \
    --input $SUBJECTS_DIR/$subject_name/mri/aseg.nii.gz \
    --output $fmriprep_DIR/$subject_name/func/${subject_name}_${task_name}_run-1_space-T1w_desc-aseg_manual.nii.gz \
    --reference-image $fmriprep_DIR/$subject_name/func/${subject_name}_${task_name}_run-1_space-T1w_boldref.nii.gz \
    --transform $fmriprep_DIR/$subject_name/anat/${subject_name}_run-1_from-fsnative_to-T1w_mode-image_xfm.txt
 antsApplyTransforms --dimensionality 3 --float 1 --interpolation genericLabel \
    --input $SUBJECTS_DIR/$subject_name/mri/aseg-in-rawavg.mgz \
    --transform identity \
    --reference-image $fmriprep_DIR/$subject_name/func/${subject_name}_${task_name}_run-1_space-T1w_boldref.nii.gz \
    --output $fmriprep_DIR/$subject_name/func/${subject_name}_${task_name}_run-1_space-T1w_desc-aseg--in-rawavg.nii.gz```

What version of fMRIPrep are you running?

20.2.4

How are you running fMRIPrep?

Singularity

Is your data BIDS valid?

Yes

Are you reusing any previously computed results?

No

Please copy and paste any relevant log output.

No response

Additional information / screenshots

Picture3 image

zhengchencai commented 2 years ago

Researching this issue a bit I could reproduce <...>_space-T1w_desc-aseg_dseg.nii.gz generated by fmriprep.

First of all, this bug https://github.com/nipreps/smriprep/pull/268 is the reason why sub-<subjID>_run-<runID>_desc-aseg_dseg.nii.gz under fmriprep/sub-<>/anat folder looks different to aseg.mgz under fmriprep/source/freesurfer/sub-<>/mri folder

Then I rerun fmriprep with 21.0.0rc2 version and did the following to reproduce <...>_space-T1w_desc-aseg_dseg.nii.gz from aseg.mgz in freefurfer folder, after reading the source code. Tested with 3 subjects and got identical segmentations resampled in BOLD resolution

# fmriprep used two steps 
# step 1, resample aseg.mgz to T1w.nii.gz and convert it to aseg.nii.gz, smriprep code below
# step 1a, use fsnative2t1w_xfm to calculate lta file () from T1.mgz to T1w.nii.gz
#fsnative2t1w_xfm = pe.Node(
#        RobustRegister(auto_sens=True, est_int_scale=True), name="fsnative2t1w_xfm"
#    )
# RobustRegister is from freesurfer mri_robust_register auto_sens is --satit est_int_scale is --iscale
mri_robust_register --mov $SUBJECTS_DIR/$subject_name/mri/T1.mgz \
    --dst $root_DIR/$subject_name/anat/${subject_name}_run-01_T1w.nii.gz \
    --lta $work_DIR/${subject_name}_run-01_fsnative2t1w_xfm.lta \
    --satit \
    --iscale 
# step 1b, use fs.ApplyVolTransform with nearest interpolation to transfer 
# aseg.mgz to T1w space and convert to aseg.nii.gz
#    # Resample from T1.mgz to T1w.nii.gz, applying any offset in fsnative2t1w_xfm,
#    # and convert to NIfTI while we're at it
#    resample = pe.Node(
#        fs.ApplyVolTransform(transformed_file="seg.nii.gz", interp="nearest"),
#        name="resample",
#    )
# fs.ApplyVolTransform is mri_vol2vol
mri_vol2vol --mov $SUBJECTS_DIR/$subject_name/mri/aseg.mgz \
     --targ $root_DIR/$subject_name/anat/${subject_name}_run-01_T1w.nii.gz \
     --lta $work_DIR/${subject_name}_run-01_fsnative2t1w_xfm.lta \
     --o $work_DIR/${subject_name}_run-01_t1_aseg.nii.gz \
     --nearest

# 2 use antsApplyTransforms resample aseg.nii.gz to bold ref 
# # Resample aseg and aparc in T1w space (no transforms needed), fmriperp code below
#        aseg_t1w_tfm = pe.Node(
#            ApplyTransforms(interpolation='MultiLabel', transforms='identity'),
#            name='aseg_t1w_tfm', mem_gb=0.1)
#        aparc_t1w_tfm = pe.Node(
#            ApplyTransforms(interpolation='MultiLabel', transforms='identity'),
#            name='aparc_t1w_tfm', mem_gb=0.1)

antsApplyTransforms --dimensionality 3 --float 1 --interpolation MultiLabel \
    --input $work_DIR/${subject_name}_run-01_t1_aseg.nii.gz \
    --reference-image $fmriprep_DIR/$subject_name/func/${subject_name}_${task_name}_run-2_space-T1w_boldref.nii.gz \
    --output $fmriprep_DIR/$subject_name/func/${subject_name}_${task_name}_run-2_space-T1w_desc-aseg-rep.nii.gz \
    --transform identity

According to https://surfer.nmr.mgh.harvard.edu/fswiki/FsAnat-to-NativeAnat freesurfer suggested using mri_label2vol to transform segmentation back to native space, so I tried the following. It could reproduce the above process, but in one subject only very few voxels are labeled differently, not sure why.

# using freesurfer mri_label2vol https://surfer.nmr.mgh.harvard.edu/fswiki/FsAnat-to-NativeAnat
mri_label2vol --seg $SUBJECTS_DIR/$subject_name/mri/aseg.mgz \
    --temp $SUBJECTS_DIR/$subject_name/mri/rawavg.mgz \
    --o $work_DIR/${subject_name}_run-01_t1_aseg_rep.mgz \
    --regheader $SUBJECTS_DIR/$subject_name/mri/aseg.mgz  

antsApplyTransforms --dimensionality 3 --float 1 --interpolation MultiLabel \
    --input $work_DIR/${subject_name}_run-01_t1_aseg_rep.mgz\
    --reference-image $fmriprep_DIR/$subject_name/func/${subject_name}_${task_name}_run-2_space-T1w_boldref.nii.gz \
    --output $fmriprep_DIR/$subject_name/func/${subject_name}_${task_name}_run-2_space-T1w_desc-aseg-rep.nii.gz \
    --transform identity

Finally, I tried the very first logic, using _run-1_from-fsnative_to-T1w_mode-image_xfm.txt provided by fmriprep to ignore recalculation of transformation from fs to T1w. The following code still can not reproduce <...>_space-T1w_desc-aseg_dseg.nii.gz. According to this bug https://github.com/nipreps/smriprep/pull/268 I was expecting the _run-1_from-fsnative_to-T1w_mode-image_xfm.txtbe different from fmriprep versions, but the one I got from 20.2.4 is identical to 21.0.0rc2, this explains why the following code gave me the same result as the snapshot showed in the first post. Maybe I am wrong so _run-1_from-fsnative_to-T1w_mode-image_xfm.txt is not the transformation from fs to T1w. Also, the above mri_robust_register and then mri_vol2vol process could not reproduce <...>_space-T1w_desc-aseg_dseg.nii.gz generated by fmriprep 20.2.4

antsApplyTransforms --dimensionality 3 --float 1 --interpolation MultiLabel \
    --input $SUBJECTS_DIR/$subject_name/mri/aseg.mgz\
    --reference-image $fmriprep_DIR/$subject_name/func/${subject_name}_${task_name}_run-2_space-T1w_boldref.nii.gz \
    --output $fmriprep_DIR/$subject_name/func/${subject_name}_${task_name}_run-2_space-T1w_desc-aseg-rep.nii.gz \
    --transform $fmriprep_DIR/$subject_name/anat/${subject_name}_run-1_from-fsnative_to-T1w_mode-image_xfm.txt