search

QSMxT / QSMxT

QSMxT is an end-to-end software toolbox for Quantitative Susceptibility Mapping.
https://qsmxt.github.io
GNU General Public License v3.0
38 stars 15 forks source link

Converting UK biobank SWI protocol to BIDS format #115

Closed karllandheer closed 6 months ago

karllandheer commented 6 months ago

Hello, thank you very much for this great package. I am interested in running it on UK biobank SWI data to generate QSM maps and I have a few questions I was hoping you could help me with.

1) Is there an easy way to convert UK biobank data into the required BIDS format for QSMxT? I tried using dicom-sort and dicom-convert, however the resulting nifti images were certainly wrong. My guess is because the individual coil data are provided by default? It's possible I did something wrong as well. A sample UK biobank SWI image dataset can be obtained here: https://biobank.ndph.ox.ac.uk/showcase/refer.cgi?id=32579

2) Is there a recommended "general purpose" algorithm that you recommend? I am very new to QSM, and as I can see there's a bunch of algorithms included in QSMxT (which is great!) which all have their own benefits and drawbacks, but in general which one would you recommend for UK biobank data?

3) How long does this recommend algorithm take to run? As you can imagine compute can become an issue when we're dealing with so many datasets.

Any help would be greatly appreciated!

karllandheer commented 6 months ago

Hello, I've been muddling my way through 1). I think the easiest way is to convert the dicoms to nifti and then grab the corresponding sum-of-squares images that are available. I then ran nifti-convert on the folder which has the two sum of squares images (both magnitude and phase). The output is a screenshot of dataset_qsmxt.csv, which I believe I have to input and then I should be able to run qsmxt on my data. I have a few questions about these input

DerivativePipeline: is this the algorithm I want to run? Sub: I assume that's the subject number, obviously 1 in my case? Ses: I assume thats session number Acq: Not sure ce: not sure rec: not sure run: not sure mod: not sure echo: is this echo number? For me this would either be 1 or 2? Or is it total number of echoes (2 for both)? flip: Excitation flip angle in degrees? inv: not sure mt: not sure part: is this magnitude or phase? It seemed to grab phase for the _ph images but not magnitude for the others. desc: description? MagneticFieldStrength: Field strength in Tesla I assume? EchoTime: Echo time in ms I assume?

Do you have a sample dataset_qsmxt.csv floating around? Sorry if this is answered somewhere else, but I'm a bit new to BIDS and it seems BIDS is less well described for SWI than T1.

Screenshot 2024-03-22 at 4 48 06 PM
astewartau commented 6 months ago

Hi @karllandheer, and thanks for your interest in using QSMxT!

1. Conversion to BIDS

You can convert UK BioBank into the BIDS format, and you are on the right track. I have found the best way to do this is by starting with the NIfTI formatted data from UK BioBank, and using nifti-convert to convert to BIDS, filling out the spreadsheet as you have started. Unfortunately, you need some expertise with the BIDS specification to know how to do this correctly - other tools can assist with BIDS conversions and are more user-friendly, but I don't know of any that reliably work with QSM data. I have documented a correct example here.

I downloaded subject 1000219 data types 20251, 20252, and 26301:

