Wrong pixdim[4] (volume acquisition time) stored for 3D fMRIs?

[baxpr]

I have some data from a fMRI with 3D acquisition where the TR (0.05 sec) was stored in the pixdim[4] field instead of the volume acquisition time (2.6 sec). I suspect this isn't compliant with the Nifti spec but not sure.

It is a Philips file and the RepetitionTime field contains the TR, 0.05 sec. There is no field containing the volume acquisition time as far as I know - it has to be computed from these:

(5200,9230)       PerFrameFunctionalGroupsSequence
  (0020,9111)     FrameContentSequence
    (0018,9074)   FrameAcquisitionDateTime
    (0020,9157)   DimensionIndexValues (0=irrelevant,1=slice,2=volume)

I can send an example file privately if you want to look.

[neurolabusc]

This is interesting. The BIDS standard explicitly states that the RepetitionTime (in seconds) should be sourced from DICOM Tag 0018, 0080 (which is in ms). By this definition, dcm2niix is correctly reporting the data provided in the DICOM header:
(0018,0080) DS [53.3363990783691] # 16, 1 RepetitionTime
Is converted to
pixdim4 0.053336

I understand this is not the desired value here, but I am very reluctant to source data from the FrameAcquisitionDateTime. In my experience, those values are often pulled from a low precision clock, sometimes are filled with time from reconstructor not true acquisition time, and are often mangled by anonymization methods.

By the way, if we do use 0018,9074, the elapsed time based on the start time of the first slice of the first volume to the first slice of the 200th volume (which should span 199*TR) is 490.8s, suggesting a TR of 2.46633166s, whereas you suggest a nominal 2.6 s:

(0018,9074) DT [20191008110752.94]                      #  18, 1 FrameAcquisitionDateTime
(0020,9157) UL 1\1\1          
...
(0018,9074) DT [20191008 11 1603.74]                      #  18, 1 FrameAcquisitionDateTime
(0020,9157) UL 1\1\200

@drmclem do you have any insight into this?

In general, the specification of 3D EPIs seems pretty unclear, as I also discuss here. Historically, 3D EPI sequences have been typically boutique schemes. However, at least on the Siemens side they are showing up, in particular for ASL and SWI. I am happy to support these, but have little insight into how to do this robustly.

@baxpr have you tried the NIfTI conversion provided by Philips? Until the community can understand these advanced sequences, using the vendor supplied tool may be the most reliable approach.

[baxpr]

Makes sense to me. Making a policy of calculating it from frame times seems iffy to me as well.

When I said 2.6 in the first place, I was remembering wrong. My actual calculation yields the 2.466... value.

The BIDS standard also explicitly states that the BIDS RepetitionTime field and pixdim[4] should be the volume acquisition time, so I guess it's not internally consistent for these files :)

For BIDS you could maybe finesse this by using the VolumeTiming field instead of RepetitionTime, and getting the info from FrameAcquisitionDateTime. But that's still complicated to do and of unknown reliability.

I haven't looked at a Philips nifti - I'll add that to my list.

Thanks for taking a look!

[drmclem]

Hi

On Philips - The TR should always be the magnetic TR of the tissue (ie the delta T between one RF excitation of the same tissue and the next). For 2D multislice single shot EPI conveniently this is the same as the inter dynamic scan time but for anything else this is not the case. So for 3D sequences this will the "little TR" ie. time between excitations and the total time per volume will be that multiplied by the number of actual slice direction phase-encodes (which may not be the final volume slice matrix due to any slice over sampling requested and or sense).

