bids-standard / legacy-validator

Validator for the Brain Imaging Data Structure
https://bids-standard.github.io/legacy-validator/
MIT License
185 stars 111 forks source link

add a check that 'PhaseEncodingDirection' is consistent between fmap and "Intended" files #341

Closed yarikoptic closed 3 hours ago

yarikoptic commented 7 years ago

For a fieldmap to be useful, its 'PhaseEncodingDirection' must match the 'PhaseEncodingDirection' in the sidecar files of "IntentedFor" files.

Here is a quick & dirty grep on openfmri which could have false positives and miss stuff... just a quick check -- whenever you see them present and not matching, might be an issue

(git)smaug:/mnt/btrfs/datasets/datalad/crawl/openfmri[master]git
$> for ds in ds*; do echo $ds; ( builtin cd $ds; s=`/bin/ls -d sub-*|head -n 1`; [ ! -z $s ] || continue; fmap=`/usr/bin/find $s -iname '*.json' | grep fmap | head -n 1`; func=`/usr/bin/find $s -iname '*.json' | grep '\(func\|dwi\)' | head -n 1`; [ -z "$fmap" ] && continue; [ -z "$func" ] && continue; grep '"PhaseEncodingDirection' $fmap $func; ); done
ds000001
ds000002
ds000003
ds000005
ds000006
ds000007
ds000008
ds000009
ds000011
ds000017
ds000017A
/bin/ls: cannot access 'sub-*': No such file or directory
ds000017B
/bin/ls: cannot access 'sub-*': No such file or directory
ds000030
ds000031
ds000051
ds000052
ds000053
ds000101
ds000102
ds000105
ds000107
ds000108
ds000109
ds000110
ds000113
/bin/ls: cannot access 'sub-*': No such file or directory
ds000113b
ds000113c
ds000113d
ds000114
ds000115
ds000116
ds000117
/bin/ls: cannot access 'sub-*': No such file or directory
ds000119
ds000120
ds000121
ds000122
ds000133
ds000138
ds000140
ds000148
ds000149
/bin/ls: cannot access 'sub-*': No such file or directory
ds000157
ds000158
ds000164
ds000168
ds000170
ds000171
ds000172
sub-control01/func/sub-control01_task-rest_acq-voxelsize112_bold.json:   "PhaseEncodingDirection": "j-", 
ds000174
ds000177
ds000200
ds000201
sub-9001/ses-1/fmap/sub-9001_ses-1_magnitude1.json:  "PhaseEncodingDirection" : "j",
sub-9001/ses-2/dwi/sub-9001_ses-2_dwi.json:  "PhaseEncodingDirection" : "j",
ds000202
ds000203
ds000204
ds000205
ds000206
ds000208
ds000210
ds000212
ds000213
ds000214
ds000216
ds000217
ds000218
ds000219
ds000220
ds000221
sub-010001/ses-02/fmap/sub-010001_ses-02_acq-GEfmap_run-01_magnitude1.json:    "PhaseEncodingDirection": "i-", 
sub-010001/ses-02/func/sub-010001_ses-02_task-rest_acq-AP_run-01_bold.json:    "PhaseEncodingDirection": "j-", 
ds000222
ds000223
sub-01/fmap/sub-01_dir-opposing_run-01_epi.json:     "PhaseEncodingDirection": "j-"
sub-01/func/sub-01_task-mag_run-01_bold.json:"PhaseEncodingDirection": "j",
ds000224
ds000228
ds000229
ds000231
ds000232
ds000233
ds000234
ds000235
ds000236
ds000237
ds000239
ds000243
chrisgorgo commented 7 years ago

I don't believe it is necessary for the fieldmap phase encoding direction to match the scan it is intended for. In other words you can estimate a fieldmap with phase encoding direcion j and apply it to a file acquired with phase encoding direction i.

kodiweera commented 7 years ago

Local magnetic field variations across the volume may change with the directions of the magnetic fields applied (AP vs RL), and hence the information in the fieldmap.

It says at https://cfmriweb.ucsd.edu/Howto/3T/fmapproc.html

"The field map script works, but the results don't look right. What's wrong?

The most common source of error is when the user does not copy the field map graphical prescription correctly between scans so check first to make sure that has been done correctly. In particular, make sure phase-encoding direction and matrix size are the same for all your fieldmap scans and your EPI scans."

So the Phase as well as Frequency enconding directions of the fieldmap is better to match to the BOLD sequence for better results.

