AFNI's volreg VS ANTs motion correction

[grandjeanlab]

Describe the bug In a dataset with func images with 5 slices, MC correction fails. Func with more slices pass. e.g. /scratch/m/mchakrav/gjoanes/multirat/ds/02008/bids/sub-0200804/ses-1/func from dataset /scratch/m/mchakrav/gjoanes/multirat/ds/02008/bids Other scans with more slices pass e.g. /scratch/m/mchakrav/gjoanes/multirat/ds/02008/bids/sub-0200805/ses-1/func

Describe RABIES call singularity run -B ${PROJECT}/template:/template -B ${bids_folder}:/rabies_in:ro \ -B ${preprocess_folder}:/rabies_out \ /project/m/mchakrav/gjoanes/rabies.0.2.0-dev.simg -p MultiProc preprocess /rabies_in /rabies_out \ --template_reg_script multiRAT \ --coreg_script Affine \ --no_STC \ --anat_template ${template} \ --brain_mask ${mask} \ --WM_mask ${wm} \ --CSF_mask ${csf} \ --vascular_mask ${csf} \ --labels ${atlas} \ --TR ${TR}s \ --commonspace_resampling 0.2x0.2x0.2

Attach log file /scratch/m/mchakrav/gjoanes/multirat/preprocess/02008/rabies_preprocess.log

Attach QC_report

Additional context we ran into this earlier and it was due to how AntsMotionCor resamples the images in the first iteration. Maybe this is still too stringent? https://github.com/CoBrALab/RABIES/blob/ebc502e69314b5b09375e3f0c2f70929c9b4baab/rabies/preprocess_pkg/utils.py#L385

[gdevenyi]

Is 5 slices even valid data?

[grandjeanlab]

not my dataset. When doing sensory-evoked fMRI, it is not uncommon to only record slices on the sensory cortex.

[gdevenyi]

Not sure if this is relevant, but both the fixed and moving images are empty:

[gdevenyi]

I'm trying to figure this out, but you don't seem to have saved your run commands in a bash script anywhere, so I don't know what any of your directories are, can you do that please.

Looks like this input files are also empty?

14:57:13 ⌂69% [gdevenyi@monster:~/niagara-mchakrav/gjoanes/multirat/ds/02008/bids/sub-0200809/ses-1] ⌀ register ./anat/sub-0200809_ses-1_T2w.nii.gz ./func/sub-0200809_ses-1_run-1_bold.nii.gz

[gdevenyi]

@Gab-D-G I think we need a pre-flight QC stage that generates an image of all inputs we can require for issue support.

[Gab-D-G]

@grandjeanlab Indeed RABIES doesn't work currently with a very small number of slices. I wasn't aware that images with such as low number of slices was common, so I'll fix this in the next updates. @gdevenyi I'm not sure what you mean by pre-flight QC. You mean make sure we generate a copy of all input images in case we need them for issue support?

[gdevenyi]

@grandjeanlab Indeed RABIES doesn't work currently with a very small number of slices. I wasn't aware that images with such as low number of slices was common, so I'll fix this in the next updates. @gdevenyi I'm not sure what you mean by pre-flight QC. You mean make sure we generate a copy of all input images in case we need them for issue support?

I'm referring to make a static jpg/PNG of the input files. In this case the failure is due to empty files, not the number of slices

[Gab-D-G]

I see, that's something else that could be added then.

[grandjeanlab]

I've re-ran the script and now it is good. Weird. I'm still lukewarm about ANTs MC. It is very slow. With some 3000 volume datasets, it can several hours to run and may lead to walltime limits. Volreg is much faster, qualitatively equal, and doesn't seem to have number of slice limits.

[Gab-D-G]

Volreg is definitely faster, but would do you have recommendations for comparing the quality of realignment? I've tried in the past to compare both ANTs MC and Volreg and they seemed to differ a fair bit in the 6 motion timecourses, but I can't figure which one is better.

[grandjeanlab]

for comparison, I would consider speed and outcome tested with a visual inspection.

If this helps, I can pick ~10 multiRAT scans with varying level of motion, compare ANTs and Volreg based on speed, motion parameters, and visual inspection. If there is no notable differences upon visual inspection, I would recommend Volreg for speed. What do you think?

[grandjeanlab]

or my awake mouse dataset. You've used it to make a point between heavy motion vs low motion awake runs

[Gab-D-G]

I guess you're right that a visual inspection probably remains our best resort for the QC of head realignment. You can conduct the inspection on your side if you have the time, and I'll do the same when I get the chance. If there's no notable difference I will implement a faster option using volreg.

