Skull Stripping in AFNI

fmrihead · June 6, 2024, 11:10am

AFNI version info (afni -ver): 24.1

Hello Afni Gurus!

Our lab uses proc.py to preprocess ME data. We use an MP2RAGE sequence for registration.
Our Pipeline is comprised of the following general steps -
(1) Dicom to AFNI
(2) Creating an anatomical mask using the INV2 series acquired during MP2RAGE, and using that mask for skull stripping:

3dUnifize with INV2 as an input
3dWarp -deoblique with both INV2 and the unified series from MP2RAGE
3dresample post deobliquing
3dskullstrip post 3dresample --> creating SSanat2
3dAutomask --> 3dcalc to create mask

I'm attaching this bit of code here for clarity -

#### correction for the inhomogeneity of the field to file anat2
3dUnifize -prefix ${subject}.7TMS.${session}.anat2.uni $path_session/'anat2'/${subject}.7TMS.${session}.anat2+orig.

#### Deoblique for files anat 2 (corrected) and anat
3dWarp -deoblique -prefix ${subject}.7TMS.${session}.anat2.Deobl ${subject}.7TMS.${session}.anat2.uni+orig

3dWarp -deoblique -prefix ${subject}.7TMS.${session}.anat.Deobl $path_session/'anat'/${subject}.7TMS.${session}.anat+orig.

#### correct the orientation (that is changed when deoblique)
3dresample -orient ASR -rmode NN -prefix ${subject}.7TMS.${session}.anat2.DeoblASR -input ${subject}.7TMS.${session}.anat2.Deobl+orig.

3dresample -orient ASR -rmode NN -prefix ${subject}.7TMS.${session}.anat.DeoblASR -input ${subject}.7TMS.${session}.anat.Deobl+orig.
 
#### 3dskullstrip for anat2 and create mask

3dSkullStrip -input ${subject}.7TMS.${session}.anat2.DeoblASR+orig -prefix SSanat2

3dAutomask -prefix Automask SSanat2+orig 

#### create anat mask out of the automask
3dcalc -a Automask+orig -b ${subject}.7TMS.${session}.anat.DeoblASR+orig -expr 'a*b' -prefix ${subject}.7TMS.${session}.anatMASK

Following these steps we use proc.py. We use SSwarper to create the inputs proc.py requires in order to run (anatQQ/SS). As input, SSwarper receives the final product of the bit of code pasted above - which is a skull-stripped MP2RAGE image. We are aware of the fact that SSwarper requires a non-skull stripped input, but unfortunately it isn't doing a good job with stripping our original MP2RAGE image :(.
This is the bit of code that follows -

cd $path

 
for file in anatQQ.*; do
    if [[ -f $file ]]; then
    	echo 'anat.QQ FILES exsist'
    	break
        else
        echo 'preforming @SSwraper'
        export OMP_NUM_THREADS 4
        @SSwarper                                                                      \
		   -input  $subject.$experiment.$session.anatMASK+orig.BRIK                       \
		   -base   $abin_path/MNI152_2009_template_SSW.nii.gz                                  \
		   -subid  $subject.$experiment.$session
    fi
done
  

for scan in ${scans[@]}
do

echo "preforming proc.py to scan $scan"


afni_proc.py \
     -subj_id $subject.$experiment.$session.$scan \
     -out_dir ${output_path}/$subject.$experiment.$session.$scan.results \
     -script ${output_path}/proc.$subject.$experiment.$session.$scan \
     -blocks despike tshift align tlrc volreg mask combine scale regress \
     -radial_correlate_blocks tcat volreg                      \
     -blip_forward_dset $subject.$experiment.$session.AP+orig.HEAD -align_unifize_epi local \
     -blip_reverse_dset $subject.$experiment.$session.PA+orig.HEAD                     \
     -copy_anat anatSS.$subject.$experiment.$session.nii \
     -anat_has_skull no \
     -dsets_me_run $subject.$experiment.$session.${scan}*0?+orig.HEAD \
     -echo_times 13.2 34.72 56.24 \
     -reg_echo 1                                               \
     -tcat_remove_first_trs 4                                  \
     -tshift_interp -wsinc9                                    \
     -align_opts_aea -cost lpc+ZZ -check_flip            \
     -tlrc_base MNI152_2009_template_SSW.nii.gz               \
     -tlrc_NL_warp                                             \
     -tlrc_NL_warped_dsets anatQQ.$subject.$experiment.$session.nii anatQQ.$subject.$experiment.$session.aff12.1D anatQQ.$subject.$experiment.${session}_WARP.nii \
     -volreg_align_to MIN_OUTLIER                              \
     -volreg_align_e2a                                         \
     -volreg_tlrc_warp                                         \
     -mask_epi_anat yes                                        \
     -mask_apply anat \
     -combine_method m_tedana                                    \
     -html_review_style pythonic \
     -scr_overwrite
     