~/neurodesktop-storage/data/2023-BioBank-QSM-Data: tree 1000219/
1000219/
├── 1000219_20251_2_0
│   └── SWI
│       ├── brain_mask.nii.gz
│       ├── filtered_phase.nii.gz
│       ├── SOS_TE1.nii.gz
│       ├── SOS_TE2.nii.gz
│       ├── SWI.nii.gz
│       ├── SWI_to_T1.mat
│       ├── SWI_TOTAL_MAG.json
│       ├── SWI_TOTAL_MAG.nii.gz
│       ├── SWI_TOTAL_MAG_orig.nii.gz
│       ├── SWI_TOTAL_MAG_TE2.json
│       ├── SWI_TOTAL_MAG_TE2_orig.nii.gz
│       ├── SWI_TOTAL_MAG_to_T1.nii.gz
│       ├── SWI_TOTAL_PHA.json
│       ├── SWI_TOTAL_PHA_TE2.json
│       ├── T1_to_SWI.mat
│       ├── T2star.nii.gz
│       └── T2star_to_T1.nii.gz
├── 1000219_20252_2_0
│   └── T1
│       ├── T1_brain_mask.nii.gz
│       ├── T1_brain.nii.gz
│       ├── T1_brain_to_MNI.nii.gz
│       ├── T1_fast
│       │   ├── T1_brain_bias.nii.gz
│       │   ├── T1_brain_pve_0.nii.gz
│       │   ├── T1_brain_pve_1.nii.gz
│       │   ├── T1_brain_pve_2.nii.gz
│       │   └── T1_brain_seg.nii.gz
│       ├── T1_first
│       │   ├── T1_first_all_fast_firstseg.nii.gz
│       │   ├── T1_first-BrStem_first.bvars
│       │   ├── T1_first-BrStem_first.vtk
│       │   ├── T1_first-L_Accu_first.bvars
│       │   ├── T1_first-L_Accu_first.vtk
│       │   ├── T1_first-L_Amyg_first.bvars
│       │   ├── T1_first-L_Amyg_first.vtk
│       │   ├── T1_first-L_Caud_first.bvars
│       │   ├── T1_first-L_Caud_first.vtk
│       │   ├── T1_first-L_Hipp_first.bvars
│       │   ├── T1_first-L_Hipp_first.vtk
│       │   ├── T1_first-L_Pall_first.bvars
│       │   ├── T1_first-L_Pall_first.vtk
│       │   ├── T1_first-L_Puta_first.bvars
│       │   ├── T1_first-L_Puta_first.vtk
│       │   ├── T1_first-L_Thal_first.bvars
│       │   ├── T1_first-L_Thal_first.vtk
│       │   ├── T1_first-R_Accu_first.bvars
│       │   ├── T1_first-R_Accu_first.vtk
│       │   ├── T1_first-R_Amyg_first.bvars
│       │   ├── T1_first-R_Amyg_first.vtk
│       │   ├── T1_first-R_Caud_first.bvars
│       │   ├── T1_first-R_Caud_first.vtk
│       │   ├── T1_first-R_Hipp_first.bvars
│       │   ├── T1_first-R_Hipp_first.vtk
│       │   ├── T1_first-R_Pall_first.bvars
│       │   ├── T1_first-R_Pall_first.vtk
│       │   ├── T1_first-R_Puta_first.bvars
│       │   ├── T1_first-R_Puta_first.vtk
│       │   ├── T1_first-R_Thal_first.bvars
│       │   └── T1_first-R_Thal_first.vtk
│       ├── T1.json
│       ├── T1.nii.gz
│       ├── T1_orig_defaced.nii.gz
│       ├── T1_unbiased_brain.nii.gz
│       └── transforms
│           ├── T1_to_MNI_linear.mat
│           └── T1_to_MNI_warp_coef.nii.gz
└── 1000219_26301_2_0
    └── SWI
        └── QSM
            ├── PHASE_TE1.nii.gz
            ├── PHASE_TE2.nii.gz
            ├── QSM_CSFref.nii.gz
            ├── QSM_CSF.txt
            ├── QSM_mask.nii.gz
            ├── QSM.nii.gz
            └── QSM_ud.nii.gz

11 directories, 69 files

Then, I started the conversion to bids using nifti-convert:

~/neurodesktop-storage/data/2023-BioBank-QSM-Data: nifti-convert 1000219 bids
[INFO]: Running QSMxT 6.4.2 (linked installation)
[INFO]: Command: /home/ashley/neurodesktop-storage/miniconda3/envs/qsmxt/bin/nifti-convert 1000219 bids
[INFO]: Python interpreter: /home/ashley/neurodesktop-storage/miniconda3/envs/qsmxt/bin/python
[INFO]: Finding NIfTI files...
[INFO]: 30 NIfTI files found.
[INFO]: Extracting details from filenames...
[INFO]: Done reading details.
[INFO]: Updating details with JSON header information...
[INFO]: Done reading JSON header files.
[INFO]: RUN AGAIN AFTER ENTERING RELEVANT BIDS INFORMATION TO /home/ashley/neurodesktop-storage/data/2023-BioBank-QSM-Data/bids/dataset_qsmxt.csv.
[INFO]: Finished