If you use this for a dynamic fMRI seuqence it is presented on out information page but it is quite hard to find it anywhere else so this could be difficult to sort out (we don't currently implement all the suggestions for recording DICOM timing - this would ideally be in "frame acquisition duration" (0018,9220) but there are arguments against using this when the "frame" is a 2d slice taken from a "3D volume" 0

At the moment - I dont see a good solution to this (and I suspect nifti the Philips NIFTI will be no better - it exports but no obvios correct timing)

Matthew

[neurolabusc]

@drmclem thanks for your feedback. For 2D EPI, we can derive the TR from 0018, 0080. Your suggestion regarding 0018,9220 is interesting. This appear to provide the time between the first slice of the series, and the final slice of the series. In @baxpr's 200 volume image we get (0018,9220) FD 495721.00830078125ms which would suggest a TR of 2.478605s. I suspect the small error between this value and the nominal 2.46633166s is a combination of measurement error and the fact that we are measuring from the start of one volume until the end of another, which does not include the temporal gap that would follow the final volume. We really want a measure that includes the start of one volume to the start of a subsequent volume.

[neurolabusc]

@drmclem I realize that (0018,9220) works if one assumes that @baxpr's dataset has 201 volumes, rather than the 200 stored. Is it possible that the time includes one initial scan that is not stored (e.g. T1 saturation effects). If so, is there a way to determine the number of initial deleted dummy volumes either by DICOM tag or by formula (e.g. for Siemens: ROUNDUP(3001/TR))? The latest development commit warns the user, but does not set pixdim[4]:

Warning: 3D EPI with FrameAcquisitionDuration = 495.721s volumes = 200s Perhaps TR = 2.46627s (assuming 1 delete volume)

[baxpr]

@neurolabusc Same is true if we look at the AcquisitionDuration 0018,9073 field:

495.7210 (AcquisitionDuration) / 201 (hypothetical volumes) = 2.4663

Which matches my estimate based on FrameAcquisitionDateTime values.

However - my exam card says we got 200 volumes with 5 initial dummy scans. So that doesn't seem to quite fit the story.

As far as I know, the dummy scan count isn't in the DICOM anywhere - at least, I couldn't find it.

[drmclem]

I had another look at this and I think we need to be a bit more careful.

Just to clear - TR means the magnetic TR (time between sucessive RF pulses on the same slice), the scan time is the time to complete one full "scan" (so, multi slice stack, or one 3D volume) and the dynamic time is the time between repetitions of sucessive scans.

For 2D slice interleaved EPI with dynamic time shortest, these are all equal to TR but need not be (we can add dead time into the dynamic time to perform "sparse" epi.

For 3D or multishot segmented you have the added step that the TR is now shorter.

Even if you are using EPI, on the Philips you have the option of adding "dead" time into the dynamics.

Acquisition duration may not reflect the dynamic scan time (for example if navigators or frequency compensation is on that adds a bit of time to the minimum dynamic time.

Matthew

[neurolabusc]

Hmm, I remember this from PAR/REC, where the requested TR would be reported in the header, but the actual TR could be different if the user included features like prospective motion correction or perhaps inserted a delay between volumes (sparse imaging). To account for this, dcm2niix reports if there is a discrepancy between the reported TR and the TR as calculated from the start of the first dynamic to the start of the last dynamic, divided by the number of dynamics. In the case of a discrepancy, we use the calculated time, and report this to the user.

  printWarning("Reported TR=%gms, measured TR=%gms (prospect. motion corr.?)\n", d.TR, TRms);

As I recall, the prospective motion correction added about 17ms per volume on the Intera (last Philips hardware I used). Therefore, over the course of a fMRI session this could lead to a considerable clock drift. I always thought this was odd, as the TR experienced by the tissue is the applied TR, not the requested TR. From a Siemens perspective, Philips was reporting the TA as the TR.

This makes me concerned that both the 2D and 3D EPI DICOMs from Philips may show the same clock drift. With the Intera, there was a handy coaxial TTL-pulse that one could record from using a computer (at the time, most computers had TTL parallel ports) or oscilloscope.

I would be surprised if this is the explanation. It would be an amazing coincidence if the drift that @baxpr measured across 200 volumes was precisely the time of one TR.

@baxpr any chance you can measure this? You would use the same EPI pulse used for fMRI studies for fMRI synchronization. It would be interesting to have a brief series of both 2D and 3D EPI where features (prospective motion correction, temporal gap for sparse imaging) are manipulated.

[mharms]

A fundamental problem here is that the vendors don't interpret the DICOM standard consistently. Is it the responsibility of dcm2niix to resolve these discrepancies? For example, as a way to operationalize "RepetitionTime", BIDS defines it as the value of DICOM Tag 0018, 0080. But Siemens and Philips use very different definitions of "TR" for their MPRAGE sequences, for example.

[baxpr]

As usual the situation is even more confusing than I initially thought :) Yes, we should be able to measure this under a few scenarios - I will look into it.

[neurolabusc]

@mharms, your comments go to the heart of my feelings regarding @fedorov and @dclunie's presentation's at MICCAI 2017. I completely support their ideal of DICOM as the grand unified format for medical imaging. However, the thriving ecosystem of NIfTI suggests it fills a niche in our community. NIfTI is very simple and very explicit and therefore while it can not contain the breadth of meta data and image types (variable slice spacing, many transfer syntaxes, etc) it is useful for scientists and big data. Likewise, since different vendors interpret DICOM differently, an intermediate representation that unifies these parameters (TR, gradient direction, etc) is really useful. From my perspective, the utility of NIfTI is clear. I think the real debate between the views of @fedorov/@dclunie and others is whether NIfTI should just be an ephemeral translation of DICOM for consistent processing, or whether curated NIfTI is a viable archival/data sharing format (e.g. BIDS).

Ideally, I would like dcm2niix to generate consistent data regardless of the vendor it came from, and I think that is the expectation of most users. Most NIfTI tools make pretty clear assumptions of what TR means for 4D datasets: the time between the start of acquisition of a 3D volume and the start of acquisition of the subsequent volume.

On the other hand, we have to be honest that not all vendors provide all the data we would like (slice times, phase encoding polarity, true time between volumes). I think research users need to be aware that the tool developers have only a limited sample of datasets, and our reverse engineering may not work for unobserved cases or for future interpretations of the DICOM format.

You provided a tremendous service to the community by conducting a deep dive into Siemens data, providing users with a much richer set of meta data. For Siemens, a tremendous amount of meta data was hidden in their proprietary CSA header (prior to XA10). At the moment, GE, XA10 and Philips provide much less data in their DICOMs (e.g. one can not 'Phoenix' a Philips dataset, since the slice timing and phase encoding polarity are unknown). However, @baxpr has provided unique Philips datasets and @drmclem has provided insight into Philips rationale. Hopefully some of this will also impact future, richer and perhaps more consistent datasets.

[mharms]

So, if the aspirational goal is for dcm2niix to provide a consistent representation, regardless of vendor, how does it treat the "RepetitionTime" value for MPRAGE from a Philips, assuming that Philips populates DICOM Tag 0018, 0080 with the "magnetic TR" (time between sucessive RF pulses) (i.e., the inner pulse sequence loop, rather than the outer pulse sequence loop)? My broader point is that users are really asking a lot of you if they expect dcm2niix to resolve these sorts of vendor-specific interpretational differences, when they don't even agree on something as basic as what to populate in the DICOM field for "Repetition Time".

[drmclem]

Hi - to be devils advocate a bit more - modern sequences have many TR's (imagine a IR prep, multi shot non-steady state gradient echo, IR-TFE in Philips speak). And as it happens, Siemens and Philips put different ones in 0018,0080 - but this is slightly off topic as actually you are assuming TR is the dynamic time (which is what fMRI wants) and which is only = TR in a very small set of sequences - My history may be inaccurate but nifti originally only handled fMRI and 3D anatomicals but mission creep means that the design no longer has the flexibility to represent all sequences without I suspect some new design.

[neurolabusc]

The concise simplicity of NIfTI means it is inherently unable to store complexity. Therefore, most tools have few expectations for the data, for fMRI/resting state they do expect pixdim[4] is the dynamic time, which we may not be able to infer from Philips 3D EPI (the focus of this issue).

I do not know of any tools that make expectations for pixdim[4] for T1 scans, but richer meta data would probably help tools do a better job. For example, having a starting estimate for the amount of flow enhancement.

For better or worse, core NIfTI just does not have the flexibility for mission creep. The same can not be said for the BIDS sidecar. I think one could raise concerns that as these sidecars grow, the performance benefit over DICOM is lost, and as meta data from DICOM is ported to BIDS we are effectively losing the simplicity of NIfTI/BIDS. On the other hand, these proposals codify the meta data needed by scientists, and much of this data is inconsistently reported by DICOM. Implementing these proposals is going to require a lot of development work, and much of this meta-data does not exist in many existing DICOMs.

[baxpr]

@neurolabusc Here are some 3D fmri scans, DICOMs exported directly on Philips scanner, with different number of dummy scans and settings for dynamic stabilization and prospective motion correction. Our usual parameters are: 5 dummy scans, TR set to shortest = 53 ms, VAT ends up ~2.46s, dynamic stabilization ENHANCED, prospective moco OFF. Hopefully clear from the filenames/protocol descriptions what has been manipulated in each case. I don't think we're able to do sparse-type scans with this pulse sequence.
https://www.dropbox.com/s/ajbbh5nqcw2y0ov/fmris.tgz?dl=0

PS Thanks to Anna Combes for putting these together.

... there's a typo in "fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy_EnhanDynStab_Mo.dcm" - this one actually had moco OFF and is identical to "fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy.dcm"

[neurolabusc]

@baxpr thanks. It looks like the number of dummy scans is included in the FrameAcquisitionDuration (0018,9220), but the number of dummies is not recorded. This defeats my earlier attempts to determine TR from 0018,9220.

I guess the remaining option would be to use FrameAcquisitionDateTime (0018,9074) but I am really reluctant to do that as in my experience these values often rely on low resolution clocks and many anonymization tools put bogus values in these fields.

My sense is that we close this issue as the current Philips images do not allow us to resolve this robustly. We can revisit this if future Philips data provides more complete data. @xiangruili do you agree with this approach, or can you think of a clever solution that I missed?

[xiangruili]

Just catch up the topic. @neurolabusc I agree it is almost game-over if the vendor does not store the inter-volume time. We may do some guess work, but as you pointed out, it may be very wrong.

I am surprised Philips does not provide a convenient way to specify/store the inter-volume time. For some stimulus presentation, this should be pre-determined. My understanding is that users should be able to specify this then render other parameters, not the other way around.

[drmclem]

Hi

So I had a bit of a dig into this - I have tested this on a release 5.3 system with a dynamic FFE sequence in enhanced dicom.

There is private tag 2005,10a0 which contains the dynamic start time
`
Line 866: (2005,10a0) FL 0 # 4, 1 Unknown Tag & Data
Line 1005: (2005,10a0) FL 29.625 # 4, 1 Unknown Tag & Data
Line 1144: (2005,10a0) FL 59.25 # 4, 1 Unknown Tag & Data
Line 1283: (2005,10a0) FL 88.875 # 4, 1 Unknown Tag & Data
Line 1422: (2005,10a0) FL 118.5 # 4, 1 Unknown Tag & Data
Line 1561: (2005,10a0) FL 148.125 # 4, 1 Unknown Tag & Data
Line 1700: (2005,10a0) FL 177.75 # 4, 1 Unknown Tag & Data
Line 1839: (2005,10a0) FL 207.375 # 4, 1 Unknown Tag & Data
Line 1978: (2005,10a0) FL 237 # 4, 1 Unknown Tag & Data
Line 2117: (2005,10a0) FL 266.625 # 4, 1 Unknown Tag & Data
'
Which should be the information you need.

It is also present in the fmri data uploaded in the previous email.

At the moment I am not aware of an offcial DICOM tag for this information (hence the misuse of TR as a surrogate) but if there are suggestions I can send the forward.

[xiangruili]

@baxpr nice dig!
(2005,10a0) is MRImageDynamicScanBeginTime
It seems it is computed (starting with 0 for first volume), rather than timestamped. This seems a reliable tag for pixdim[4]. How do you think, @neurolabusc ?

[baxpr]

Not me, the detective was @drmclem :) Thanks, this is really helpful!

[neurolabusc]

Well spotted @drmclem - I missed that! @baxpr can you confirm that these TRs look correct, if so we have a solution everyone is happy with:

Warning: Assuming TR = 2481.44ms, not 0018,0080 = 53.6647ms
 201_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_54TR_0dummy
Warning: Assuming TR = 2481.44ms, not 0018,0080 = 53.6647ms
 301_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_54TR_1dummy
Warning: Assuming TR = 2481.33ms, not 0018,0080 = 53.6647ms
 401_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy
Warning: Assuming TR = 2481.33ms, not 0018,0080 = 53.6647ms
 501_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_54TR_10dummy
Warning: Assuming TR = 3079.44ms, not 0018,0080 = 66.6647ms
 601_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_80TR_0dummy
Warning: Assuming TR = 3079.44ms, not 0018,0080 = 66.6647ms
 701_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_80TR_1dummy
Warning: Assuming TR = 3079.44ms, not 0018,0080 = 66.6647ms
 801_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_80TR_5dummy
Warning: Assuming TR = 3079.44ms, not 0018,0080 = 66.6647ms
 901_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_80TR_10dummy
Warning: Assuming TR = 2430ms, not 0018,0080 = 54ms
 1001_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy_noDynStab_noMoc
Warning: Assuming TR = 2442.78ms, not 0018,0080 = 54ms
 1101_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy_noDynStab_Moco
Warning: Assuming TR = 2442.78ms, not 0018,0080 = 54ms
 1201_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy_RegDynStab_noMo
Warning: Assuming TR = 2481.33ms, not 0018,0080 = 53.6647ms
 1301_WIP_fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy_EnhanDynStab_Mo

[baxpr]

@neurolabusc These look right to me. They match my own computation, and are within a ms of the values computed from FrameAcquisitionDateTime (0018,9074) as well. We don't have an easy way to get an extra cross-check though - the values reported at the scanner console are truncated or rounded at the nearest 10 ms.

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_54TR_0dummy.dcm
VAT from frame timestamps: 2.481111
VAT from private field:    2.481444

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_54TR_10dummy.dcm
VAT from frame timestamps: 2.481111
VAT from private field:    2.481333

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_54TR_1dummy.dcm
VAT from frame timestamps: 2.481111
VAT from private field:    2.481444

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy.dcm
VAT from frame timestamps: 2.482222
VAT from private field:    2.481333

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy_EnhanDynStab_Mo.dcm
VAT from frame timestamps: 2.482222
VAT from private field:    2.481333

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy_RegDynStab_noMo.dcm
VAT from frame timestamps: 2.443333
VAT from private field:    2.442778

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy_noDynStab_Moco.dcm
VAT from frame timestamps: 2.443333
VAT from private field:    2.442778

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_54TR_5dummy_noDynStab_noMoc.dcm
VAT from frame timestamps: 2.430000
VAT from private field:    2.430000

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_80TR_0dummy.dcm
VAT from frame timestamps: 3.078889
VAT from private field:    3.079444

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_80TR_10dummy.dcm
VAT from frame timestamps: 3.078889
VAT from private field:    3.079444

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_80TR_1dummy.dcm
VAT from frame timestamps: 3.078889
VAT from private field:    3.079444

Getting volume acquisition time from fMRI_3DFFE_MultSht_SENSE_14slice_80TR_5dummy.dcm
VAT from frame timestamps: 3.080000
VAT from private field:    3.079444

[neurolabusc]

@baxpr please close this issue if you are happy with the performance of the latest commit to the development branch. Thanks @baxpr and @drmclem nice to have an elegant solution.

[baxpr]

Looks great and works for me. Thanks.