tcsh -xef $output_path/proc.$subject.$experiment.$session.$scan 2>&1 | tee $output_path/output.proc.$subject.$experiment.$session.$scan

## copy final output to shared folder
echo "copying final errts and anat files to shared folder: $shared_output"

cp ${output_path}/$subject.$experiment.$session.$scan.results/errts* $shared_output

done

cp ${output_path}/$subject.$experiment.$session.$scan.results/anat_final* $shared_output/anat_final.$subject.$experiment.$session

done
done

The problem is that no matter what options we try, SSwarper (and also SSwarper2) sometimes takes a bite out of the skull stripped input it receives, taking out part of the cortex / cerebellum.
For instance -

Any suggestions on how we can work this out? We have a method for skull stripping without SSwarper. We want to move to MNI, and use that product with proc.py, but to avoid skull stripping any further.

Is there a way to add a flag to SSwarper to prevent it from Skull Stripping, and just warping to MNI?
Any way creating the files proc.py needs without using SSwarper?

Thank you in advance!!!!!

ptaylor · June 6, 2024, 5:53pm

Howdy-

If you already have a skullstripped anatomical, then you could just use 3dQwarp directly, instead of @SSwarper or sswarper. Those latter programs are wrappers for 3dQwarp, and they do work assuming that they need to pair skull removal with alignment. It wouldn't surprise me that results might be a bit odd sometimes if the properties of the inputs are changed notably, as here.

There are a couple different ways to do it---call 3dAllineate before 3dQwarp, or just use the latter with the -allineate option. In the end, to provide the pre-calculated warp to afni_proc.py, you still need to provide 3 datasets: the warped dataset in template space; the affine part of the alignment (*aff12.1D file); and the nonlinear part of the warp (*_WARP* dset).

I would recommend perhaps starting with this, which will assume that your initial skullstripping is good, and we'll use the same *SSW* reference template base that has multiple bricks:

set Basedset = MNI152_2009_template_SSW.nii.gz

 3dQwarp        \
        -allineate \
        -allineate_opts "-twopass -cmass" \
        -base     "${Basedset}[0]"            \
        -source   DSET_IN    \
        -lpa                           \
        -pblur 0 1 \
        -weight "${Basedset}[2]"              \
        -inedge                  \
        -prefix DSET_OUT.nii

Your input anatomical is DSET_IN, and the result aligned to the template will be DSET_OUT.nii, while there will also be a DSET_OUT_WARP.nii.

NB: There will be an *.aff12.1D matrix file also output, but that that file is just intermediate and it shouldn't be used with the warp here when passed to afni_proc.py. It can be a little confusing, but when 3dQwarp does an initial affine estimate, it will already include that information in the output *_WARP* dset; but when @SSwarper or sswarper2 uses affine early on, it does not include it directly in the warp, so both have to be given. Therefore, in this scenario, we have to hand afni_proc.py a placeholder identity matrix file in the set of 3 inputs. You can make one called mat_identity.aff12.1D whose contents are just this (I put spaces to break up the quartets visually):

1 0 0 0   0 1 0 0   0 0 1 0

Then, if alignment is good, you can provide the results to afni_proc.py like:

 -tlrc_NL_warped_dsets  DSET_OUT.nii      \
       mat_identity.aff12.1D   \
        DSET_OUT_WARP.nii   \

How does that sound?

--pt

dglen · June 7, 2024, 12:41am

If the input anatomical dataset, then you should provide this option (same as before):

           -anat_has_skull           no                                  \