This generated a spreadsheet:

image

Which I filled out as follows:

image

Running again completes the conversion:

(qsmxt) ~/neurodesktop-storage/data/2023-BioBank-QSM-Data: nifti-convert 1000219/ bids
[INFO]: Running QSMxT 6.4.2 (linked installation)
[INFO]: Command: /home/ashley/neurodesktop-storage/miniconda3/envs/qsmxt/bin/nifti-convert 1000219/ bids
[INFO]: Python interpreter: /home/ashley/neurodesktop-storage/miniconda3/envs/qsmxt/bin/python
[INFO]: CSV spreadsheet '/home/ashley/neurodesktop-storage/data/2023-BioBank-QSM-Data/bids/dataset_qsmxt.csv' found! Reading...
[INFO]: CSV spreadsheet loaded.
[INFO]: Computing new NIfTI file names and locations...
[INFO]: New NIfTI file names and locations determined.
Summary of identified files and proposed new names (following BIDS standard):
SWI_TOTAL_MAG_TE2_orig.nii.gz 
     -> sub-1000219_ses-20190507_acq-swi_echo-2_part-mag_MEGRE.nii.gz
SWI_TOTAL_MAG_orig.nii.gz 
     -> sub-1000219_ses-20190507_acq-swi_echo-1_part-mag_MEGRE.nii.gz
T1_orig_defaced.nii.gz 
     -> sub-1000219_ses-20190507_acq-t1w_T1w.nii.gz
PHASE_TE1.nii.gz 
     -> sub-1000219_ses-20190507_acq-swi_echo-1_part-phase_MEGRE.nii.gz
PHASE_TE2.nii.gz 
     -> sub-1000219_ses-20190507_acq-swi_echo-2_part-phase_MEGRE.nii.gz
QSM_mask.nii.gz 
     -> sub-1000219_ses-20190507_acq-swi_mask.nii.gz
Confirm copy + renames? (n for no): 

[INFO]: Copying NIfTI files to new locations with new names...
[INFO]: Done copying NIfTI files.
[INFO]: Copying JSON header files if present and generating them if needed...
[INFO]: Automatically generated JSON header file '/home/ashley/neurodesktop-storage/data/2023-BioBank-QSM-Data/bids/sub-1000219/ses-20190507/anat/sub-1000219_ses-20190507_acq-swi_echo-2_part-mag_MEGRE.json'
[INFO]: Automatically generated JSON header file '/home/ashley/neurodesktop-storage/data/2023-BioBank-QSM-Data/bids/sub-1000219/ses-20190507/anat/sub-1000219_ses-20190507_acq-swi_echo-1_part-mag_MEGRE.json'
[INFO]: Automatically generated JSON header file '/home/ashley/neurodesktop-storage/data/2023-BioBank-QSM-Data/bids/sub-1000219/ses-20190507/anat/sub-1000219_ses-20190507_acq-t1w_T1w.json'
[INFO]: Automatically generated JSON header file '/home/ashley/neurodesktop-storage/data/2023-BioBank-QSM-Data/bids/sub-1000219/ses-20190507/anat/sub-1000219_ses-20190507_acq-swi_echo-1_part-phase_MEGRE.json'
[INFO]: Automatically generated JSON header file '/home/ashley/neurodesktop-storage/data/2023-BioBank-QSM-Data/bids/sub-1000219/ses-20190507/anat/sub-1000219_ses-20190507_acq-swi_echo-2_part-phase_MEGRE.json'
[INFO]: Automatically generated JSON header file '/home/ashley/neurodesktop-storage/data/2023-BioBank-QSM-Data/bids/derivatives/biobank/sub-1000219/ses-20190507/anat/sub-1000219_ses-20190507_acq-swi_mask.json'
[INFO]: Done copying and generating JSON header files.
[INFO]: Generating details for BIDS datset_description.json...
[INFO]: Writing BIDS dataset_description.json...
[INFO]: Writing BIDS .bidsignore file...
[INFO]: Writing BIDS dataset README...
[INFO]: Finished