chrisgorgo commented 7 years ago

This sentence could be interpreted as a heuristic to avoid accidental misspecification of phase encoding direction (something we could adopt as a heuristic).

Let's get second and third opinion: @oesteban @mharms can I use a fieldmap estimated using data acquired with one phase encoding direction (say "i") to correct data with a different phase encoding direction (say "j")?

yarikoptic commented 7 years ago

another piece: https://www.fmrib.ox.ac.uk/datasets/techrep/tr04ss2/tr04ss2/node17.html "One approach to solving this problem is to use a measured field-map to ``unwarp'' the distorted images by performing pixel shifts in the phase-encode direction [25] "

So, as I read it, you can compensate if the distortion introduced in the direction of the phase-encoding since that is what fieldmap estimates... if your func/dwi phase encoding direction is different -- distortion will be in different direction, thus correction for some other direction wouldn't be applicable.

yarikoptic commented 7 years ago

@neurolabusc , how much do you know about fieldmap-based distortion correction(s)? ;-)

chrisgorgo commented 7 years ago

@yarikoptic When you apply fieldmap correction you specify which direction to unwarp the input data. That direction in theory could be different from the phase encoding direction of the data you used to estimate the fieldmap.

mharms commented 7 years ago

In theory, SE EPI volumes collected with "i" and "i-" polarity could indeed be used to correct a GE EPI scan collected with "j" (or "j-") polarity. [In fact, I think that the susceptibility field could be estimated in 'topup' if the two SE EPI volumes were "i" and "j", which is to say that it isn't even a requirement that the SE EPI volumes have 180 degree opposing polarity]. Further, it also isn't a theoretical requirement that the SE EPI have the same FOV and matrix (or voxel size) as the GE EPI that is to be corrected.

That said, it is important to be aware, when setting up acquisition protocols, that many existing software tools/wrappers may not be written to support the fully generalized possibilities that are theoretically allowed on the acquisition side -- i.e., for convenience, the software assumes some constraints on how the SE EPI was acquired. As one example with which I'm familiar, the HCP Pipelines (@glasserm @tbbrown) currently assume that the SE EPI volumes are collected with 180 deg opposing polarity, and along the same PE axis as the GE EPI that is to be corrected. They also assume that the SE EPI have the same FOV and matrix as the GE EPI.

chrisgorgo commented 7 years ago

Thank you for a prompt reply @mharms. I think this makes it clear that such inconsistency should not be a validator error.

However, such inconsistency could be a hint of a misspecified phase encoding direction. In all cases I can remember the phase encoding direction is consistent between the data used to estimate the fieldmap and the data the fieldmap was applied to. We might consider this a warning then.

oesteban commented 7 years ago

Maybe I am reading all your responses too quickly to grasp all details, but first we need to clearly state the use-cases.

All @mharms indicated seems correct to me, but please note that all of that is indicated for SE EPI fieldmaps (meaning, you use topup or a tool as such to estimate the actual distortion). So, section 8.9.3 of the BIDS spec.

Then, we have sections 8.9.1 and 8.9.2, using actual field-mapping sequences. In this case, again, the encoding does not need to match. I agree with the first statement from @chrisfilo:

I don't believe it is necessary for the fieldmap phase encoding direction to match the scan it is intended for. In other words you can estimate a fieldmap with phase encoding direcion j and apply it to a file acquired with phase encoding direction i.

It is worth checking, as @yarikoptic suggested, that these fieldmaps really need the "PhaseEncodingDirection" field. Most of the times, GRE sequences that do not show strong susceptibility-derived distortions are preferred because otherwise the distortion map needs to be inverted first before applying to the target image. So, if the fieldmap is in undistorted space, then the field "PhaseEncodingDirection" of the fieldmap can be dismissed. When the fieldmap is acquired in distorted space, then different "PhaseEncodingDirection" for fieldmap and "IntendedFor" images might be suspicious (but not forcefully wrong).

EDIT: so, in case of mismatch, a warning is probably well justified.

kodiweera commented 7 years ago

Yes, there are two different fieldmap correction methods. One (FSL - TOPUP correction) requires two scans with oposite PE directions acquired with a spin echo based sequence. In this topic , we are talking about GE- based fieldmap correction.

oesteban commented 7 years ago

AFAIK, the "IntendedFor" field applies to both types of field map estimation methods.

Anyways, I think in both situations the "PhaseEncodingDirection" of the field map does not need to forcefully match the one of the file(s) pointed to by "IntendedFor", right?

yarikoptic commented 7 years ago

I might need to read more on fieldmap corrections, but in my understanding -- the idea was to map the distortion (either via fieldmap or via TOPUP pair) which is introduced in the phase encoding direction (A1), and then to correct for that by mapping values back (along that direction only, A2) accounting for that. Knowing how distortion is introduced in some other encoding direction, wouldn't quite be useful since it wouldn't even know how distortion is introduced in the direction of interest (A1).

Please correct me if A1 and/or A2 are faulty assumption(s)? I could see that overall distortion might slightly leak into other directions (e.g. in case of AP phase encoding, voxel might need to be mapped not strictly within j axis, but also leak a bit along i and k), but as far as I see it, distortion primarily happens in the direction of the phase encoding (thus along j in AP), and would be different if we choose another axis (e.g. LR).

yarikoptic commented 7 years ago

if my logic above is more or less correct, let me make a possibly bad analogy toward "the distortion map needs to be inverted first before applying to the target image." comment. If I smash someones forehead with a hammer, and we learn how skull and matter would displace into the head, would the inverse of it be the same if I took the same action if it was done from inside the head pointing outside? my understanding is "no" since gradients of tissue density one way is different from the gradient "the other way" (e.g. skull pieces could extrude much more outside since skin/air provides less resistance), so distortion would be different. Would it be a reasonable approximation (better than none) -- may be. ok, I will stop and let experts figure out the truth... if no conclusive answer, I might scan myself in few days with varying PE directions to empirically test the issue.

chrisgorgo commented 7 years ago

@yarikoptic a reconstructed fieldmap is not a warp field, but a map of field susceptibility. The field susceptibility is not related to phase encoding direction. Image distortion is a combination of the susceptibility field and the phase encoding direction which leads to a warp field that can be used to correct for the distortion.

Thus you can use a fieldmap to correct for any phase encoding direction. Please note that tools for doing field unwarping require you to specify which direction to unwarp in.

Having said that - this is just how I understand the theory. Your experiment would be very useful. Sharing that dataset would also be valuable.

oesteban commented 7 years ago

Both are correct.

The culprit is that with whatever mapping technique, you estimate the inhomogeneity of the B0 field inside the scanner. That map is assumed to be (almost exactly) equal to the inhomogeneity that happens during the acquisition of the EPI scans.

Once you measure the map, you can use it to correct for either (i, j) distortion that happens. The PhaseEncodingDirection you use to acquire the images you'll use to estimate the field map can be different to the PhaseEncodingDirection of your EPIs.

The comment about inverting the distortion map only affects the fieldmap estimation. You can measure the variations of susceptibility across tissue using techniques affected or not by distortion. If you use EPI to measure the fieldmap (which is not very common), then your measurements are also warped.

oesteban commented 7 years ago

Regarding fieldmaps and their inversion, I find (Hutton et al., 2002) particularly clear:

If EPI is used to acquire the phase information, the map of voxel shifts will be in distorted space such that a voxel with coordinates (x, y, z) in the distorted image has been shifted by dy(x, y, z). The map of one-dimensional shifts, dy(x, y, z) can be considered as the warp field in distorted image space that describes the mapping from points in an undistorted image to points in the distorted image. To estimate the undistorted image from the distorted image, we must find the warp field that maps from points in the distorted space to points in the undistorted space. This requires the calculation of the inverse transformation dy'(x, y', z), where y' is the y-coordinate in undistorted space corresponding to y in the distorted image space.

yarikoptic commented 7 years ago

thanks you @oesteban -- this quote applies to just how to correct for the distortion, not how it could apply if the distortion/warping has happened in the opposite direction (i.e. e.g. if PE was i, while EPI was i-)

@chrisfilo thanks for the explanation! re The field susceptibility is not related to phase encoding direction -- ignorant me could only trust and check ;) will do (unless you has happened to have a data which would show that)

oesteban commented 7 years ago

Coming back to the original issue, I think that a warning could be a great idea for the validator.

glasserm commented 7 years ago

I'll note that for the topup-based form of distortion correction, it is convenient to have the GRE and SE images match exactly so that the same phase encoding direction SE as GRE can be rigidly registered in distorted space. For a regular field map, this is accomplished usually by taking the magnitude image and distorting it so that its distortions match those of the image to be corrected before registration. In theory this could be done as Mike says with field maps computed using other SE approaches, though I wonder if topup would handle 90 difference in phase encoding direction...