Sign inSign up

diannepat/dwiqc

By diannepat

•Updated over 6 years ago

FSL 6 dwi processing through eddy, followed by qc metrics. ~10 GB based on fsl6-core-nogpu

Image
0

10K+

diannepat/dwiqc repository overview

⁠DWIQC PROJECT

Dianne Patterson, University of Arizona, SLHS Dept.
Created April 20, 2019
Updated October 10, 2020

Special thanks to Patricia Klobusiakova without whose eagle eye and sharp mind, there would be many more errors in this set of scripts.

WARNING: This works for non-session data. It requires creating two text files for the dwi acq parameters. These need to be placed at the top level of the BIDS directory so that they apply to all subjects.

TO DO: Calculate all text files automatically so the user does not have to create them manually. Handle session data.

⁠REVISION CONTROL

The dwiqc project is under version control using git.
git clone https://[email protected]/dpat/dwiqc.git⁠ to create a copy of the dwiqc directory.

⁠WHAT IS DWIQC?

The scripts primarily run FSL 6 commands (these should be backward compatible to ~5.0.7). Everything was written to run in an Ubuntu 18 docker container using BIDS compliant file structure and naming: http://bids.neuroimaging.io/⁠

  • DWIQC creates separate derivatives directories for:
    • dwi (fsl_dwi_proc),
    • anatomical processing (fsl_anat_proc),
    • If you use diannepat/bip2 to run further processing on your dwi data, you'll find it uses the same directory structure and interoperates with dwiqc to avoid redundant processing or files.

⁠DOCKER

The dwiqc directory contains a Dockerfile and run.py at the top level. If you have docker running on Mac or Linux, the easiest way to get the current container is to pull it from dockerhub like this:

docker pull diannepat/dwiqc

⁠dwiqc_wrap.sh

To facilitate running dwiqc on your local machine, use dwiqc_wrap.sh.

dwiqc_wrap.sh

If dwiqc_wrap.sh is in your path, the above call will show you help.

By examining the contents of dwiqc_wrap.sh, you can see examples of the different calls that can be made to the docker container.

You must be in your main bids directory for dwiqc_wrap.sh to work correctly.

⁠SINGULARITY

To run dwiqc on the HPC, you must build a singularity container. Note that FSL 6.X is ~18 GB, so you'll want to make sure you have enough space to build it.

singularity build dwiqc.sif docker://diannepat/dwiqc

singularity run ./dwiqc.sif ${PWD}/Data ${PWD}/Data/derivatives participant --participant_label 327

This command says to run dwiqc.sif (the singularity container) on subjects under the present working dir in Data and write to Data/derivatives. We are running for one subject, sub-327. Look at the example run commands in dwiqc_wrap.sh for other examples. Run commands are a bit different for Singularity and Docker, but a lot of the BIDS structures are the same.

⁠Sample Dataset

https://arizona.box.com/v/bip⁠

⁠PARTICIPANT LEVEL STEPS

Run for selected subjects (e.g., 327 is sub-327 etc.):

dwiqc_wrap.sh QC 327 328

Run for all available subjects:

dwiqc_wrap.sh QC_All

⁠GROUP LEVEL STEP

Run group level reports (depends on individual reports)

dwiqc_wrap.sh QC_Group

Outputs will be placed under derivatives/fsl_dwi_proc

⁠GETTING READY

Your file naming and directory structure are presumed to be bids compliant and to include fieldmaps. Fieldmaps are notoriously difficult to use consistently. These examples work with our Siemens Skyra at U of A:

| acqp_AP.txt
| acqparams.txt
| index_AP.txt
| pe.txt
|-- sub-001
|   |-- anat
|   |   |-- sub-001_T1w.json
|   |   `-- sub-001_T1w.nii.gz
|   |-- dwi
|   |   |-- sub-001_acq-AP_dwi.bval
|   |   |-- sub-001_acq-AP_dwi.bvec
|   |   |-- sub-001_acq-AP_dwi.json
|   |   `-- sub-001_acq-AP_dwi.nii.gz
|   `-- fmap
|       |-- sub-001_dir-PA_epi.json
|       |-- sub-001_dir-PA_epi.nii.gz
|       |-- sub-001_magnitude1.json
|       |-- sub-001_magnitude1.nii.gz
|       |-- sub-001_magnitude2.json
|       |-- sub-001_magnitude2.nii.gz
|       |-- sub-001_phasediff.json
|       `-- sub-001_phasediff.nii.gz
⁠dwiqc runs topup, applytopup, eddy and eddyqc.

The function Top runs topup and applytopup. This takes ~25 minutes to run on a 2013 mac pro. topup calculates the distortion and movement in the DWI images. Top looks for ${bids_dir}/acqparams.txt