With the final folder appearing as follows:

(qsmxt) ~/neurodesktop-storage/data/2023-BioBank-QSM-Data: tree bids
bids
├── dataset_description.json
├── dataset_qsmxt.csv
├── derivatives
│   └── biobank
│       └── sub-1000219
│           └── ses-20190507
│               └── anat
│                   ├── sub-1000219_ses-20190507_acq-swi_mask.json
│                   └── sub-1000219_ses-20190507_acq-swi_mask.nii.gz
├── log_2024-03-25_16-53-20248899.txt
├── log_2024-03-25_16-53-33292075.txt
├── README
├── references.txt
└── sub-1000219
    └── ses-20190507
        └── anat
            ├── sub-1000219_ses-20190507_acq-swi_echo-1_part-mag_MEGRE.json
            ├── sub-1000219_ses-20190507_acq-swi_echo-1_part-mag_MEGRE.nii.gz
            ├── sub-1000219_ses-20190507_acq-swi_echo-1_part-phase_MEGRE.json
            ├── sub-1000219_ses-20190507_acq-swi_echo-1_part-phase_MEGRE.nii.gz
            ├── sub-1000219_ses-20190507_acq-swi_echo-2_part-mag_MEGRE.json
            ├── sub-1000219_ses-20190507_acq-swi_echo-2_part-mag_MEGRE.nii.gz
            ├── sub-1000219_ses-20190507_acq-swi_echo-2_part-phase_MEGRE.json
            ├── sub-1000219_ses-20190507_acq-swi_echo-2_part-phase_MEGRE.nii.gz
            ├── sub-1000219_ses-20190507_acq-t1w_T1w.json
            └── sub-1000219_ses-20190507_acq-t1w_T1w.nii.gz

9 directories, 18 files

2. Recommended algorithm

QSMxT has several pipelines you can run, each applying different collections of algorithms and parameter values. I would usually recommend using the default pipeline (GRE) to start with and branching out depending on the data - e.g. the EPI pipeline tends to deal with the lower SNR of accelerated acquisitions better than the GRE pipeline; the body pipeline was optimised for a prostate imaging application and doesn't assume brain anatomy etc.

After some limited testing using the BioBank data, I found that the result improved when we added a BET mask to the default threshold-based masking algorithm, which you can do by applying the --add_bet flag (or enabling the equivalent setting using the interactive menus). I used the following command:

qsmxt bids qsm --add_bet --auto_yes

image

This result is comparable to the BioBank images, with fewer artefacts around veins. Importantly, the scale of QSMxT's images is parts-per-million (ppm), whereas the BioBank data uses parts-per-billion (ppb).

3. Runtime

For one subject, I just timed this at 11 minutes on my workstation (11th Gen Intel(R) Core(TM) i5-11500 @ 2.70GHz). However, if you have access to a high-performance cluster (HPC) with SLURM or PBS, QSMxT can submit each algorithmic step as a job, parallelising all of the processing and making this much more feasible to run in a short amount of time. To do this, follow the HPC install instructions here, then when you run, use the --slurm or --pbs flag. For SLURM, use:

qsmxt bids qsm --add_bet --auto_yes --slurm ACCOUNT_STRING SLURM_QUEUE

Replacing ACCOUNT_STRING with the account string or group with SLURM permissions enabled and SLURM_QUEUE with the queue name to which jobs should be submitted. For PBS, it is similar:

qsmxt bids qsm --add_bet --auto_yes --pbs ACCOUNT_STRING

Hope that helps! Please let me know how you go so I can close the issue once resolved. :smile:

karllandheer commented 6 months ago

Hello, thanks for your response! We do use a cluster, however compute costs can become quite relevant depending on length of modules. 10 minutes isn't such a major concern. One quick thing, I am getting this error:

