GE Direct Fieldmaps - how to debug PE direction?

Hi all,
First off, thanks so much for all you do.

I’m feeling pretty stuck trying to preprocess my data and any advice is welcome. My two main questions are:
1. How should one debug just one step of the fmriprep pipeline without running the whole thing? I have tried to figure out the function it uses for fieldmap distortion correction but not found a clear answer.
2. In a case where fmri and B0 maps were collected with different PE directions, which should be used? Or should we just not do SDC?

Here are more details about my situation:

I have a data set of resting state fMRI and B0 maps from GE direct fieldmapping. Some of the data were collected with different PE directions for the fMRI and the fieldmaps- A/P for fMRI and R/L for fieldmaps. After this was corrected, the rest of the data is all A/P. Checking the json sidecars for the fMRI data, the “PhaseEncodingDirection” is either “i” or “j-”. However, when I run fmriprep, neither seems correct; the output of SDC looks works than before.

It’s not clear to me how dcm2niix infers PE direction for GE dicoms. In my case, the “PhaseEncodingDirection” maps to the “InPlanePhaseEncodingDirectionDICOM” variable: “i” corresponds with “ROW” and “j-” corresponds with “COL”. I’ve come to understand that i,j,k refer to the first, second, and thirds dimensions of the image, and when I check the nifti header info, these generally do correspond with what I know of the acquisition parameters. However, I’m wondering if the “PhaseEncodingDirection” in the func jsons is derived from the fmap files?

To give a concrete example, I have a func nifti that says:
qform_xorient Right-to-Left
qform_yorient Posterior-to-Anterior
qform_zorient Inferior-to-Superior
the corresponding json says:
"PhaseEncodingDirection": "i",
"InPlanePhaseEncodingDirectionDICOM": "ROW",

Assuming the “i” refers to the “xorient”, it would be R/L, which is how the fieldmap was acquired. Or at least that’s how I’m understanding it.

I used a dcm2bids conda environment to convert the dicoms to bids. The info is as follows:
“ConversionSoftware”: “dcm2niix”,
“ConversionSoftwareVersion”: “v1.0.20211006”,
“Dcm2bidsVersion”: “2.1.6”,

I’m using a fmriprep singularity container, version 20.2.3

Hi @jenncummings,

I cannot help with #1, but I can try to help with #2.

A few questions:
When you say “GE direct fieldmapping”, do you mean that you have a pair of magnitude/phase (or magnitude/field) images? In that case, the PE direction is not relevant, since the images are acquired with a Gradient Echo sequence (non EPI), and they are nominally undistorted.

For the fMRI images, there are two aspects: the direction and the sign. As you were able to infer, i,j,k refer to the first, second, and thirds dimensions of the image. To double-check that the sidecar has the correct direction, you should try to eyeball in which of the two inplane directions (i or j) the distortions occur. A relatively easy way is to go to a slice above orbito-frontal cortex and see if the brain is more or less symmetric. If it is, then the PE would be along the A/P direction (which, according to your qform, should be j or j-). If, on the contrary, the distortions are clearly asymmetric, the PE would be along the R/L direction (i or i-). (In my very limited experience, because of the peripheral nerve stimulation, most of the functional data acquired with whole body scanners have PE in the A/P direction; but there is nothing preventing you to run it in the other direction.)

Now, to determine the sign it is trickier. In general, I would recommend that you run the SDC with each of the signs, and see which one works. One sign should work and the other should make the distortions worse.

Finally, another thing to consider is making sure that your fmap images are indeed direct fieldmap images (_fieldmap.nii.gz) and that they have the correct units. It is possible that you have the direction and the sign correct, but if your units are wrong (e.g, if the values go from 0 to 2pi or from -pi to pi it is not a direct fieldmap), the SDC will not work.

I hope this helps!

-Pablo

As @pvelasco noted, the phase encoding polarity is not important for gradient echo field maps (e.g. for FSL FUGUE). However, different phase encoding polarity is required for spatial undistortion using FSL TOPUP (ideally using spin echo sequences to avoid signal dropout).

Your question is underspecified, as you do not report the GE system software or acquisition sequence used. Perhaps you could include the entire BIDS file created by dcm2niix for a better feedback.

The method used by dcm2niix is provided in the dcm_qa_ge validation dataset and in the dcm2niix notes for GE.

Hopefully, the direction and polarity are reported by the DICOM public tags In-plane Phase Encoding Direction (0018,1312) and Rectilinear Phase Encode Reordering (0018,9034):

(0018,1312) CS [ROW]
(0018,9034) CS [REVERSE_LINEAR]

However, when this fails you can inspect the private tag User Define Data GE (0043,102A) which uses a proprietary format. To understand dcm2niix’s behavior you can run it in logorrheic mode (-v 2). This generates an overwhelming amount of information, so it is typically wise to isolate a single slice from your series for conversion with this mode.

For completeness, some GE sites acquire polarity reversed images using the epi_pepolar sequence. This kludge generates DICOM images that are not completely truthful, so dcm2niix needs some logic to convert these correctly. A validation dataset with sample images and description is available.

I think this is an exhaustive answer, but perhaps @mr-jaemin has further thoughts.