applytopup applies distortion and movement corrections to the B0 images. applytopup uses two images: blip1 consists of two B0 blip-up volumes. Call these volumes 1 and 2. blip2 consists of two B0 blip-down volumes. Call these volumes 3 and 4. applytopup combines the first blip1 volume (1) with the first blip2 volume (3) to create a volume (1+3). Then the 2nd blip1 image (2) is combined with the second blip2 image (4) to create a volume (2+4). These two new volumes are stored in ${output_dir}/${subj}/dwi/b0_hifi.nii.gz.

eddy expects four parameter files in the bids_dir: acqp_AP.txt acqparams.txt index_AP.txt pe.txt These are explained in more detail below.

You will need a structural and a dwi image with bvals and bvecs. They should be named and organized in the BIDS data structure.

You also need standard Siemens fieldmaps (magnitude and phase images and reverse phase encoded B0 image), place these in the fmap directory. These fieldmaps are used with topup and eddy to improve distortion correction.

Reverse Phase encode B0 images: My experiments suggest it is better if the main dwi image is A-P data rather than P-A data (less noise, fewer spikes and holes in weird places, better eyeball reconstruction). Thus, the reverse phase encode B0 images should be P-A.

You will need the following four text files in your bids_dir:

  • acqp_??.txt (e.g., acqp_AP.txt or acqp.txt)
  • acqparams.txt
  • index_??.txt (e.g., index_AP.txt or index.txt)
  • pe.txt

dwiqc will calculate pe if there is no pe.txt and it will create index.txt if it fails to find index*.txt.

For my 32 direction dataset, the files look like this:

⁠acq_AP.txt:

0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602
0 -1 0 0.0602

⁠acqparams.txt:

0 -1 0 0.0602
0 -1 0 0.0602
0 1 0 0.0602
0 1 0 0.0602

⁠index_AP.txt (row or column):

1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1

⁠pe.txt:

-y

⁠CREATING THE 2 Text FILES

These should be placed at the top-level of the BIDS directory to be accessible to all subject runs.

⁠acq_AP.txt and acqparams.txt:
⁠TOTAL READOUT TIME OF DWI (0.0602)

Total Readout time is used by topup, and appears in the acqp and acqparams files. There are several ways to compute total readout time.

⁠The Topup Users guide says:

"If your readout time is identical for all acquisitions you don't necessarily have to specify a valid value in this column (you can e.g. just set it to 1), but if you do specify correct values the estimated field will be correctly scaled in Hz, which may be a useful sanity check." https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/topup/TopupUsersGuide/⁠

FSL calculates total readout time as the echo spacing (dwell time) * the number of phase encoding steps. If you use grappa (may be called "ParallelReductionFactorInPlane" or "AccelFactPE" in the JSON file), divide the number of phase encoding steps (128 for sub-001) by the grappa factor (2 for sub-001). "PhaseEncodingLines": 128 "ParallelReductionFactorInPlane": 2 (n.b. sometimes called "AccelFactPE": 2) "DerivedVendorReportedEchoSpacing": 0.000940019,

e.g., (128/2)*.000940019 sec =0.0602 sec

From the JSON file, we see that dcm2niix has performed a similar (but apparently slightly different) calculation for sub-001: "TotalReadoutTime": 0.0596912

I use the calculated value (but either is probably fine)

0 -1 0 0.0602 (acqp_A-P.txt...this line is repeated 32 times, once for each DWI volume)

acqparams.txt: Used for topup. Specifies the values for the 4 blip_down and blip_up B0 volumes

acqp_??.txt (e.g., acqp_AP.txt or acqp.txt): The acqp file is used by eddy, as described here: http://fsl.fmrib.ox.ac.uk/fsl/fslwiki/EDDY/UsersGuide⁠. The format of this file is identical to that used by topup (though the parameter is called --datain there). The value is repeated 32 times for the 32 volumes. The -1 indicated the phase encode direction for y in A-P. The last value, Total Readout time of DWI, is described above.

Keep in mind, that if you have 2 BO images in your main acquisition, these will be concatenated with the B0 images in your reverse-phase encode image to sum to four B0 images.

⁠The following files no longer have to be created

⁠index_A-P.txt (Now created automatically, not explicitly needed)

Used by eddy to index all the values in acqp_A-P.txt (all 1's)

⁠pe.txt (Now calculated automatically from the JSON file, not explicitly needed)
⁠PHASE ENCODE DIRECTION OF DWI

epi_reg requires a phase encode direction (--pedir) for the primary dwi image if you have fmap. We read ${bids_dir}/pe.txt to get this value (e.g., -y). Unfortunately, the JSON files use i,j,k and FSL uses x,y,z to record this value. "PhaseEncodingDirection": "j-" in the json file corresponds to -y in FSL's terminology (-y=AP).

Tag summary

Content type

Image

Digest

Size

3.2 GB

Last updated

over 6 years ago

docker pull diannepat/dwiqc