Screenshot 2024-03-25 at 8 50 12 AM

my guess is you are parsing by something like string.split('.')[0]? My user name has a . in the path

Screenshot 2024-03-25 at 8 51 09 AM

anyways I am going to try a workaround, but this may come up for others.. i'll get back to you

karllandheer commented 6 months ago

Hello, yes that was indeed the issue. I am having another issue now. Working on figuring out what the issue is... Here's the error log:

(UKB) karl.landheer@C02F203JMD6R UKB_Brain_Pipeline % qsmxt /Users/karl.landheer/Downloads/bids_output/ /Users/karl.landheer/Downloads/qsm_output --auto_yes
[INFO]: Loading previous QSMxT settings... [INFO]: QSMxT v6.4.2 [INFO]: Python interpreter: /Users/karl.landheer/opt/anaconda3/envs/ukb/bin/python [INFO]: Command: qsmxt /Users/karl.landheer/Downloads/bids_output /Users/karl.landheer/Downloads/qsm_output --auto_yes [INFO]: Creating QSMxT workflow for sub-1.ses-1.acq-swi... [INFO]: Running using MultiProc plugin with n_procs=12 240325-10:12:29,589 nipype.workflow INFO: Workflow workflow settings: ['check', 'execution', 'logging', 'monitoring'] 240325-10:12:29,644 nipype.workflow INFO: Running in parallel. 240325-10:12:29,648 nipype.workflow INFO: [MultiProc] Running 0 tasks, and 4 jobs ready. Free memory (GB): 14.40/14.40, Free processors: 12/12. 240325-10:12:31,652 nipype.workflow INFO: [MultiProc] Running 3 tasks, and 2 jobs ready. Free memory (GB): 13.80/14.40, Free processors: 9/12. Currently running:

Cmdline: /Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/qsmxt/scripts/romeo_unwrapping.jl --no-rescale --compute-B0 B0.nii --phase /Users/karl.landheer/Downloads/qsm_output/workflow/sub-1/ses-1/qsmxt_acq-swi/mrt_romeo_combine-phase/multi-echo-phase.nii --mag /Users/karl.landheer/Downloads/qsm_output/workflow/sub-1/ses-1/qsmxt_acq-swi/mrt_romeo_combine-phase/multi-echo-mag.nii -t '[9.42 19.7]' --phase-offset-correction monopolar Stdout:

Stderr: env: julia: No such file or directory Traceback: Traceback (most recent call last): File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nibabel/loadsave.py", line 100, in load stat_result = os.stat(filename) FileNotFoundError: [Errno 2] No such file or directory: 'unwrapped.nii'

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/interfaces/base/core.py", line 400, in run
    outputs = self.aggregate_outputs(runtime)
  File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/interfaces/base/core.py", line 429, in aggregate_outputs
    predicted_outputs = self._list_outputs()  # Predictions from _list_outputs
  File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/qsmxt/interfaces/nipype_interface_romeo.py", line 97, in _list_outputs
    outputs['phase_unwrapped'] = split_multi_echo("unwrapped.nii", [extend_fname(f, "_romeo-unwrapped", ext="nii", out_dir=os.getcwd()) for f in self.inputs.phase])
  File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/qsmxt/interfaces/nipype_interface_romeo.py", line 20, in split_multi_echo
    image4d = nib.load(in_path).get_fdata()
  File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nibabel/loadsave.py", line 102, in load
    raise FileNotFoundError(f"No such file or no access: '{filename}'")
FileNotFoundError: No such file or no access: 'unwrapped.nii'

