valosekj / dcm-brno

GNU General Public License v3.0
0 stars 0 forks source link

PAM50 template registration - what `-param` to use? #18

Open valosekj opened 2 months ago

valosekj commented 2 months ago

I need to register C spine T2w images (either 0.28 x 0.28 x 1.3 mm or 0.8 iso mm) of subjects with cord compression to the PAM50 template. I have SC seg and disc labels (1 to 10) available.

[!NOTE]
The primary usage of the registration to the template is to use its warping fields as initialization for template-DWI registration

[!NOTE]
I am interested in extracting DTI metrics from the C3 level (i.e., level above the compression), meaning that the registration does not have to be perfect at the lower levels with compression sites.

I compared the following -param configurations.

Link to a presentation comparing different configurations. Link to the data used (without registration results, as there would be too many files due to different configurations).

Option 1: three steps

-param step=1,type=seg,algo=centermassrot:step=2,type=seg,algo=syn,slicewise=1,smooth=0,iter=5:step=3,type=im,algo=syn,slicewise=1,smooth=0,iter=3

https://github.com/valosekj/dcm-brno/blob/866e408c01ddc5aac46048ce0bde062cf6786505/02_processing_scripts/02_process_data.sh#L242-L244

Inspired by https://github.com/spine-generic/spine-generic/blob/master/process_data.sh#L165-L166 from spine-generic. Note that the spine-generic processing was done on healthy subjects (i.e., without cord compression).

For Option 1, I also tried tweaking the template->dwi registration; see this commit and slide 16 in the attached slide deck. UPDATE: related issue: https://github.com/valosekj/dcm-brno/issues/19

Option 2: two steps

-param step=1,type=seg,algo=centermassrot:step=2,type=seg,algo=columnwise

https://github.com/valosekj/dcm-brno/blob/6af1607a6a0772d4addc865d4f1cc9c2a0d48780/02_processing_scripts/02_process_data.sh#L244-L246

Inspired by the SCT course: slide 81. Note that algo=columnwise is used to deal with the compressed cord.

Option 3: two steps

-param step=1,type=seg,algo=centermassrot:step=2,type=im,algo=syn,iter=5,slicewise=1,metric=CC,smooth=0

https://github.com/valosekj/dcm-brno/blob/ba23a6d006f9c4ffe4a342cb3b1abfec6b161682/02_processing_scripts/02_process_data.sh#L243-L245

Inspired by https://github.com/sct-pipeline/spine-park/blob/main/batch_processing.sh#L143-L146 from spine-park. Note that the spine-park processing is done on PD patients; probably without severe cord compression.


@jcohenadad, I would greatly appreciate your opinion.

jcohenadad commented 2 months ago

@valosekj i cannot assess via keynote slides-- i need to see the qc report-- i'm going to regenerate it but it would save me some time if you could just write the syntaxes of the parameters you try so i can just copy/paste the code, run it, and open qc-- thx

valosekj commented 2 months ago

Sure!

Link to the data used.

Go to the data directory after download:

~/Downloads/sub-1860B6472B/

Option 1:

cd anat
# Register T2w to PAM50 template using C3 and C5 mid-vertebral levels
sct_register_to_template -i sub-1860B6472B_ses-1860B_T2w.nii.gz -s sub-1860B6472B_ses-1860B_T2w_seg.nii.gz -l sub-1860B6472B_ses-1860B_T2w_label-disc_c3c5.nii.gz -c t2 -param step=1,type=seg,algo=centermassrot:step=2,type=seg,algo=syn,slicewise=1,smooth=0,iter=5:step=3,type=im,algo=syn,slicewise=1,smooth=0,iter=3 -qc ../qc -qc-subject sub-1860B6472B_ses-1860B
# Rename warping fields for clarity
mv warp_template2anat.nii.gz warp_template2T2w.nii.gz
mv warp_anat2template.nii.gz warp_T2w2template.nii.gz

cd ../dwi
# Register template->dwi (using T2w-to-template as initial transformation)
sct_register_multimodal -i $SCT_DIR/data/PAM50/template/PAM50_t1.nii.gz -iseg $SCT_DIR/data/PAM50/template/PAM50_cord.nii.gz -d sub-1860B6472B_ses-1860B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean.nii.gz -dseg sub-1860B6472B_ses-1860B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean_seg.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3 -initwarp ../anat/warp_template2T2w.nii.gz -initwarpinv ../anat/warp_T2w2template.nii.gz -owarp warp_template2dwi.nii.gz -owarpinv warp_dwi2template.nii.gz -qc ../qc -qc-subject sub-1860B6472B/ses-1860B

