Mouse-Imaging-Centre / pydpiper

Python code for flexible pipeline control
Other
25 stars 10 forks source link

Addressing nu_correct issues #307

Open mcvaneede opened 7 years ago

mcvaneede commented 7 years ago

We recently had a file fail during the non uniformity correction stage. This happened because of some combination between a mask that included the background and some particular combination of intensities in that entire area.

After talking to John, it turns out that it's very important to exclude all background when running non uniformity correction. There is an option to use both a mask produced using a bi modal threshold and the user supplied mask. We should implement this in our pipelines, but then explicitly. We should produce a mask based solely on the bimodlt value from mincstats and the mask gotten from the initial model (or from a better masking procedure) and intersect them to produce a more refined mask. This mask should only be used for the non uniformity correction. (It's possible that some dark areas inside the brain will be excluded in this mask, which is good for nu_correct, but which we don't want for the rest of the pipeline.)

bcdarwin commented 7 years ago

I wonder if better masking would also avoid the issue in which inormalize flips image intensities -- do we know anything about why this occurs?

gdevenyi commented 7 years ago

If you're continuing the transition to minc-toolkit-v2, you could try N4BiasFieldCorrection which offers both weighting and masking options. By default it constructs an Otsu mask to consider only those voxels during correction.

bcdarwin commented 7 years ago

Interesting.

Also on the topic of masking, have you had good results using BEAST or something to mask mouse brains? We now have a MAGeT-based masking step, but this is fairly expensive in terms of time before the lsq12/nlin stages can even start ...

gdevenyi commented 7 years ago

Interesting thought. I haven't tried it, but it should work quite well.

It would definitely be a lot less computationally expensive than MAGeT masking.

gdevenyi commented 7 years ago

The command headmask can make an otsu mask for a minc volume

gdevenyi commented 7 years ago

Here's how I'm now handling preprocessing:

#!/bin/bash
#Script to do optimal preprocessing on in-vivo/ex-vivo structural scans
#Taken using the CIC Bruker 7T
#usage:
#mouse-preprocessing-v2.sh input.mnc output.mnc

#Operations
# swaps zy to re-orient mouse
# flips x to fix left-right mixup
# centers brain in space
# denoises
# registers to DSUR_2016 atlas
# gets approximate brain mask from atlas
# expands mask
# Bias field correction with N4

set -euo pipefail
set -v
REGTARGET=${QUARANTINE_PATH}/resources/DSUR_2016/DSUR_40micron_average.mnc
REGMASK=${QUARANTINE_PATH}/resources/DSUR_2016/DSUR_40micron_mask_version2.mnc

tmpdir=$(mktemp -d)

cp $1 $tmpdir/input.mnc

input=$tmpdir/input.mnc
output=$2

#Fix "too long history" from some dcm2mnc conversion
minc_modify_header $input -sinsert :history=‘’

#Swap orientations so that scans start oriented the same way as MiCE atlases
volmash -swap zy $input $tmpdir/mash.mnc
volflip $tmpdir/mash.mnc $tmpdir/flip.mnc

#Center in space and cleanup slicing
clean_and_center_minc.pl $tmpdir/flip.mnc $tmpdir/centered.mnc

#Denoise
minc_anlm $tmpdir/centered.mnc $tmpdir/denoise.mnc

#Make a full-scan mask to short-circuit some N4 sillyness
minccalc -byte -unsigned -expression 'A[0]?1:1' $tmpdir/denoise.mnc $tmpdir/fullmask.mnc

#Register to atlas
antsRegistration --dimensionality 3 --float 0 --collapse-output-transforms 1 --verbose --minc \
--output $tmpdir/trans \
--use-histogram-matching 0 \
--transform Rigid[0.1] --metric Mattes[$tmpdir/denoise.mnc,${REGTARGET},1] --convergence [2000x2000,1e-6,10,1] --shrink-factors 6x4 --smoothing-sigmas 3x2 --masks [NULL,NULL] \
--transform Similarity[0.1] --metric Mattes[$tmpdir/denoise.mnc,${REGTARGET},1] --convergence [2000x2000,1e-6,10,1] --shrink-factors 4x2 --smoothing-sigmas 2x1 --masks [NULL,NULL] \
--transform Affine[0.1]     --metric Mattes[$tmpdir/denoise.mnc,${REGTARGET},1] --convergence [2000x2000,1e-6,10,1] --shrink-factors 2x1 --smoothing-sigmas 2x1 --masks [NULL,NULL]  \
--transform Affine[0.1]     --metric Mattes[$tmpdir/denoise.mnc,${REGTARGET},1] --convergence [2000x2000,1e-6,10,1] --shrink-factors 2x1 --smoothing-sigmas 1x0.5 --masks [NULL,${REGMASK}]