240325-10:13:01,705 nipype.workflow ERROR: Node romeo-voxelquality failed to run on host C02F203JMD6R. 240325-10:13:01,706 nipype.workflow ERROR: Saving crash info to /Users/karl.landheer/PythonProjects/RGCi_Repo/rgc-imaging/other_packages/UKB_Brain_Pipeline/crash-20240325-101301-karl.landheer-romeo-voxelquality-db1a0fac-8489-4718-b395-6d4fd0b6e1ad.pklz Traceback (most recent call last): File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/pipeline/plugins/multiproc.py", line 67, in run_node result["result"] = node.run(updatehash=updatehash) File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/pipeline/engine/nodes.py", line 527, in run result = self._run_interface(execute=True) File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/pipeline/engine/nodes.py", line 645, in _run_interface return self._run_command(execute) File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/pipeline/engine/nodes.py", line 771, in _run_command raise NodeExecutionError(msg) nipype.pipeline.engine.nodes.NodeExecutionError: Exception raised while executing Node romeo-voxelquality.

Cmdline: /Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/qsmxt/scripts/romeo_voxelquality.jl --TEs '[0.00942]' --mag /Users/karl.landheer/Downloads/qsm_output/workflow/sub-1/ses-1/qsmxt_acq-swi/mask_workflow/nibabal-numpy_combine-magnitude/sub-1_ses-1_acq-swi_echo-1_part-mag_MEGRE_resampled_combined.nii.gz --phase /Users/karl.landheer/Downloads/qsm_output/workflow/sub-1/ses-1/qsmxt_acq-swi/nibabel_numpy_nilearn_axial-resampling/mapflow/_nibabel_numpy_nilearn_axial-resampling0/sub-1_ses-1_acq-swi_echo-1_part-phase_MEGRE_scaled_resampled.nii.gz --output sub-1_ses-1_acq-swi_echo-1_part-phase_MEGRE_scaled_resampled_romeo_voxelquality.nii --type grad+second+mag_weight+mag_coherence Stdout:

Stderr: env: julia: No such file or directory Traceback: RuntimeError: subprocess exited with code 127.

240325-10:13:01,709 nipype.workflow INFO: [MultiProc] Running 0 tasks, and 0 jobs ready. Free memory (GB): 14.40/14.40, Free processors: 12/12. 240325-10:13:03,700 nipype.workflow INFO:

240325-10:13:03,701 nipype.workflow ERROR: could not run node: workflow.sub-1.ses-1.qsmxt_acq-swi.mrt_romeo_combine-phase 240325-10:13:03,701 nipype.workflow INFO: crashfile: /Users/karl.landheer/PythonProjects/RGCi_Repo/rgc-imaging/other_packages/UKB_Brain_Pipeline/crash-20240325-101301-karl.landheer-mrt_romeo_combine-phase-ef840f93-8243-4e1c-b6b1-48b222c0e46c.pklz 240325-10:13:03,702 nipype.workflow ERROR: could not run node: workflow.sub-1.ses-1.qsmxt_acq-swi.mask_workflow.romeo-voxelquality 240325-10:13:03,703 nipype.workflow INFO: crashfile: /Users/karl.landheer/PythonProjects/RGCi_Repo/rgc-imaging/other_packages/UKB_Brain_Pipeline/crash-20240325-101301-karl.landheer-romeo-voxelquality-db1a0fac-8489-4718-b395-6d4fd0b6e1ad.pklz 240325-10:13:03,704 nipype.workflow INFO:

RuntimeError: Traceback (most recent call last): File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/pipeline/plugins/multiproc.py", line 67, in run_node result["result"] = node.run(updatehash=updatehash) File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/pipeline/engine/nodes.py", line 527, in run result = self._run_interface(execute=True) File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/pipeline/engine/nodes.py", line 645, in _run_interface return self._run_command(execute) File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/pipeline/engine/nodes.py", line 771, in _run_command raise NodeExecutionError(msg) nipype.pipeline.engine.nodes.NodeExecutionError: Exception raised while executing Node mrt_romeo_combine-phase.

Cmdline: /Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/qsmxt/scripts/romeo_unwrapping.jl --no-rescale --compute-B0 B0.nii --phase /Users/karl.landheer/Downloads/qsm_output/workflow/sub-1/ses-1/qsmxt_acq-swi/mrt_romeo_combine-phase/multi-echo-phase.nii --mag /Users/karl.landheer/Downloads/qsm_output/workflow/sub-1/ses-1/qsmxt_acq-swi/mrt_romeo_combine-phase/multi-echo-mag.nii -t '[9.42 19.7]' --phase-offset-correction monopolar Stdout:

