With an EPI sequence, a 3D volume is usually not acquired at once but rather in a sequence of 2D slices, obtained at different times (within one time of repetition (TR)). In order to fix this issue, the voxels activations of every slices of a volume can be interpolated to the same timepoint (the reference slice) if we know when each slice was acquired: this is the role of slice timing correction[1]. This is why one must ensure to know the precise TR and slice order acquisition.
Slice timing correction has been shown to reliably increase sensitivity and effect power without any adverse effect[2].
There is some debate about whether one should use slice timing correction first or motion correction (realignment) first. Indeed, unbiasing one first will ensure maximum precision for this step but will add some additional bias for the subsequent correction step, as errors accumulate. Usually, it is advised to use slice timing correction first if you use either a complex slice order sequence (i.e., anything that is not sequential such as interleaved, central, etc.) or with a sequential slice order if significant head movement is expected, or motion correction (realignment) first if you use a sequential slice order and you expect only slight head movement[3]. See also Neuroimaging_Data_Processing/Slice_Timing. There are some methods which were devised to perform joint optimization of both objectives in 4D simultaneously, such as there is no error accumulation, as implemented in nipy.SpaceTimeRealigner() of the Python library nipy[4].
For multiband acquisitions, it is necessary to use slice timings instead of slice order. Also, with multiband and a small TR, the slice timing correction can usually be skipped without much impact, but using the correction will always be beneficial. See for more information the subsection about multiband EPI acquisitions.
If you plan to use Dynamic Causal Modeling (DCM), it is mandatory to use slice timing correction.
As an alternative, if the TR < 2s, it is possible to skip slice timing correction and instead use a temporal derivative, which can account for +/- 1 second of changes in timing[5], however it was shown that slice timing correction is always beneficial, even for TR < 2s[2], and that slice timing correction increases sensitivity and power compared to temporal derivatives[2]. In addition, using both slice timing correction and temporal derivatives only decreases power compared to using slice timing correction alone[2].
SPM makes two assumptions about the convention to specify slice order: 1- the temporal order of slices would be coded by the left-to-right order of slices in a vector ; 2- The 1st slice is the most bottom slice, following the Analyze convention[9].
There are currently a few known ways to access the slice order information:
Excerpt page from an EP2D BOLD resting-state Siemens scanner sequence printout. Note the "Series" parameter, which states that the slice order is ascending sequential (and not interleaved). The most secure way is to access the MRI scanner console, which will display the exact parameters used for the EPI BOLD sequence acquisition. It is possible to save a "printout" of the sequence's parameters, which is a PDF file summarizing the scanner parameters. In general, having a copy of the sequence printout is always a good idea to check the parameters. These printout are PDF files that can be generated from the scanner software interface.
From the DICOM files for some recent scanners (see Siemens below).
If you created a BIDS sidecar when converting your DICOM images to NIfTI (e.g. you used dcm2niix or dicm2nii to convert your images), the slice order can be inferred by the "SliceTiming" tag in the JSON-format BIDS file. The "SliceTiming" tag lists the time each slice was acquired in seconds[7].
Some NIFTI viewers (here Mricron) can display the slice order if the information is available. Even without a BIDS sidecar, some DICOM to NIFTI converters (e.g. dcm2niix or MRIconvert) can keep the slice order (and TR) informations inside the generated nifti files, and some nifti viewers (e.g. mricron) can display it or you can then access this information programmatically like this[8]:
fMRIname = 'path_to_nifti_file.nii';
slice_order = 0; % Set 0 to autodetect
%these are the possible slice_orders http://nifti.nimh.nih.gov/pub/dist/src/niftilib/nifti1.h
kNIFTI_SLICE_UNKNOWN = 0; %AUTO DETECT
kNIFTI_SLICE_SEQ_INC = 1; %1,2,3,4
kNIFTI_SLICE_SEQ_DEC = 2; %4,3,2,1
kNIFTI_SLICE_ALT_INC = 3; %1,3,2,4 Siemens: interleaved with odd number of slices, interleaved for other vendors
kNIFTI_SLICE_ALT_DEC = 4; %4,2,3,1 descending interleaved
kNIFTI_SLICE_ALT_INC2 = 5; %2,4,1,3 Siemens interleaved with even number of slices
kNIFTI_SLICE_ALT_DEC2 = 6; %3,1,4,2 Siemens interleaved descending with even number of slices
[pth,nam,ext,vol] = spm_fileparts( deblank(fMRIname(1,:)));
fMRIname1 = fullfile(pth,[ nam, ext]); %'img.nii,1' -> 'img.nii'
if slice_order == 0 %attempt to autodetect slice order
fid = fopen(fMRIname1);
fseek(fid,122,'bof');
slice_order = fread(fid,1,'uint8')
fclose(fid);
if (slice_order > kNIFTI_SLICE_UNKNOWN) && (slice_order <= kNIFTI_SLICE_ALT_DEC2)
fprintf('Auto-detected slice order as %d\n',slice_order);
else
fprintf('%s error: unable to auto-detect slice order. Please manually specify slice order or use recent versions of dcm2nii.\n');
return;
end;
end
Here are a few other notes one should pay attention to:
The scanner head position for the first slice relatively to the patient can be found in the DICOM fields Image Position Patient (0020, 0032) together with Image Orientation Patient (0020, 0037) and Patient Position (0018, 5100). The latter is a required DICOM field and is necessary for softwares to interpret the orientation of the images, so you can reliably use this field. See DICOM reference C.7.3.1.1.2 for more informations on the meaning of the different values[10].
SPM expects the first (spatial) slice to be the bottom slice (orientation acquisition Transversal and from Inferior to Superior). In other words, the slice 1 as in your slice order must represent the bottom-most slice of the brain. If not, you need to change the reference slice and use as slice order TR - INTRASCAN_TIME - SLICE_TIMING_VECTOR[11].
If possible, prefer to use the slice timing instead of the slice order, as some scanners will round off some of the slice timing in order to more simplify and more reliably acquire these slices over all volumes (because it is nearly impossible for a machine to acquire at a time that is defined with a very long floating point number, if it is too precise the mechanical parts might not follow, hence the rounding, which logic differs depending on the brand and machine). Using the slice order will always assume a perfectly spaced acquisition, whereas using the slice timing, as found in the DICOMs, BIDS or Scanner console log, will account for these specific roundings.
When you do slice timing correction, all slices of one volume are interpolated in time to one slice of reference. This reference slice becomes the most accurate slice since it gets no interpolation, only other slices are interpolated.
In SPM, you can change this parameter in Reference Slice (refslice). The slice of reference is 1 by default in SPM, which corresponds to a slice order sequential ascending or interleaved ascending. If your slice order is different, then you need to change this value. In other words, if the first slice in the slice order is not 1, you need to change the reference slice.
If you want to choose the first acquired slice as a reference slice, then you need to use the slice spatial number as reference slice[11]. For example, with a slice order sequential descending 4 3 2 1, to set the first temporal slice as the reference slice, you need to set the reference slice to 4. If instead of slice order you use slice timing in seconds (e.g., 2.0 1.0 0.0 1.5 0.5), then you need to also specify the reference slice in seconds (e.g., to set the first slice as the reference, use '0').
However, you can also choose to use another reference slice as you wish. Although setting the reference slice to the first is the easiest and the most common, another common reference is to use the middle slice for more precision with sequential slice orders. Indeed, using the middle slice theoretically guarantees that we minimize the amount of temporal interpolation error, because then the maximum interpolation will be of TR/2 (in negative and positive shifts). Also, since we assume we are in a sequential slice order, the middle slice in time is also the middle slice in space, thus we also minimize the spatial interpolation error and push the accumulating errors to the top and bottom slices, where there is usually less tissues of interest.
There is however some debate about the benefits of using any reference slice other than the first[12]:
« As Rik Henson noted, it doesn't matter what the slice actually contains as regards its use as a reference (the slice could cover only air!). I would disagree, however, that it is optimal to use the middle slice of the acquisition series or a slice about which one might have an a priori hypothesis. Instead, I think the most appropriate reference slice is the first slice acquired in a given TR.
I would argue that this is because the goal of slice acquisition correction is to create a data set that corresponds to what one would have obtained if all of the slices were acquired simultaneously at the moment in time that the TR was initiated. If one models in the G matrix the pattern of neural activity (and induced hemodynamic response) with regard to the onsets of the TRs, any other reference slice would result in a fixed offset (e.g., 1/2 the TR) between the model and the actual data. Of course, one could also shift the covariates of the G matrix by this fixed offset, but this seems needlessly complicated!
Also, I don't believe that it is the case that using the middle slice affords any greater "accuracy" in making the correction. The slice correction routine makes the assumption that no meaningful power is present in the data above the Nyquist frequency. (For which there is some empirical support). Given this assumption, all shifts in time using sinc-interpolation are equally valid. »
It is important to note that whatever reference slice you use, you should check that the Microtime Onset (fMRI_T0) in the statistical analysis module corresponds to your Reference Slice (refslice) of the slice timing module to ensure that your onsets during the statistical test are shifted appropriately (else they might be too early or too late!)[3] ·[13] ·[14] ·[15]. Note that by default, the microtime onset is set to 8 over 16 in SPM12, which corresponds to a middle reference slice, thus if you use the first temporal slice as the reference slice, you should change the microtime onset to 1.
In any case, one must ensure that the slice timing correction reference slice is the same as the statistical test microtime onset. To make it easier, you can change the microtime resolution to the number of slices there are in each EPI volume, so that the number of bins in the statistical test will be the same as the number of slices used for slice timing correction, hence the microtime onset will correspond to the slice position in time[16]
Note that the microtime onset is to be set in the temporal convention (the number is the slice position in time) and scaled to the microtime resolution (you can also set microtime resolution to the number of slices of your EPI volumes, which will make things easier), whereas the Reference Slice is in the spatial convention (the number is the slice position in space).
The following applies to presumably all Siemens scanners using the Syngo system, but might or might not apply to other scanners. Please note that most of the points below are only correct if the "Orientation" scanner parameter is "Transversal" [17]. Also please note that we use here the 1-based indexing convention (i.e., slices start at 1 then 2, 3, ...).
Slice times (in milliseconds) can sometimes be read directly from the header of DICOM files (from any volume of an EPI BOLD sequence after the first one) using the private vendor field MosaicRefAcqTimes (0019, 1029) below the field (0019, 0010) SIEMENS MR HEADER:
This will give you the relative time of acquisition of each slice, which you can combine with (0008, 0033) Acquisition Time to compute the absolute time of each slice. Note that (0019, 1029) is a private field, so it might not always be present on Siemens machines (also your sequence definition might change this).
If you want the slice order offsets instead of timing:
[~, slice_order] = sort(slice_times);
Note that with the slice timing, you don't need the slice order, as you can provide the slice timing in milliseconds directly to SPM instead of the slice order.
The slice order is defined by the parameter "Series", and not "Multi-slice mode"[17] (although both present the option "interleaved", only the option in "Series" is related to the slice order, the "Multi-slice mode" being interleaved for an EPI sequence just means that the acquisition is simultaneous and is equivalent to single-shot[17]). There are usually three options: ascending, interleaved (ascending) and descending. Interleaved is particular on Siemens machines as this mode always acquire in ascending fashion, but the starting slice will change depending on the number of slices[17]: if the total number of slices is even, the slice order will be even-first, otherwise with an odd number of slices the slice order will be odd-first. Note that this odd-first vs even-first interleave behavior that is specific to Siemens machines is only true if using a standard sequence definition, as a user-defined sequence (such as CMRR's[18]) can change this to always be odd-first as with other MRI machines[19]. See the table below for a summary of the possible slice orders for Siemens machine and the sample codes to supply to SPM. If you are wondering which of interleaved or sequential/contiguous acquisition is better, a brief answer is that interleaved is better to reduce slice cross-talk but slower at recovering from T1 signal variations due to motion noise than sequential acquisition[20].
Siemens advises the scanner parameter Transversal to always be set to F >> H, which means "from foot to head", in order to simplify viewing images. The Transversal parameter should be described in the protocol printout that the machine can output. Else if you have H >> F, the mosaic display (ie, how the slices numbering will be stored) will be reversed[17]. Although the official documentation specifies that the slice order is reversed[17], other studies found in practice that only the mosaic display is affected, and not the slice order acquisition nor storage, thus this parameter does not affect the slice order[21] ·[22] ·[23] ·[19], but you should be aware that the mosaic display might be confusing and not reflect the slice order. Note also that since H >> F always reverses the mosaic display, an interleave will then have a mosaic display (but not slice order) that will always be odd-first, whether the number of slices is odd or even[17] ·[23]. However, on some newer machines, such as the Magnetom, the slice timing is indeed affected, thus not only the mosaic display is affected but also it seems the slice timing as can be found in the related DICOM field (see above), as is stated in the Siemens documentation. In any case, if H >> F, one should check the slice timing from the DICOM files to ensure if the slice timing got reversed or not. See the table below for a summary.
Ascending vs descending meaning: to know what exactly ascending and descending refer to, please refer to your printout: it will detail for each dimension what is the ascending direction. By default on Siemens machines, Transversal should be F >> H which is foot to head, Sagittal should be R >> L which is right to left, Coronal should be A >> P which is anterior to posterior[19]. Finally, the main slice acquisition dimension, and the one that must be accounted for slice timing correction, is given by the Orientation parameter, which can be either Transversal, Sagittal or Coronal[17]. Thus, if your main orientation is Sagittal with direction R >> L, it means that the slices will be acquired right-first and left-last. However, according to practical observations, the direction set in the scanner parameters appear to not matter much except for mosaic display, as the slices are still acquired in the default direction[21] ·[23]. SPM expects the main acquisition dimension to be Transversal, so slices are acquired from foot to head[11].
Here is a table summarizing the various slice orders that are possible with Siemens machines, with sample code to supply to SPM for correct slice timing correction[21] ·[23]:
Note that the SPM code is affected by the acquisition mode, the number of slices and potentially the transversal direction (depending on the machine model, this can impact only the mosaic display order or also the slice timing).
Depending on the TR, Siemens machines will acquire additional "dummy scans", which are volumes that are discarded prior to beginning the real sequence acquisition in order to reduce magnetic saturation. The number of dummy scans is chosen to guarantee at least 3 seconds of stabilization: Dummy scans = ROUNDUP(3001/TR), and is enabled depending on the TR: if there is only one volume in the sequence, there is no dummy scan; if 1501 < TR <= 1001ms: 3 dummy scans; 3001 < TR <= 1501ms: 2 dummy scans; else if TR > = 3001ms: 1 dummy scan[24].
The Siemens machines can be subdivided into two groups depending on the underlying system: Numaris 3/3.5 and Syngo. Whereas the slice order also changes the visualization of slices for Numaris-based machines (hence you can easily deduce the slice order by just looking at the slices with any DICOM viewer), the Syngo-based machines differenciate the slices storage with the slice acquisition order, what is termed "mosaic". Indeed, in a mosaic, the slice storage order is dependent on the slice spatial number, instead of the slice temporal acquisition number. Numaris machines are using SUN OS as the operating system, whereas Syngo machines use Windows. Numaris machines include Open, Impact, Vision and older Harmony, Symphony scanners ; Syngo machines include Concerto, Harmony, Symphony, Trio, Allegra, Magnetom[25].
Multi-band acquisition, also known as POMP (GE), Simultaneous excitation or SMS for simultaneous multi-slice (Siemens), Multi-slice (Philips), Dual-slice (Hitachi), and QuadScan (Toshiba), is a technique that allows the acquisition of a few slices simultaneously[27]. In practice, if you look at the slice timing, you will see that multiple slices have the same slice timing, which you will never see with the slice timing of a non-simultaneous EPI acquisition.
It is very difficult to check if you have multi-band acquisition enabled on your sequence, as the printout may not reveal this information. One reliable way is to check the slice timing in the scanner console, or also in the DICOM if the private vendor field is available (see above)[28] ·[29]: if two values are the same in the slice timing, then your acquisition is multi-band. For newer machines such as Siemens Magnetom VIDA, the DICOM can include the private fields (0018,9077) Parallel Acquisition and (0018,9078) Parallel Acquisition Technique to describe the use of multiband. Note that these two tags will only report if multiband is enabled ("(0018,9077) CS [YES]" and "(0018,9078) CS [SMS]") not the acceleration factor. Alternatively, the field (0021,1009) will report the in-plane (e.g. iPAT, SENSE, GRAPPA) and between slice (e.g. SMS) acceleration (for example "(0021,1009) LO [p2 s4]" would suggest multiband 4).
If you use multi-band acquisition, you cannot use the slice order as an input to slice timing correction, since a slice order cannot represent multiple slices acquired at the same time (if it was a matrix it would be possible, but SPM only accepts a vector). However, you can use the slice timing instead of slice order when using a multi-band EPI acquisition[30] ·[31]. If you do know your slice order but not your slice timing, you can artificially create a slice timing manually, by generating artificial values from the slice order with equal temporal spacing, and then scale the numbers on the TR, so that the last temporal slices timings = TR - TR/(nslices/multiband_channels).
Is slice timing correction necessary? (for small TR, multiband, etc)
With very small TRs such as multiband, it becomes very complicated to slice time correct because there are usually a high number of slices and any slight movement might make slice time correction to correct the wrong slice, and this correction is less necessary as the TR is small and the temporal difference between slices is thus reduced. In such cases, you can choose to skip slice timing correction (or use a temporal derivative, which can account for +/- 1 second of changes in timing[5]).
However, it was shown that slice timing correction is always beneficial, even for fast TR < 2s such as obtained with multiband[2]. Indeed, it was observed that slice timing correction helps more with the slices that are more delayed (ie, towards the end of the slice order) and that slice timing correction increases the statistical effect size the longer the TR, with diminishing returns as the TR is shorter, down to 0.5s where the common methods of slice timing correction do not help anymore[32]. However, with new slice timing correction methods such as upsampling and lowpass filtering, there is still a significant gain, that is constant whatever the TR[32].
The time of repetition can be retrieved from DICOM files in the mandatory field (0018,0080) Repetition Time, which will give you the value of the TR setting set in the scanner. It can also be found in NIFTI files headers except if an anomyzation process stripped it.
However, you might find that the TR value might be a floating value different from the TR set in the scanner settings. It seems some scanners such as the Philips Healthcare Ingenia include in the Repetition Time field the "real" time of repetition (the time it really took to acquire one volume of this sequence). You will hence get a floating point value, which might be a bit off from the TR you set in the settings (eg, if TR is 2 then you can have a "real" TR of 2.00392).
It would be interesting to be able to use this precise TR timing per sequence to more precisely correct the TR at the subject level, unfortunately most (all?) softwares do not currently support this feature.
At the statistical analysis step, after the "fMRI model specification" module, thus after preprocessing and slice timing correction, SPM does not work anymore in the original series domain but in its own: the (hopefully) slice time corrected volumes are split in "microtime bins". This microtime design is made to super sample the timeseries, in other word to increase the resolution of the BOLD signal, so that more precise data points can be sampled at a much greater time resolution than the original timeseries.
This design is configured by two parameters in the "fMRI model specification" module: microtime resolution (fmri_T) and microtime onset (fmri_T0). The first sets the total number of bins per volume (in other words the "magnification" factor, e.g., if you set 16 you will get 16 points per TR and voxel instead of one), whereas the second is used to shift the microtime design to correspond to the slice timing correction parameters (i.e., if you set the reference slice in the slice timing correction module, you need to change the microtime onset in the statistical test).
The HRF Time Derivative can be enabled during the fmri model specification to allow for tolerance of small shifts in time (< 1s), and can thus be used as an alternative to slice timing correction (even though slice timing correction increase sensitivisy more than HRF Time Derivative[2], the latter is simpler and more generic to use).
See HRF 'Informed' Basis Set for more information.
↑ abcdefSladky, R., Friston, K. J., Tröstl, J., Cunnington, R., Moser, E., & Windischberger, C. (2011), "Slice-timing effects and their correction in functional MRI", Neuroimage, 58(2): 588–594{{citation}}: CS1 maint: multiple names: authors list (link)
↑ abRik Henson, Slice timing questions, archived from the original on February 08, 2013, retrieved 2017-10-23{{citation}}: Check date values in: |archivedate= (help); Unknown parameter |dead-url= ignored (|url-status= suggested) (help)
↑Roche, A. (2011), "A four-dimensional registration algorithm with application to joint correction of motion and slice timing in fMRI", IEEE transactions on medical imaging, 30(8): 1546–1554
↑ abSPM12 Manual(PDF), archived from the original(PDF) on 2017-07-11, retrieved 2017-11-08{{citation}}: Unknown parameter |dead-url= ignored (|url-status= suggested) (help)