#Transform brain mask back
antsApplyTransforms -d 3 -i ${REGMASK} -o $tmpdir/mask.mnc -t $tmpdir/trans0_GenericAffine.xfm -r $tmpdir/denoise.mnc -n NearestNeighbor

#Bias field first round
N4BiasFieldCorrection -d 3 -i $tmpdir/denoise.mnc -b [20] -c [200x200x200,0.0] -w $tmpdir/mask.mnc -r 1 -x $tmpdir/fullmask.mnc -o $tmpdir/denoise.N4.mnc

#Re-regiser to atlas
antsRegistration --dimensionality 3 --float 0 --collapse-output-transforms 1 --verbose --minc \
--output $tmpdir/trans \
--use-histogram-matching 0 \
--transform Rigid[0.1] --metric Mattes[$tmpdir/denoise.N4.mnc,${REGTARGET},1] --convergence [2000x2000,1e-6,10,1] --shrink-factors 6x4 --smoothing-sigmas 3x2 --masks [NULL,NULL] \
--transform Similarity[0.1] --metric Mattes[$tmpdir/denoise.N4.mnc,${REGTARGET},1] --convergence [2000x2000,1e-6,10,1] --shrink-factors 4x2 --smoothing-sigmas 2x1 --masks [NULL,NULL] \
--transform Affine[0.1]     --metric Mattes[$tmpdir/denoise.N4.mnc,${REGTARGET},1] --convergence [2000x2000,1e-6,10,1] --shrink-factors 2x1 --smoothing-sigmas 2x1 --masks [NULL,NULL]  \
--transform Affine[0.1]     --metric Mattes[$tmpdir/denoise.N4.mnc,${REGTARGET},1] --convergence [2000x2000,1e-6,10,1] --shrink-factors 2x1 --smoothing-sigmas 1x0.5 --masks [NULL,${REGMASK}]

#Transform new mask
antsApplyTransforms -d 3 -i ${REGMASK} -o $tmpdir/mask.mnc -t $tmpdir/trans0_GenericAffine.xfm -r $tmpdir/denoise.N4.mnc -n NearestNeighbor

#Correct again
N4BiasFieldCorrection -d 3 -i $tmpdir/denoise.mnc -b [20] -c [200x200x200,0.0] -w $tmpdir/mask.mnc -r 1 -x $tmpdir/fullmask.mnc -o $output

#Remove scale/shear parts to make an lsq6 transform
xfminvert $tmpdir/trans0_GenericAffine.xfm $tmpdir/trans0_GenericAffine_inverse.xfm
param2xfm $(xfm2param $tmpdir/trans0_GenericAffine_inverse.xfm | grep -E 'scale|shear') $tmpdir/scale.xfm
xfminvert $tmpdir/scale.xfm $tmpdir/unscale.xfm
xfmconcat $tmpdir/trans0_GenericAffine_inverse.xfm $tmpdir/unscale.xfm $tmpdir/lsq6.xfm

#Output masks
cp $tmpdir/mask.mnc $(dirname $output)/$(basename $output .mnc)_mask.mnc

#Provide lsq6 aligned scans in addition to native-space scans
mincresample -tfm_input_sampling -invert_transform -transform $tmpdir/lsq6.xfm $output $(dirname $output)/$(basename $output .mnc)_lsq6.mnc
mincresample -tfm_input_sampling -invert_transform -transform $tmpdir/lsq6.xfm $tmpdir/mask.mnc $(dirname $output)/$(basename $output .mnc)_lsq6_mask.mnc

rm -rf $tmpdir