Stderr: env: julia: No such file or directory Traceback: Traceback (most recent call last): File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nibabel/loadsave.py", line 100, in load stat_result = os.stat(filename) FileNotFoundError: [Errno 2] No such file or directory: 'unwrapped.nii'

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/interfaces/base/core.py", line 400, in run
    outputs = self.aggregate_outputs(runtime)
  File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/interfaces/base/core.py", line 429, in aggregate_outputs
    predicted_outputs = self._list_outputs()  # Predictions from _list_outputs
  File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/qsmxt/interfaces/nipype_interface_romeo.py", line 97, in _list_outputs
    outputs['phase_unwrapped'] = split_multi_echo("unwrapped.nii", [extend_fname(f, "_romeo-unwrapped", ext="nii", out_dir=os.getcwd()) for f in self.inputs.phase])
  File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/qsmxt/interfaces/nipype_interface_romeo.py", line 20, in split_multi_echo
    image4d = nib.load(in_path).get_fdata()
  File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nibabel/loadsave.py", line 102, in load
    raise FileNotFoundError(f"No such file or no access: '{filename}'")
FileNotFoundError: No such file or no access: 'unwrapped.nii'

The above exception was the direct cause of the following exception:

Traceback (most recent call last): File "/Users/karl.landheer/opt/anaconda3/envs/UKB/bin/qsmxt", line 8, in sys.exit(main()) File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/qsmxt/cli/main.py", line 1561, in main wf.run( File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/pipeline/engine/workflows.py", line 638, in run runner.run(execgraph, updatehash=updatehash, config=self.config) File "/Users/karl.landheer/opt/anaconda3/envs/ukb/lib/python3.8/site-packages/nipype/pipeline/plugins/base.py", line 224, in run raise error from cause RuntimeError: 2 raised. Re-raising first.

karllandheer commented 6 months ago

I think I am just missing Julia...

karllandheer commented 6 months ago

Hello, I had issues with it even after I installed Julia 1.9, but anyways it worked fine from the docker, I think. This is an example output, does this look right?

image

It looks like there may have been an issue with the add_bet part? I have the brain masks from a separate part of the pipeline, is there a way to add it in without relying on the bet tool? Would manually masking the mag/phase images by the brain mask work, or is it not that simple?

As an aside, will you be in Singapore? I would be interested to learn more about qsmxt

karllandheer commented 6 months ago

I'm not so sure this is right. One thing is that I am working from the raw dicom data, not the data you suggested (the reason being those are already processed and many are missing). The raw data is available here: https://biobank.ndph.ox.ac.uk/showcase/refer.cgi?id=32579

If you unpack it it produces an image (magnitude and phase) for each of the 32 coils. This data is then processed into the PHASE_TE1 / PHASE_TE2, etc, that you used by a pipeline. This pipeline does things like coil combination, phase unwrapping via fsl, etc. I could implement this pipeline, but I am wondering if qsmxt is able to handle multi coil data, or if you have such a tool available? Or is it expected that the images should be in this coil combined and phase-unwrapped form?

astewartau commented 6 months ago

[Errno 13] Permission denied: 'Users/karl.json

Thanks very much for flagging this! It is a bug in QSMxT, which I have just resolved in the main branch. The fix will be included in v6.4.3, which will be released soon. I'm glad you found a workaround for now.

env: julia: No such file or directory

As you have discovered, it looks like you were missing Julia. Indeed, I recommend using the containerised version of QSMxT to avoid these dependency problems. In case you are interested, there is a detailed list of requirements in the container build script, including a series of Julia packages described alongside it, which may explain whatever the additional problem you encountered was (possibly you were missing these Julia packages).

Output