[grandjeanlab]

So, I've ran a comparison (script below). For 180 volume mouse fMRI, 3dvolreg took 31 +/- 8 sec while ANTs took 1556 +/- 386 sec, or 50x more time.

Visual inspection shows that in low movement subjects, both perform similarly (see attached image), and in heavy motion cases, neither do miracles (e.g. image distortions and muscle activity due to movement cannot be corrected).

Overall, that is a 25min preprocessing time for a very simple case. Some users advocate for >>1000 volumes for proper FC estimates (I got datasets with 3000 volumes). This could lead to max wall-time issues.

Here is the code I ran. It should use same ANTs MoCo parameters as in RABIES.

qsub -I -l 'walltime=24:00:00,mem=10gb, procs=1'

change basefolder to where you want to do the analysis

basefolder=/project/4180000.19/rabies_moco/

afni_runtime=$basefolder'afni.txt' ANTs_runtime=$basefolder'ANTs.txt'

module load ANTs module load afni

cd ${basefolder}

ls bids/ | while read subject do cd bids/$subject

ls . | while read session do cd $session/func

ls *.nii | while read func do

Aout=$(remove_ext $func) 3dresample -inset ${func}[0] -prefix tmp_ref.nii.gz

start=date +%s antsMotionCorr -d 3 -o [${Aout}, ${Aout}_ANTs.nii.gz] -m MI[ tmp_ref.nii.gz , ${func} , 1 , 20 , Regular, 0.2 ] -t Rigid[ 0.1 ] -i 100x50x30 -u 1 -e 1 -l 1 -s 2x1x0 -f 4x2x1 -n 10 end=date +%s

echo $((end-start)) >> $ANTs_runtime

rm tmp_ref.nii.gz

start=date +%s 3dvolreg -prefix ${Aout}_anfi.nii.gz -1Dfile ${Aout}_anfi.txt -float ${func} end=date +%s

echo $((end-start)) >> $afni_runtime

cp ${Aout}_anfi.txt a.1D

1dcat ${Aout}MOCOparams.csv > tmp.txt 1dcat -nonconst tmp.txt > b.1D

1ddot -terse a.1D b.1D >> ${Aout}_cormat.txt rm b.1D rm a.1D rm tmp.txt rm tmp_ref.nii.gz

done cd ../.. done cd ${basefolder} done

[gdevenyi]

This should be in a different issue

[Gab-D-G]

Great! Thanks for carrying this out @grandjeanlab . It's fair to at least provide the option for AFNI's volreg. I'm not sure how complicated it will be to integrate, it will depend on its outputs format, but I'll work on it.

[gdevenyi]

@grandjeanlab I think the ANTs devs would be very interesting in this report of a 50x performance difference. I encourage you to report these findings there, https://github.com/ANTsX/ANTs

[grandjeanlab]

@gdevenyi , I suppose the 50x performance diff depends on the ANTs MoCo parameters. Specifically the number of iterations by the -i argument and the degree of data reduction with the -s argument. To have a proper comparison, we'd need to know the ground truth (maybe synthetic data?), and compare across different parameters. Not an undertaking I am willing to take :-)

[grandjeanlab]

maybe a simpler solution is to offer a fast ANTs MoCo option with fewer iterations, and no data reduction on the -s argument?

[gdevenyi]

Its very hard to tell exactly how any of the AFNI tools work, I don't know the metric, or registration pyramid framework (if any) but it appears:

1. There is only the base-resolution used (-s 0 -f 1).
2. Max iterations is 23 (-i 23)
3. Cost function looks to be raw mean squared error (which they call least squared). This metric isn't available in antsMotionCorr, the closest comparison here would be the GC metric.

[gdevenyi]

antsMotionCorr -d 3 -o [ants_mc_tmp/motcorr,ants_mc_tmp/motcorr.nii.gz,ants_mc_tmp/motcorr_avg.nii.gz] \ -m GC[ fixed.nii.gz , moving.nii.gz , 1 , NA , None, 1 ] -t Rigid[ 0.1 ] -i 23 -u 1 -e 1 -l 1 -s 0 -f 1

[grandjeanlab]

I'll try running the comparison with those parameters see how comparable it is to 3dVolreg

[Gab-D-G]

I've run the comparison on a subject with high motion and I get fairly different outputs somehow between the two. I'm not sure if AFNI's maxdisp.1D is different from ANTs' measure of maximal positioning/framewise displacement. ANTs' call:

