Summary of what happened:
Hello, I have successfully preprocessed my data using QSIPrep. Next, I would like to proceed with reconstruction. Since my data is from infants, I am concerned about the accuracy of the segmentation. Therefore, instead of using the method of generating the 5TT file with Freesurfer as done in the QSIPrep HSVS workflow, I used iBEAT to generate the segmentation maps and assembled the 5TT file myself. Modifying the input data in QSIPrep’s source code for this workflow was too challenging for me, so I decided to extract the MRtrix commands from QSIPrep’s HSVS workflow and write a small script for the reconstruction. (In this script, I only kept the key commands, omitting the quality control steps and file path adjustments.) The input files required for the script are all sourced from QSIPrep’s output directories, along with the 5TT file I generated.
I would like to ask if my script below aligns with the actual commands QSIPrep executes in the HSVS workflow(except for 5ttgen). The results of the fiber tracking I obtained appear unusual, seemingly due to a mismatch between the input DWI image and the bvec file. Could you help me understand why this is happening?
Thank you very much for your time and help!
Command used (and if a helper script was used, a link to the helper script or the command generated):
#!/bin/bash
# test.sh
resmodel=$1
maskpath=$2
dwinii=$3
dwipath=$(echo ${dwinii} | cut -d '.' -f 1)
dwimif=${dwipath}.mif
dwidtmif=${dwipath}_dt.mif
dwipddmif=${dwipath}_pdd.mif
trackpath=$4
ttnii=$5
wmmask=$6
ttpath="${ttnii%.*}"
ttmif=${ttpath}.mif
mrconvert -debug ${dwinii} ${dwipath}.mif -fslgrad ${dwipath}.bval ${dwipath}.bvec -force
mrconvert ${ttnii} ${ttmif} -force
dwi2response dhollander -mask ${maskpath} -nthreads 16 -force ${dwimif} ${trackpath}/wm_${resmodel}.txt ${trackpath}/gm_${resmodel}.txt ${trackpath}/csf_${resmodel}.txt
dwi2fod -mask ${maskpath} -lmax 8,8,8 -nthreads 16 -force msmt_csd ${dwimif} ${trackpath}/wm_${resmodel}.txt ${trackpath}/wm_${resmodel}.mif ${trackpath}/gm_${resmodel}.txt ${trackpath}/gm_${resmodel}.mif ${trackpath}/csf_${resmodel}.txt ${trackpath}/csf_${resmodel}.mif
mtnormalise ${trackpath}/wm_${resmodel}.mif ${trackpath}/wm_${resmodel}_mtnorm.mif ${trackpath}/gm_${resmodel}.mif ${trackpath}/gm_${resmodel}_mtnorm.mif ${trackpath}/csf_${resmodel}.mif ${trackpath}/csf_${resmodel}_mtnorm.mif -check_mask ${trackpath}/inliers_${resmodel}.nii.gz -mask ${maskpath} -check_norm ${trackpath}/norm_${resmodel}.nii.gz -force -nthreads 16
tckgen -algorithm ${resmodel} -maxlength 250.000000 -minlength 10.000000 -samples 4 -nthreads 16 -act ${ttmif} -force -output_seeds ${trackpath}/out_seeds_${resmodel}.nii.gz -power 0.330000 -seed_dynamic ${trackpath}/wm_${resmodel}_mtnorm.mif -select 10000000 ${trackpath}/wm_${resmodel}_mtnorm.mif ${trackpath}/tracked_${resmodel}_act_wmmask.tck
/path/to/python /path/to/convert_tck_to_trk.py ${dwinii} ${trackpath}/tracked_${resmodel}_act_wmmask.tck
tcksift2 -nthreads 16 -out_mu ${trackpath}/wm_${resmodel}_mtnorm_mu.txt ${trackpath}/tracked_${resmodel}_act_wmmask.tck ${trackpath}/wm_${resmodel}_mtnorm.mif ${trackpath}/wm_${resmodel}_mtnorm_weights.csv -act ${ttmif} -force
This is the script for iterating over each subject’s folder to call the fiber tracking script:
#!/bin/bash
# Set the QSIPrep output directory
QSIPREP_OUTPUT_DIR="/path/to/qsiprep"
SCRIPT_PATH="/path/to/mrtrix_recon.sh"
# Create the qsub scripts directory if it doesn't exist
QSUB_DIR="/path/to/qsub_scripts"
mkdir -p "${QSUB_DIR}"
# Loop over each subject directory in the QSIPrep output directory
for SUBJ_DIR in "${QSIPREP_OUTPUT_DIR}"/sub-*; do
# Extract the subject ID
SUBJ_ID=$(basename "${SUBJ_DIR}")
SUBJ_SES_DIR="${SUBJ_DIR}/ses-baseline/dwi"
# Define variables
resmodel="iFOD2"
maskpath="${SUBJ_SES_DIR}/${SUBJ_ID}_ses-baseline_acq-multishell_space-T1w_desc-brain_mask.nii.gz"
dwinii="${SUBJ_SES_DIR}/${SUBJ_ID}_ses-baseline_acq-multishell_space-T1w_desc-preproc_dwi.nii.gz"
ttnii="${SUBJ_DIR}/ses-baseline/anat/${SUBJ_ID}_5tt_useT2.nii.gz"
wmmask="${SUBJ_DIR}/ses-baseline/anat/${SUBJ_ID}_anat2dwi_affine_wmmask.nii.gz"
# Check if ttnii exists, try useT1 if useT2 is not found
if [ ! -f "${ttnii}" ]; then
ttnii="${SUBJ_DIR}/ses-baseline/anat/${SUBJ_ID}_5tt_useT1.nii.gz"
fi
# Skip if ttnii still does not exist
if [ ! -f "${ttnii}" ]; then
echo "ttnii file not found for ${SUBJ_ID}: ${ttnii}"
continue
fi
# Generate wmmask by extracting the third volume from ttnii
echo "Generating wmmask for ${SUBJ_ID}"
mrconvert "${ttnii}" -coord 3 2 "${wmmask}" -force
# Check if wmmask was successfully generated
if [ ! -f "${wmmask}" ]; then
echo "Failed to generate wmmask for ${SUBJ_ID}"
continue
fi
# Prepare the command to run
CMD="${SCRIPT_PATH} ${resmodel} ${maskpath} ${dwinii} ${SUBJ_SES_DIR} ${ttnii} ${wmmask}"
# Create a job submission script
JOB_SCRIPT="${QSUB_DIR}/mrtrix_recon_${SUBJ_ID}.sh"
cat <<EOF > "${JOB_SCRIPT}"
#!/bin/bash
#PBS -N mrtrix_recon_${SUBJ_ID}
#PBS -o ${SUBJ_DIR}/mrtrix_recon.out
#PBS -e ${SUBJ_DIR}/mrtrix_recon.err
#PBS -l nodes=1:ppn=16
#PBS -l walltime=24:00:00
#PBS -V
cd ${SUBJ_DIR}
${CMD}
EOF
# Make the job script executable and submit the job
chmod +x "${JOB_SCRIPT}"
qsub "${JOB_SCRIPT}"
done
Version:
0.22.0
Environment (Docker, Singularity / Apptainer, custom installation):
singularity
Data formatted according to a validatable standard? Please provide the output of the validator:
PASTE VALIDATOR OUTPUT HERE
Relevant log outputs (up to 20 lines):
PASTE LOG OUTPUT HERE