Your output indeed does not look perfect and has some intrusive artefacts, unlike the example that I processed earlier. What was the exact command that you used to generate it? Could you provide some of the output logs to confirm the version etc.? If it is exactly the same as how I processed the other example, this dataset may be an outlier or have different characteristics. I would consider trying --premade bet to see if the result is any better - this pipeline uses a more straightforward masking algorithm that depends less on threshold parameter tuning but does not include the two-pass artefact reduction algorithm.

In case it provides any further insight; this is what my input magnitude and phase images look like from my BIDS directory:

image

Coil combination

Unfortunately, while the DICOM images include combined and uncombined magnitude and phase data, the phase data is incorrectly combined. So you are right that you would need to run a coil combination algorithm separately on the uncombined images since QSMxT expects already-combined images.

Including coil combination in QSMxT is something I have considered, but there is also no standard way to store uncombined coil data in BIDS as inputs, although an extension to support it was proposed which we could use.

If you want to run a separate coil combination, you could also consider ASPIRE/MCPC-3D-S, which also conveniently runs in Julia (see korbinian90/MriResearchTools.jl). You could then convert the output to BIDS to run in QSMxT. This is what I am considering integrating into QSMxT at some stage.

About the NIfTI data I suggested - when you say that many are missing, did you look into both the 26301 and 20251 data fields (since both are required)? Or do many subjects not have this?

astewartau commented 6 months ago

Oh, and I forgot to mention that I will be in Singapore! I am presenting work about evaluating QSM and ideas for a future 'always-online' QSM challenge (prototype here). It would be great to meet you. :smile:

karllandheer commented 6 months ago

Hello again, I think I solved and found the issue. I was using the phase maps from the scanner, but I reconstructed them using the MCPC-3D-S coil combination method I found in the literature. This is an example of the resulting QSM maps now:

Screenshot 2024-03-26 at 2 28 58 PM

Does this look right? I think so but not sure. A few remaining quick questions:

1) The scale seems to be substantially different from the values from the UKB pipeline. Perhaps this is just a difference in the units. For example the maps I have from qsmxt seem have a value in the csf of ~0.015. The UKB QSM pipeline has a value of ~90. Is the values on my end correct or am I still incorrectly specifying something?

2) Do you have a poster or talk at ISMRM? If so please let me know the time and I will do my best to make it.

astewartau commented 6 months ago

That updated reconstruction looks a lot better! This subject appears to have larger ventricles - the larger amount of free fluid may lead to noisier phase values which different algorithms may respond to differently depending on settings such as regularisation. Overall, this reconstruction looks much more reasonable.

  1. The scale of QSMxT's images is parts-per-million (ppm), whereas the BioBank data uses parts-per-billion (ppb). You can see how they scale it here. In the same units, QSMxT is producing 15 vs. UKB's 90 for CSF, so there is still some difference. By default, QSMxT also references QSM values with respect to the mean susceptibility value across the brain. Interestingly, UKB is supposed to have referenced QSM values to CSF, so we would expect the CSF value to be closer to zero in theirs... If you used QSMxT's pipeline for segmentation and analysis, which is produced using FastSurfer for segmentation and ANTs for registration, this may also have an effect since UKB used FSL's FAST for segmentation and FIRST for registration. That is not to say that any particular result is more correct, but that this could explain some of the differences in addition to the trade-offs that different QSM algorithms make.

  2. I have a poster presentation at ISMRM titled "QSM-CI: An automated continuous QSM challenge" (abstract #6388). It will be during the 'Quantitative Susceptibility Mapping' poster session on the 6th of May at 9:15 AM. I also plan to attend the EMTP study group session, which I may have the opportunity to present at, depending on available spots.

karllandheer commented 6 months ago
  1. Hello, OK great. Regarding the UKB pipeline I was looking at their QSM image, not QSM_CSFref (which subtracts the value form CSF to make susceptibility in CSF = 0). And good to know about the ppm vs ppb. I am using the same FAST and FIRST registration as UKB (at least initially), but I do have FastSurfer segmentations so will probably use those as well later.
  2. OK great I will do my best to be there! 9:15 on singapore time will be hard though :)
astewartau commented 6 months ago

Wonderful! 😊 I will close this issue for now and reopen if you have further concerns on this topic.