If there are no warps provided, then afni_proc.py will call auto_warp.py to do the affine and nonlinear alignment to a template. auto_warp.py calls @auto_tlrc and 3dQwarp. It's a bit different than the 3dQwarp-allineate method because @auto_tlrc uses 3dWarpDrive by default and not 3dAllineate for its affine alignment. For afni_proc.py, you would still need to provide a template, of course, along with the option to do the nonlinear warping and be sure to include the tlrc block in the processing blocks.

           ... \
           -tlrc_base   MNI152_2009_template.nii.gz   \
           -tlrc_NL_warp   \

In the end, there are fewer options for afni_proc.py because you won't provide the warp and the affine transformation.

fmrihead · June 16, 2024, 12:43pm

Hi Paul,

Thanks for responding - I tried your suggestion, and I'm afraid the end product wasn't very successful.

Code -

echo "Running 3dQwarp..."
3dQwarp \
    -allineate \
    -allineate_opts "-twopass -cmass" \
    -base     "${Basedset}[0]" \
    -source   "$source_file" \
    -lpa \
    -pblur 0 1 \
    -weight "${Basedset}[2]" \
    -inedge \
    -prefix "${subject}.7TMS.${session}.3dQwarp.nii"

echo "3dQwarp completed."

As stated in your suggestion - I see 4 output files -
(1) 3dQwarp_Allin.aff12.1D
(2) 3dQwarp_Allin.nii
(3) 3dQwarp_WARP.nii
(4) 3dQwarp.nii

Looking at file (4) - it seems the output is quite smeared..
I suspect this is because our anatomical input is not perfectly skull stripped -
Screenshot 2024-06-16 at 15.41.53

Do you think this might be it?

I'm trying other methods to skull strip and will update if 3dQwarp performs better with those inputs.

Thanks again, will update!

fmrihead · June 17, 2024, 7:31am

Hi guys,

I'm still getting weird results with 3dQwarp (trying out multiple inputs). The warp to MNI isn't working out, and the result Warp is distorted.
I'm attaching my code and an example of the results - would love to hear your input on what you think is going on.

Thanks again.

(1) Code -

Basedset=$abin_path/'MNI152_2009_template_SSW.nii.gz'
data='/Volumes/Labs/ramot/All_data/AFNI_raw'/$experiment/$subject/'epi'/$session/'anatmask'
mkdir -p '/Volumes/Labs/ramot/All_data/AFNI_raw'/$experiment/$subject/'epi'/$session/'3dQwarp'
warp_path='/Volumes/Labs/ramot/All_data/AFNI_raw'/$experiment/$subject/'epi'/$session/'3dQwarp'



echo 'Creating directory:' $warp_path
cd $warp_path

source_file=$data/brainmask.nii
if [ ! -f $source_file ]; then
    echo "Error: Source file not found: $source_file"
    exit 1
fi

echo "Running 3dQwarp..."
3dQwarp \
    -allineate \
    -allineate_opts "-twopass -cmass" \
    -base     "${Basedset}[0]" \
    -source   "$source_file" \
    -lpa \
    -pblur 0 1 \
    -weight "${Basedset}[2]" \
    -inedge \
    -prefix "${subject}.7TMS.${session}.3dQwarp_brainmask_freesurfer.nii"

echo "3dQwarp completed."

(2) Example result -

Thanks again.

ptaylor · June 17, 2024, 2:45pm

What does the output of the Allineate part look like, which I guess is 3dQwarp_Allin.nii in the file naming above?

The first/solo image above didn't look that bad, it is hard to tell from a single slice, but it looked symmetric and undistorted. The pair of images above (sagittal+coronal views) do look distorted).

--pt

fmrihead · June 18, 2024, 10:34am

Hi,
This is the Allineate part, I think it looks OK.

fmrihead · June 24, 2024, 6:11am

Hi again, just wondering if you had any other suggestions to try.
I've been playing with the code but still with no success with 3dQwarp.

Thanks again!

ptaylor · June 25, 2024, 9:06pm

Howdy-

Sorry for the delayed reply (at OHBM conference+brainhack).

I'm not exactly sure. It might be easiest if I try looking at the dataset itself and run this on my computer, if you can share the data. I will send a PM.

-pt