--distortion-group-merge concat fails with DimensionError while average works

Summary

When using two BIDS-formatted DWI acquisitions with opposite phase encoding directions (AP and PA), I encounter an error whenever I set --distortion-group-merge concat. However, if I set it to average, QSIPrep runs successfully. This seems unexpected and may represent a bug.

Additional details

  • QSIPrep version: 1.0.1
  • Docker version: 28.2.2, build e6534b4

What were you trying to do?

I have BIDS data structured as follows:

sub-01/
├── anat/
│   ├── sub-01_ses-anat_T1w.nii.gz
│   ├── sub-01_ses-anat_T1w.json
│   ├── sub-01_ses-anat_T2w.nii.gz
│   └── sub-01_ses-anat_T2w.json
├── dwi/
│   ├── sub-01_ses-anat_acq-a_dir-AP_dwi.nii.gz
│   ├── sub-01_ses-anat_acq-a_dir-AP_dwi.bval
│   ├── sub-01_ses-anat_acq-a_dir-AP_dwi.bvec
│   ├── sub-01_ses-anat_acq-a_dir-AP_dwi.json
│   ├── sub-01_ses-anat_acq-a_dir-PA_dwi.nii.gz
│   ├── sub-01_ses-anat_acq-a_dir-PA_dwi.bval
│   ├── sub-01_ses-anat_acq-a_dir-PA_dwi.bvec
│   ├── sub-01_ses-anat_acq-a_dir-PA_dwi.json
│   ├── sub-01_ses-anat_acq-b_dir-AP_dwi.nii.gz
│   ├── sub-01_ses-anat_acq-b_dir-AP_dwi.bval
│   ├── sub-01_ses-anat_acq-b_dir-AP_dwi.bvec
│   ├── sub-01_ses-anat_acq-b_dir-AP_dwi.json
│   ├── sub-01_ses-anat_acq-b_dir-PA_dwi.nii.gz
│   ├── sub-01_ses-anat_acq-b_dir-PA_dwi.bval
│   ├── sub-01_ses-anat_acq-b_dir-PA_dwi.bvec
│   └── sub-01_ses-anat_acq-b_dir-PA_dwi.json
└── sub-01_ses-anat_scans.tsv

Using --bids-filter-file, I selected only acq-b. Acquisition details:
• 160 volumes, 84 slices
• 1.798 mm isotropic voxels
• TR/TE = 3.2 s / 81 ms
• Multiband MB=4
• EffectiveEchoSpacing = 0.650 ms
• TotalReadoutTime = 73.45 ms
• Both AP and PA acquisitions use the same b-value tables:

0 2000 3000 2000 1000 3000 3000 2000 3000 1000 2000 3000 2000 1000 500 3000 3000 0 2000 3000 1000 2000 3000 3000 1000 2000 3000 2000 3000 1000 2000 3000 3000 1000 0 2000 3000 2000 3000 1000 3000 500 2000 3000 1000 2000 3000 2000 3000 1000 3000 0 2000 3000 1000 2000 3000 3000 2000 1000 3000 2000 3000 1000 2000 3000 500 3000 0 2000 1000 3000 2000 3000 1000 3000 2000 3000 2000 1000 3000 2000 3000 1000 3000 0 2000 3000 2000 1000 3000 3000 2000 500 1000 3000 2000 3000 2000 1000 3000 3000 0 2000 1000 3000 2000 3000 3000 1000 2000 3000 2000 1000 3000 2000 3000 3000 500 0 1000 2000 3000 2000 1000 3000 3000 2000 3000 1000 2000 3000 2000 1000 3000 3000 0 2000 3000 1000 2000 3000 3000 1000 2000 500 3000 2000 3000 1000 2000 3000 2000 0 1000 3000 2000 3000 1000 2000

Command:

docker run --rm \
  --platform linux/amd64 \
  --privileged \
  --gpus all \
  -v "${bids_root_dir}/Rawdata":/input:ro \
  -v "${bids_root_dir}/Rawdata/derivatives/qsiprep":/output \
  -v "${bids_root_dir}/work":/work \
  -v "${bids_root_dir}/qsiprep/license.txt":/license.txt:ro \
  -v "${bids_root_dir}/qsiprep/filters_b.json":/filters_b.json:ro \
  -v "${bids_root_dir}/qsiprep/eddy_params.json":/eddy_params.json:ro \
  pennlinc/qsiprep:1.0.1 \
    /input /output participant \
    --fs-license-file /license.txt \
    --skip-bids-validation \
    --stop-on-first-crash \
    --participant-label 01 \
    --session-id anat \
    --distortion-group-merge concat \
    --bids-filter-file /filters_b.json \
    --eddy-config /eddy_params.json \
    --output-resolution 1.8 \
    --pepolar-method TOPUP+DRBUDDI \
    --unringing-method rpg \
    -w /work