antsMotionCorr -d 3 -o [ants_mc_tmp/motcorr,ants_mc_tmp/motcorr.nii.gz,ants_mc_tmp/motcorr_avg.nii.gz] -m MI[ mean.nii.gz , ${func} , 1 , 20 , Regular, 0.2 ] -t Rigid[ 0.1 ] -i 100x50x30 -u 1 -e 1 -l 1 -s 2x1x0 -f 4x2x1 -n 10 -v

AFNI's call:

3dvolreg -prefix prefiltered_func_data_mcf.nii.gz -1Dfile motion.1D -maxdisp1D maxdisp.1D -float ${func}

I'm not sure which one is better, it isn't obvious when looking at the output timeseries visually. It shoudn't be too complicated to integrate AFNI's volreg as an option so it will be offered soon. I tried

antsMotionCorr -d 3 -o [ants_mc_tmp/motcorr,ants_mc_tmp/motcorr.nii.gz,ants_mc_tmp/motcorr_avg.nii.gz] \ -m GC[ fixed.nii.gz , moving.nii.gz , 1 , NA , None, 1 ] -t Rigid[ 0.1 ] -i 23 -u 1 -e 1 -l 1 -s 0 -f 1

but it didn't perform well, and was actually much slower than the normal ANTs script.

[gdevenyi]

, and was actually much slower than the normal ANTs script.

That's probably because I dropped the 20% random sampling for dense sampling for the metric.

As far as I'm concerned, I can't be sure how to make ANTs and AFNI comparable because AFNI doesn't document how the code actually works.

[gdevenyi]

Sidenote, this is what the ANTs guys coded in the ANTsR code for "motion correction", with 0-3 being progressively "more accurate"

https://github.com/ANTsX/ANTsR/blob/60eefd96fedd16bceb4703ecd2cd5730e6843807/R/ants_motion_estimation.R#L126-L164

Edit: they also default to an Affine registration for motion it seems.

[Gab-D-G]

Also if looking at the rigid body parameters it's not 100% obvious how they compare. AFNI's parameters: ANTs' parameters:

It isn't clear which parameter correspond to which movement with ANTs. The range of values isn't even comparable, so I am unsure how I could actually convert AFNI's output to fit with the rest of the pipeline. If we can figure what these outputs are, then maybe it will be possible. It seems like the first three roll parameters of AFNI's correspond to the first three of ANTs, with inverted signs.

[Gab-D-G]

Sidenote, this is what the ANTs guys coded in the ANTsR code for "motion correction", with 0-3 being progressively "more accurate"

https://github.com/ANTsX/ANTsR/blob/60eefd96fedd16bceb4703ecd2cd5730e6843807/R/ants_motion_estimation.R#L126-L164

Great @gdevenyi , that may be the simplest way to go here, just to offer a 'quick and dirty' option from ANTs parameters.

[gdevenyi]

Great @gdevenyi , that may be the simplest way to go here, just to offer a 'quick and dirty' option from ANTs parameters.

Yeah, its probably just worth copying their versions and offering as options, 0-3 at least.

[gdevenyi]

It isn't clear which parameter correspond to which movement with ANTs.

This page describes an AFNI command which will plot the results, https://andysbrainblog.blogspot.com/2012/10/fmri-motion-correction-afnis-3dvolreg.html

[gdevenyi]

One more note, it appears the "default" afni pipeline which does everything adds the option "-cubic" to 3dvolreg

[Gab-D-G]

I implemented the proposed options and tried them on the same subject. In all cases I was using a Rigid transform. Here's the comparison for framewise displacement: Here's for positioning:

The runtimes for each of them was approximately 6sec for 0, 3min for 1, 4min for 2, about 16min for 3 and 1:30min for intrabold. It seems that only option 0 is significantly different from the others, and is likely insufficient for a proper correction.

The intrabold option looks like a good tradeoff between time and performance, so that will become the default, and all options will be also available. Runtime for AFNI was 20sec, so it is still significanlty lower than intrabold, but hopefully the >2-fold reduction in runtime will be helpful @grandjeanlab . Otherwise, there should be some way to tweak the parameters to obtain something faster in between option 0 and intrabold, probably by reducing -m MI[ %s , %s , 1 , %s , regular, 0.25 ] further down from 0.25 to 0.5.

[grandjeanlab]

hey, 2-fold reduction sounds great. Thanks for looking into this. :-)

[Gab-D-G]

@grandjeanlab The latest updates were added to the version 0.2.1 and the container was updated, so you should be able to have the new options available when you re-download it.

[grandjeanlab]

seems to work 100%