Option 2:

# Register T2w to PAM50 template using C3 and C5 mid-vertebral levels
sct_register_to_template -i sub-1860B6472B_ses-1860B_T2w.nii.gz -s sub-1860B6472B_ses-1860B_T2w_seg.nii.gz -l sub-1860B6472B_ses-1860B_T2w_label-disc_c3c5.nii.gz -c t2 -param step=1,type=seg,algo=centermassrot:step=2,type=seg,algo=columnwise -qc ../qc -qc-subject sub-1860B6472B_ses-1860B
# Rename warping fields for clarity
mv warp_template2anat.nii.gz warp_template2T2w.nii.gz
mv warp_anat2template.nii.gz warp_T2w2template.nii.gz

cd ../dwi
# Register template->dwi (using T2w-to-template as initial transformation)
sct_register_multimodal -i $SCT_DIR/data/PAM50/template/PAM50_t1.nii.gz -iseg $SCT_DIR/data/PAM50/template/PAM50_cord.nii.gz -d sub-1860B6472B_ses-1860B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean.nii.gz -dseg sub-1860B6472B_ses-1860B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean_seg.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3 -initwarp ../anat/warp_template2T2w.nii.gz -initwarpinv ../anat/warp_T2w2template.nii.gz -owarp warp_template2dwi.nii.gz -owarpinv warp_dwi2template.nii.gz -qc ../qc -qc-subject sub-1860B6472B/ses-1860B

Option 3:

# Register T2w to PAM50 template using C3 and C5 mid-vertebral levels
sct_register_to_template -i sub-1860B6472B_ses-1860B_T2w.nii.gz -s sub-1860B6472B_ses-1860B_T2w_seg.nii.gz -l sub-1860B6472B_ses-1860B_T2w_label-disc_c3c5.nii.gz -c t2 -param step=1,type=seg,algo=centermassrot:step=2,type=im,algo=syn,iter=5,slicewise=1,metric=CC,smooth=0 -qc ../qc -qc-subject sub-1860B6472B_ses-1860B
# Rename warping fields for clarity
mv warp_template2anat.nii.gz warp_template2T2w.nii.gz
mv warp_anat2template.nii.gz warp_T2w2template.nii.gz

cd ../dwi
# Register template->dwi (using T2w-to-template as initial transformation)
sct_register_multimodal -i $SCT_DIR/data/PAM50/template/PAM50_t1.nii.gz -iseg $SCT_DIR/data/PAM50/template/PAM50_cord.nii.gz -d sub-1860B6472B_ses-1860B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean.nii.gz -dseg sub-1860B6472B_ses-1860B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean_seg.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=seg,algo=bsplinesyn,slicewise=1,iter=3 -initwarp ../anat/warp_template2T2w.nii.gz -initwarpinv ../anat/warp_T2w2template.nii.gz -owarp warp_template2dwi.nii.gz -owarpinv warp_dwi2template.nii.gz -qc ../qc -qc-subject sub-1860B6472B/ses-1860B
jcohenadad commented 2 months ago

First observation, the DWI segmentation is slightly over-segmented, therefore a seg-based registration will reflect the over-segmentation. Also: there is a slight AL shift in the top slices:

ezgif-5-d045598072

Then, for the pipeline, i would go with a 2-step approach as you did, but with very little DOFs for the first step (T2w reg), and leave the high DOFs for the DWI reg. Something like this:

cd anat
# Register T2w to PAM50 template using C3 and C5 mid-vertebral levels
sct_register_to_template -i sub-1860B6472B_ses-1860B_T2w.nii.gz -s sub-1860B6472B_ses-1860B_T2w_seg.nii.gz -l sub-1860B6472B_ses-1860B_T2w_label-disc_c3c5.nii.gz -c t2 -param step=1,type=seg,algo=centermassrot:step=2,type=seg,algo=syn,slicewise=1,smooth=0,iter=5 -qc ../qc -qc-subject sub-1860B6472B_ses-1860B
# Rename warping fields for clarity
mv warp_template2anat.nii.gz warp_template2T2w.nii.gz
mv warp_anat2template.nii.gz warp_T2w2template.nii.gz
cd ../dwi
# Register template->dwi (using T2w-to-template as initial transformation)
sct_register_multimodal -i /Users/julien/code/spinalcordtoolbox/data/PAM50/template/PAM50_t1.nii.gz -iseg /Users/julien/code/spinalcordtoolbox/data/PAM50/template/PAM50_cord.nii.gz -d sub-1860B6472B_ses-1860B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean.nii.gz -dseg sub-1860B6472B_ses-1860B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean_seg.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=im,algo=syn,slicewise=1,iter=5 -initwarp ../anat/warp_template2T2w.nii.gz -initwarpinv ../anat/warp_T2w2template.nii.gz -owarp warp_template2dwi.nii.gz -owarpinv warp_dwi2template.nii.gz -qc ../qc -qc-subject sub-1860B6472B/ses-1860B