What did you expect to happen?

I expected qsiprep to concatenate sub-01_ses-anat_acq-b_dir-AP_dwi.nii.gz and sub-01_ses-anat_acq-b_dir-PA_dwi.nii.gz, so that the reconstruction reflects the entire AP/PA series.

What actually happened?

Error log (excerpt):

250819-12:58:43,306 nipype.workflow INFO:
         [Node] Setting-up "qsiprep_1_0_wf.sub_03_ses_anat_wf.sub_03_ses_anat_acq_b_final_merge_wf.distortion_merger" in "/work/qsiprep_1_0_wf/sub_03_ses_anat_wf/sub_03_ses_anat_acq_b_final_merge_wf/distortion_merger".
250819-12:58:43,313 nipype.workflow INFO:
         [Node] Executing "distortion_merger" <qsiprep.interfaces.dwi_merge.MergeDWIs>
250819-12:58:43,381 nipype.workflow INFO:
         [Node] Finished "distortion_merger", elapsed time 0.067441s.
250819-12:58:43,382 nipype.workflow WARNING:
         Storing result file without outputs
250819-12:58:43,382 nipype.workflow WARNING:
         [Node] Error on "qsiprep_1_0_wf.sub_03_ses_anat_wf.sub_03_ses_anat_acq_b_final_merge_wf.distortion_merger" (/work/qsiprep_1_0_wf/sub_03_ses_anat_wf/sub_03_ses_anat_acq_b_final_merge_wf/distortion_merger)
250819-12:58:44,473 nipype.workflow ERROR:
         Node distortion_merger failed to run on host 0d5009a5e82a.

        Traceback (most recent call last):
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nipype/interfaces/base/core.py", line 401, in run
            runtime = self._run_interface(runtime)
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/qsiprep/interfaces/dwi_merge.py", line 76, in _run_interface
            to_concat, b0_means, corrections = harmonize_b0s(
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/qsiprep/interfaces/dwi_merge.py", line 700, in harmonize_b0s
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nipype/interfaces/base/core.py", line 401, in run
            runtime = self._run_interface(runtime)
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/qsiprep/interfaces/dwi_merge.py", line 76, in _run_interface
            to_concat, b0_means, corrections = harmonize_b0s(
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/qsiprep/interfaces/dwi_merge.py", line 700, in harmonize_b0s
            b0_mean = index_img(dwi_nii, b0_indices).get_fdata().mean()
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nilearn/image/image.py", line 665, in index_img
            imgs = check_niimg_4d(imgs)
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nilearn/_utils/niimg_conversions.py", line 412, in check_niimg_4d
            return check_niimg(
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nilearn/_utils/niimg_conversions.py", line 329, in check_niimg
            raise DimensionError(len(niimg.shape), ensure_ndim)
        nilearn._utils.exceptions.DimensionError: Input data has incompatible dimensionality: Expected dimension is 4D and you provided a 3D image. See https://nilearn.github.io/stable/manipulating_images/input_output.html.

250819-12:58:44,944 nipype.workflow CRITICAL:
         QSIPrep failed: Traceback (most recent call last):
  File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nipype/pipeline/plugins/multiproc.py", line 66, in run_node
    result["result"] = node.run(updatehash=updatehash)
  File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nipype/pipeline/engine/nodes.py", line 525, in run
    result = self._run_interface(execute=True)
  File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nipype/pipeline/engine/nodes.py", line 643, in _run_interface
    return self._run_command(execute)
  File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nipype/pipeline/engine/nodes.py", line 769, in _run_command
    raise NodeExecutionError(msg)
nipype.pipeline.engine.nodes.NodeExecutionError: Exception raised while executing Node distortion_merger.

Traceback:
        Traceback (most recent call last):
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nipype/interfaces/base/core.py", line 401, in run
            runtime = self._run_interface(runtime)
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/qsiprep/interfaces/dwi_merge.py", line 76, in _run_interface
            to_concat, b0_means, corrections = harmonize_b0s(
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/qsiprep/interfaces/dwi_merge.py", line 700, in harmonize_b0s
            b0_mean = index_img(dwi_nii, b0_indices).get_fdata().mean()
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nilearn/image/image.py", line 665, in index_img
            imgs = check_niimg_4d(imgs)
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nilearn/_utils/niimg_conversions.py", line 412, in check_niimg_4d
            return check_niimg(
          File "/opt/conda/envs/qsiprep/lib/python3.10/site-packages/nilearn/_utils/niimg_conversions.py", line 329, in check_niimg
            raise DimensionError(len(niimg.shape), ensure_ndim)
        nilearn._utils.exceptions.DimensionError: Input data has incompatible dimensionality: Expected dimension is 4D and you provided a 3D image. See https://nilearn.github.io/stable/manipulating_images/input_output.html.