@pvelasco you provided an excellent answer. However assuming the typical case that the only acquisition parameter that changed is the phase encoding polarity, your comment One sign should work and the other should make the distortions worse. This is not correct: if both your acquisitions perform phase-encoding along the same axis, but with different signs/readout-times, topup will still be able to calculate and correct for the distortions, though the estimated off-resonance field may be sign reversed. Of course, if you estimate the distortion from one pair of series and apply it to a different series, you need to make sure the direction of your polarity matches between the paired image and the series used to apply the transform.

Hi @neurolabusc and @jenncummings ,

Regarding point #2, it is hard to say without a better understanding of the pulse sequence used for B0 field mapping, as already pointed out.

I guess the “GE Direct Fieldmaps” refer to the sequence that outputs the fieldmap in Hz and a Magnitude image as in here https://github.com/mr-jaemin/ge-mri/blob/main/doc/B0_fieldmap.pdf.

We have seen a similar post not long ago:
fMRIPrep 20.2.3 directly measured B0 map SDC correction created worse functional images

Can it be solved with a similar solution (reorient or flip the sign of the fieldmap in some cases) ?

Thanks

@neurolabusc You are correct, for the pepolar case. However, it seems (still to be confirmed) that the original post refers to a fmap collected with a GRE sequence (magnitude/phase-diff). In that case, then the sign of the field map will be correct, and if you were to use fugue or topup to correct for the field map, you will obtain opposite results depending on the sign of the PE direction in your sidecar.

Thank you for all of your responses! This is very helpful.

Yes, when I say direct field mapping, I mean that I have one fieldmap in Hz and one magnitude image. My current plan is to change the sign of the PE direction in the func json and rerun fmriprep- hopefully that works.

Below is the whole json of one of the fieldmap files (with some info removed), in case it’s helpful:
“Modality”: “MR”,
“MagneticFieldStrength”: 3,
“ImagingFrequency”: 127.77,
“Manufacturer”: “GE”,
“ManufacturersModelName”: “SIGNA Premier”,
“BodyPartExamined”: “BRAIN”,
“PatientPosition”: “HFS”,
“SoftwareVersions”: “28\LX\MR Software release:RX28.0_R05_2035.a”,
“MRAcquisitionType”: “3D”,
“SeriesDescription”: “fMRI B0map”,
“ScanningSequence”: “RM”,
“SequenceVariant”: “NONE”,
“ScanOptions”: “PFP”,
“ImageType”: [
“ORIGINAL”,
“PRIMARY”,
“OTHER”,
“REAL”,
“FIELDMAPHZ”
],
“RawImage”: false,
“SeriesNumber”: 5,
“AcquisitionTime”: “10:45:3.000000”,
“AcquisitionNumber”: 1,
“Units”: “Hz”,
“SliceThickness”: 2.4,
“SpacingBetweenSlices”: 2.4,
“SAR”: 0.291624,
“EchoTime”: 0.000952,
“RepetitionTime”: 0.0071,
“FlipAngle”: 9,
“CoilString”: “RM:32NovaHead”,
“PercentPhaseFOV”: 90,
“PercentSampling”: 77.4187,
“AcquisitionMatrixPE”: 90,
“ReconMatrixPE”: 256,
“PixelBandwidth”: 1302.11,
“ImageOrientationPatientDICOM”: [
1,
0,
0,
0,
1,
0
],
“InPlanePhaseEncodingDirectionDICOM”: “ROW”,
“ConversionSoftware”: “dcm2niix”,
“ConversionSoftwareVersion”: “v1.0.20211006”,
“Dcm2bidsVersion”: “2.1.6”,
“IntendedFor”: “ses-01/func/sub-0106003_ses-01_task-rest_bold.nii.gz”

This clarifies things. The image is a gradient echo (GRE) fieldmap. Spatial distortion will use FUGUE, so the phase encoding polarity is not required. However, I note that the image is stored as FIELDMAPHZ so I suspect you need to scale this to rad/s.

I am puzzled that the PhaseEncodingPolarityGE and PhaseEncodingDirection were not populated in this image. Perhaps it is because the field map is derived from other images, or perhaps some anonymization tool stripped the relevant tags from your DICOMs.

As an aside, I am surprised that the AcquisitionMatrixPE does not match the ReconMatrixPE. Does this reflect some in-plane acceleration that dcm2niix did not detect, or is the data being interpolated. In general, I would avoid interpolation at the acquisition stage an you need to be aware that it can disrupt some tools.

Interesting! So just to be clear on next steps, you would suggest scaling the fmap images to rad/s and just using those as the input to fmriprep? I assume I should also be making sure everything is in LPS orientation as was suggested here? And the PhaseEncodingDirection entry of the func json files is irrelevant?

I’ve been told that AcquisitionMatrixPE and ReconMatrixPE are normally different for GE DICOMs - acquisition matrix size is 90x90 and GE normally provides output DICOM matrix size as 128x128 or 256x256 through interpolation, by default. We could probably disable this if necessary

I did try running fmriprep after scaling the fieldmap to rad/s (and changing the Units in the corresponding json) and this created very strange output (image below) in addition to several errors. Looking at this page of the sdcflows documentation, it seems to me that sdcflows, and thus fmriprep, would want the units to be in Hz.

It’s probably important to note that I’m using fmriprep version 21.0.2 - I’ve upgraded since my original post - but all of the problems remain. I’m currently experimenting with SyN-SDC to get around this