QC report: qc.zip

valosekj commented 2 months ago

First observation, the DWI segmentation is slightly over-segmented, therefore a seg-based registration will reflect the over-segmentation. Also: there is a slight AL shift in the top slices:

Yes! I also noticed that later yesterday! After some investigation (https://github.com/valosekj/dcm-brno/issues/19), I switched from contrast-agnostic to sct_deepseg_sc.

Then, for the pipeline, i would go with a 2-step approach as you did, but with very little DOFs for the first step (T2w reg), and leave the high DOFs for the DWI reg. Something like this:

cd anat
# Register T2w to PAM50 template using C3 and C5 mid-vertebral levels
sct_register_to_template -i sub-1860B6472B_ses-1860B_T2w.nii.gz -s sub-1860B6472B_ses-1860B_T2w_seg.nii.gz -l sub-1860B6472B_ses-1860B_T2w_label-disc_c3c5.nii.gz -c t2 -param step=1,type=seg,algo=centermassrot:step=2,type=seg,algo=syn,slicewise=1,smooth=0,iter=5 -qc ../qc -qc-subject sub-1860B6472B_ses-1860B
# Rename warping fields for clarity
mv warp_template2anat.nii.gz warp_template2T2w.nii.gz
mv warp_anat2template.nii.gz warp_T2w2template.nii.gz
cd ../dwi
# Register template->dwi (using T2w-to-template as initial transformation)
sct_register_multimodal -i /Users/julien/code/spinalcordtoolbox/data/PAM50/template/PAM50_t1.nii.gz -iseg /Users/julien/code/spinalcordtoolbox/data/PAM50/template/PAM50_cord.nii.gz -d sub-1860B6472B_ses-1860B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean.nii.gz -dseg sub-1860B6472B_ses-1860B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean_seg.nii.gz -param step=1,type=seg,algo=centermass:step=2,type=im,algo=syn,slicewise=1,iter=5 -initwarp ../anat/warp_template2T2w.nii.gz -initwarpinv ../anat/warp_T2w2template.nii.gz -owarp warp_template2dwi.nii.gz -owarpinv warp_dwi2template.nii.gz -qc ../qc -qc-subject sub-1860B6472B/ses-1860B

Thank you! Pipeline updated in https://github.com/valosekj/dcm-brno/commit/fd1ee8b66a022901ca511f5d32c37259864aebb2. Trying it now across multiple subjects.

valosekj commented 2 months ago

WM atlas warped to DWI space seems to be "undersegmented" for some subjects despite the fact that the DWI mean seg seems reasonable. Two examples:

sub-3013B6405B_ses-6405B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean.nii.gz: - focus mainly on 4th and 5th top slices ![Kapture 2024-07-24 at 06 40 31](https://github.com/user-attachments/assets/ac74a434-3348-4db2-8884-c32a8a1433db)
sub-3013B6405B_ses-3013B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean.nii.gz: - focus mainly on 3rd, 4th, and 5th top slices ![Kapture 2024-07-24 at 06 39 18](https://github.com/user-attachments/assets/247c2bd2-26cd-4a1e-8aa1-34eeeae22cb7)

Other subjects: sub-3060B6384B_ses-6384B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean.nii.gz sub-3120B6558B_ses-3120B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean.nii.gz


Suggestions from the 2024-07-23 DCM norm meeting:

UPDATE:

iter=3 vs iter=5 -- iter=3 seems to be less undersegmented:

sub-3013B6405B_ses-3013B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean.nii.gz:

Kapture 2024-07-24 at 08 37 42

sub-3060B6384B_ses-6384B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean.nii.gz:

Kapture 2024-07-24 at 08 36 04

jcohenadad commented 2 months ago

that looks good to me. Note that the warped atlas is thresholded on the QC report so in reality it is 'bigger'

valosekj commented 2 months ago

Trying also bsplinesyn vs syn, so we now have three variants:

sub-3013B6405B_ses-3013B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean: gif

sub-3060B6384B_ses-6384B_acq-ZOOMit_dir-AP_dwi_crop_crop_moco_dwi_mean: gif

(for some reason the GIFs now animate only once...but when downloaded, they work fine)

algo=bsplinesyn,iter=5 looks even better than algo=syn,iter=3 to me