|Home | About | Journals | Submit | Contact Us | Français|
A variable-density spiral k-space trajectory is introduced for brain functional magnetic resonance imaging (fMRI). The proposed spiral trajectory consists of an Archimedean spiral from the k-space origin to an arbitrary fraction r of the maximum k-space radius, extending beyond this point with a variable-density spiral in which the sampling density decreases as the k-space radius increases. It therefore permits a reduction in readout time at the expense of undersampling only the high spatial frequencies, in which the energy in T2*-weighted brain images is low. The trajectory was implemented in a 2D spiral-in/out sequence, and single-shot high-resolution (1.72 ×1.72 mm2 in-plane) fMRI data were acquired from human volunteers. Compared to a two-shot fully Archimedean spiral sequence with the same spatial coverage and total scan time, the variable-density sequence yielded greater activation magnitudes with improved temporal efficiency and minor artifacts.
Functional magnetic resonance imaging (fMRI) using blood-oxygen level dependent (BOLD) contrast (1-3) has become a powerful and widely-used method for advancing fundamental understanding of the brain. BOLD contrast is based on local changes in the oxygenation of cortical vasculature consequent to neural activation, and can be mapped across the brain using T2- or T2*-weighted rapid imaging sequences.
While fMRI data are conventionally acquired at in-plane voxel dimensions of ≥ 3 × 3 mm2, there is a growing interest in imaging at higher spatial resolutions. A key motivation stems from the investigation of neuroscience questions that depend critically on the detailed spatial mapping of functional architecture (e.g. (4,5)). In addition, although imaging at higher resolution decreases the signal-to-noise ratio (SNR) in magnetically uniform brain regions, reduced voxel sizes can provide increased SNR in regions near air-tissue interfaces that are compromised by susceptibility-related dephasing (6), and furthermore may improve the contrast-to-noise ratio by reducing partial-volume effects. However, for single-shot acquisitions, the long readout duration necessary for high resolution is more sensitive to motion and off-resonance effects; thus, multi-shot (segmented or interleaved) techniques (7), partial k-space methods (8), or parallel imaging (9,10) are required to sufficiently shorten the acquisition window. Unfortunately, multi-shot acquisitions can exhibit artifacts from motion- or pulsatility-induced inter-shot inconsistencies, and furthermore suffer a loss of temporal resolution; acceleration via parallel imaging reduces SNR in the thermal noise-dominated limit encountered at high resolution. Thus, there is a need for high-resolution techniques that maintain high SNR.
In the present work, we propose to increase the efficiency of single-shot fMRI by using a variable-density spiral trajectory. Variable-density spiral trajectories have been previously proposed (11,12) and applied to MR fluoroscopy (11), spectroscopy (13), and for reducing motion artifacts in cine and diffusion-weighted imaging (14,15). Lee et al. applied variable-density sampling to a 3D stack-of-spirals trajectory and demonstrated its utility in limb perfusion imaging (16). The use of a composite constant/variable-density (logarithmic) spiral trajectory was proposed by Cline et al. (17) and was shown to reduce artifacts in cardiac imaging. Here, the advantages of variable-density sampling are demonstrated for brain fMRI. Since the spectral energy of brain images is concentrated at the center of k-space and decays markedly toward the outer edges, we propose to use an Archimedean spiral that critically samples k-space between the origin and a user-specified k-space radius, and extend the trajectory to the desired resolution with a judiciously undersampled variable-density spiral, similar to designs suggested in (12,17). Because the high spatial frequencies contain little energy, the aliasing artifacts contributed by these components will be small and perhaps outweighted by the benefits conveyed by a more rapid readout (12,16). Unlike some implementations of the variable-density spiral (15), we choose not to use additional oversampling near the k-space origin; the Archimedean spiral already possesses a degree of oversampling at the origin because the trajectory starts with zero velocity due to gradient slew rate constraints (18,19).
This paper first describes a procedure for designing the gradient waveforms for the proposed variable-density spiral trajectory. We implement the variable-density spiral in a two-dimensional (2D) single-shot spiral-in/out (20) pulse sequence, and demonstrate its performance in detecting BOLD functional activation in visual, auditory, and sensorimotor regions at high (1.72×1.72mm2 in-plane) resolution. A comparison with a conventional 2-shot interleaved, Archimedean spiral-in/out sequence is presented, as the latter constitutes a fully-sampled trajectory without the use of parallel imaging and utilizes two shots to obtain acceptable image quality at the desired resolution; multi-shot spiral sequences are employed successfully in current high-resolution fMRI studies (21-23). This work was presented in preliminary form at the 2010 ISMRM annual meeting (24).
The proposed variable-density spiral trajectory, parameterized by 0<r<1 and α>1, is given by
where k(t) = kx (t) + iky (t) is the complex k-space location in radians, τ(t) is a function of time, λ1 = 1/D, λ2 = kmax = πN/D, and ω = πN, where N and D are the desired matrix size and field of view (FOV), respectively. The trajectory consists of a critically sampled Archimedean section, kA, and a variable-density section, kV, that join at a fraction r of the maximum k-space radius kmax. The constants r0 and τ0 are chosen so that the k-space radius and its first derivative (with respect to τ) are continuous at the boundary between the two sections (i.e. at |k| = rkmax), yielding
Note that choosing α close to 1 results in a trajectory that approaches uniform spacing between the turns of the spiral, and larger values of α correspond to sparser k-space sampling in the range |k|> rkmax.
The function τ(t) can be described in terms of its values during the Archimedean ( τA (t)) and variable-density ( τV (t)) sections of the trajectory, noting from Eq.  that τ = r at the boundary between the two sections and that τmax τ0 + (1− r0)1/α is the value of τ at the end of the trajectory (when |k|= kmax; Eq. ). Similarly, the gradient waveform is defined as
The functions τA (t) and τV (t) can be derived such that the gradient waveforms satisfy slew-rate and amplitude-limited hardware constraints, thereby minimizing the readout duration. The Archimedean section of the trajectory begins in the slew rate-limited regime, and transitions to the amplitude-limited regime if the gradient amplitude attains the maximum allowed value before reaching the variable-density section of the trajectory. If the latter occurs, the variable-density section begins in the amplitude-limited regime; otherwise, it begins in the slew rate-limited regime, and transitions to the amplitude-limited regime if the gradient amplitude attains its maximum value prior to completion of the trajectory.
For both the slew rate- and amplitude-limited regimes of the Archimedean section of the trajectory, the function τA (t) and associated gradient waveform gA (t) are precisely as given in Eqs. 7, 12, and 15 of (18), using the relationship ωτ = θ since the equations in (18) are specified in terms of θ. For the variable-density section of the trajectory, the slew rate- and amplitude-limited forms of τV (t) and gV (t) are derived in a manner analogous to (19), accounting for the offsets r0 and τ0. Figure 1 illustrates a trajectory design that includes both slew rate- and amplitude-limited cases for a single-shot readout (N=128, D=22 cm).
The variable-density spiral parameters (r, α) dictate the readout duration, the k-space undersampling pattern, and – in the case of spiral-in or in/out trajectories – constrain the minimum TE. Thus, different choices of (r, α) yield images with different tissue contrast, BOLD sensitivity, blurring, and aliasing properties (Fig. 2a).
The aliasing patterns associated with different choices of (r, α) can be characterized by the point spread function (PSF), calculated from the inverse Fourier transform of the density-weighted sampling function of the trajectory. Figure 2b shows the two-dimensional PSFs for the indicated set of parameters, and Fig. 3 depicts their corresponding radial energy. The location of the aliased energy reflects the decrease in effective FOV associated with variable-density sampling. Note that the parameter values chosen in columns 2-4 of Fig. 2 illustrate the effects of varying (r, α) while maintaining a constant readout time. For a fixed readout duration, higher values of r (which necessitate correspondingly larger values of α in order to complete the trajectory in the allotted time) lead to greater apodization of the high spatial frequencies, resulting in increased blurring but smaller aliased energy. Lower values of r (permitting smaller values of α) lead to larger aliased energy, but situated farther from the origin of the PSF, thereby expanding the effective FOV. Decreasing the total readout duration leads to an increase of total aliased energy.
Variable-density spirals suffer from a slight loss of SNR efficiency due to the non-uniformity of the trajectory. With the same spatial resolution and scan time, the SNR loss can be expressed by the multiplicative factor
as defined in (12), where AK is the area of the sampled k-space, D(k) = D(kx,ky) is the local sampling density, T is the total readout duration, and Δt is the gradient sampling interval. This factor tends to vary only slightly from unity, ranging from 0.960 to 0.972 for the variable-density parameters shown in Fig. 2 (with η=0.976 for an Archimedean spiral; it is < 1 due to slightly non-uniform sampling around the origin with the current design (18)). However, the SNR of the single-shot variable-density sequence will additionally deviate from that of single- or multi-shot Archimedean spiral sequences due to differences in the achievable TE and TR. Table 1 compares the calculated SNR for a 2-shot Archimedean trajectory, a variable-density trajectory (with r=0.5 and α=3.6, as used in the fMRI experiments below), and a single-shot Archimedean trajectory, assuming identical spatial resolution and scan time and using their respective minimum TE/TR values in a spiral-in/out sequence with N=128, FOV = 22 cm (T1=1350 ms, T2* = 50 ms). The relative SNR was also measured experimentally in a subset of 4 subjects (see section Data Analysis: SNR).
The proposed variable-density spiral trajectory was designed according to the method described above, and was implemented in a 2D gradient echo spiral-in/out pulse sequence using the scanner’s standard pulse sequence programming language, with user control variables for easily modifying r and α. A single-shot variable-density spiral-in/out sequence was compared with a 2-shot interleaved Archimedean spiral-in/out sequence (with respect to measures of BOLD functional activation and signal-to-noise ratio (SNR) at high (128×128) resolution, using experiments described below.
BOLD functional data from human volunteers were obtaining using a sensory task. One version of the task (“SM-block”) was a block design in which each 30-s block contained a 15 s “on” period and a 15 s “off” period, repeated 8 times for a total duration of 4 min. “On” periods consisted of simultaneous visual, auditory, and sensorimotor stimulation: the visual stimulus was a circular checkerboard of alternating black and white contrast that reversed at 4 Hz; auditory stimuli consisted of a randomized binaural tone sequence (4 Hz, synchronized with the visual stimuli); sensorimotor stimulation was delivered by pneumatically-driven plungers arranged in left- and right-hand gloves so as to push the fingers up and down in a randomized pattern at 4 Hz, in synchrony with the checkerboard reversals and auditory tones. During the “off” periods, subjects viewed a fixation cross at the center of the screen. A second version of the task (“SM-event”) was an event-related design, in which each “on” event lasted for 1 s and whose onset times were randomly jittered with inter-trial-intervals ranging from 20 s to 36 s. The total run length of this task was 5 min. Stimuli for the event-related design were identical to those of the block design.
For each task (SM-block, SM-event), two separate runs were acquired in order to compare the following two pulse sequences (Fig. 4): (1) a two-shot Archimedean spiral-in/out sequence (henceforth abbreviated AR2), and (2) a single-shot variable-density spiral-in/out sequence (r=0.5, α=3.6, abbreviated VD). For both sequences, the TE was set to the minimum value necessary to accommodate their respective spiral-in readout durations as well as the RF pulse and slice-select gradients (TE = 37 ms for each interleave of the AR2 sequence, and TE = 48 ms for the VD sequence), and identical spatial coverage (slice prescription) was used.
Two comparisons were performed (Table 2). In the first comparison, the time per volume acquisition (frame), and total number of time frames were set to be identical across the two pulse sequences, providing conditions of equal temporal resolution, inflow characteristics, and statistical degrees of freedom. For this comparison, the SM-block task was used, and 18 oblique axial slices (4mm slice thickness, 1mm gap) were acquired with the following parameters: for AR2, TE/TR=37ms/1.5s, flip angle=71°, #frames=88; for VD, TE/TR=48ms/3s, flip angle=84°, #frames=88. Note that the minimum time per volume, 3 s for 18 slices, was determined by the (longer) AR2 sequence.
In a second comparison, we take advantage of the fact that the time per frame can be shorter for the VD compared to the AR2, allowing for improved temporal resolution and statistical power with the same spatial coverage. We compared the 2 sequences using the SM-event task, and 24 axial slices (4mm slice thickness, 0mm gap) were acquired with the following parameters: for AR2, TE/TR=37ms/2s, flip angle=77°, #frames=74; for VD, TE/TR=48ms/2.5s, flip angle=81°, #frames=120. The event-related task was utilized for this comparison since the detection power in this paradigm is more sensitive to temporal resolution than are block designs. Parameters of the two comparisons are summarized in Table 2.
A subset of subjects were also scanned with the two sequences while resting with their eyes closed for 5 min. These data were acquired for purposes of temporal SNR comparison, described below. Parameters of this acquisition were identical to those used in the SM-event task.
Twelve healthy adult volunteers (mean age = 30.3 years) were scanned after providing informed consent in accordance with a protocol approved by the Stanford Institutional Review Board. Nine volunteers performed the block design version of the task, and 5 performed the event-related version (2 performed both). A subset of 4 subjects were scanned at resting state for the SNR comparison. For each task, the order in which the two sequences were run was counterbalanced across subjects. All imaging data were acquired at 3T (DVMR 750, rev 20, GE Healthcare Systems, Milwaukee, WI) using a GE 8-channel head coil. Head movement was minimized with a bite bar. T2-weighted fast spin echo structural images (TR = 3000 ms, TE = 68 ms, ETL = 12, FOV = 22 cm, matrix 192 × 256) were acquired for anatomical reference. A higher-order shimming procedure (25) was employed to reduce B0 heterogeneity prior to the functional scans. All functional scans used the following parameters: FOV = 22 cm, matrix size 128 × 128, bandwidth = ±125 kHz. Cardiac and respiratory processes were monitored using the scanner’s built-in photoplethysmograph placed on the right toe and a pneumatic belt strapped around the upper abdomen, respectively.
Functional images were reconstructed off-line using density compensation, gridding (with the calculated trajectories) and fast Fourier transforms. A Kaiser-Bessel kernel with an oversampling ratio of 2 was used for gridding (26). Density compensation was performed by (1) gridding the sampling function to obtain the weighting function, as described in Jackson et al. (26), and (2) using the approximation of Meyer et al. (7) as a second-order correction to (1). Spiral-in and spiral-out images were separately reconstructed, and combined using a weighted average based on their time-series mean values (20). Linear shim corrections for each slice were applied using individual field maps obtained during the scan (27), and a frequency navigator correction was employed to reduce blurring due to respiration-related dynamic off-resonance (“DORK”; (28)). Corrections were also performed for concomitant field effects (29).
Functional images were pre-processed using custom C and Matlab routines. Pre-processing included slice-timing correction using sinc interpolation, removal of linear and quadratic temporal trends, and corrections for physiological noise by applying RETROICOR (30) for the removal of time-locked cardiac and respiratory artifacts, and filtering out low-frequency respiratory volume and heart rate effects using methods (RVHRCOR) described in (31). No spatial smoothing was performed.
For the SM-block and SM-event tasks, the time series from each voxel was analyzed for functional activation using linear regression against the modeled task waveform, which was formed by convolving the binary stimulus vector with a canonical hemodynamic response function (32). A sigma filter was used to cluster pixels in a 3×3 region, as described previously (20,27). Activation maps were overlaid on the T2 anatomic images for visual inspection.
Activation was quantified using the method proposed by Kleinschmidt (20,33). Images were masked to exclude non-brain voxels as well as excessively noisy voxels, i.e. those whose percent signal change across the functional run exceeded 6%, which most often included the eyes and major sinuses. The histogram of t-scores across the remaining voxels was calculated and smoothed, and a normal distribution was fit to the central peak to determine the background noise distribution, which was subtracted from the smoothed histogram to generate a new histogram of activated voxels. An activation threshold was defined as t>2, and the activation volume (spatial extent) and magnitude were calculated as the number of suprathreshold voxels in the activated distribution and their average t-score value, respectively. Significant differences in activation volume and magnitude between the VD and AR2 sequences were examined using paired, two-tailed t-tests.
To compare the noise sensitivities of the trajectories, a temporal signal-to-fluctuation noise ratio (SFNR, (27)) calculation was performed on the resting-state data. For a given subject, an ROI was defined by taking the union of activated (t>2) voxels across the AR2 and VD runs of the sensory task; the time series of all voxels within this ROI were extracted from the subject’s resting-state scan, and SFNR values were computed by dividing the mean of each voxel’s time series by its standard deviation. For comparison, the SNR for a single-shot spiral-in/out sequence with Archimedean sampling (TE/TR=67ms/3.4s, flip angle=86°) was calculated and measured experimentally for one of the 4 subjects.
The activation volumes and magnitudes for both comparisons are provided for all subjects in Table 3. Figure 5 shows thresholded t-maps of 3 individual subjects for the SM-block task (comparison #1), and Figure 6 shows thresholded t-maps of 4 individual subjects for the SM-event task (comparison #2). In comparison #1, in which the VD and AR2 sequences had identical temporal resolution and number of time frames, there was no significant difference in mean t-score magnitudes within the activated voxels of the two sequences (t(8)=0.07, p=0.95), nor was there a difference in the spatial extent of activation (t(8)=1.99, p=0.08). As is apparent from Fig. 5, some subjects exhibited higher t-scores for the VD sequence at a given spatial location (Fig. 5a,b), while others exhibited the opposite effect (Fig. 5c).
Comparison #2, in which the two sequences were set to their respective maximum possible number of time frames per scan (as well as maximum temporal resolution), yielded significantly greater activation magnitudes for the VD sequence compared to the AR2 sequence (t(4)=4.12, p=0.015; Table 3, Fig. 6). There was no significant difference in the spatial extent of activation between the two sequences (t(4)=1.32, p=0.26).
The measured SNR values within the functional ROIs for each of the 4 subjects who performed the resting-state scan are shown in Table 1.
The present work introduces the use of variable-density spiral k-space trajectories for fMRI, and specifically examines their utility in high resolution imaging, where the long readout duration of conventional trajectories often dictate multi-shot acquisitions. In the proposed design, k-space is critically sampled with an Archimedean spiral between the center and a user-defined radius, and extended to the desired spatial resolution with variable-density undersampling. fMRI experiments using sensory stimuli indicate that a single-shot spiral-in/out sequence with this variable-density trajectory is more efficient, and at least as effective, as a conventional 2-shot spiral-in/out sequence in detecting neural activation at high (128×128) resolution.
By reducing the readout duration, the variable-density trajectory allows for reduced sensitivity to motion and off-resonance, improved temporal resolution (and more time frames) for a given spatial coverage, and/or better spatial coverage for a given temporal resolution. Another key advantage of the shortened readout duration of the variable-density spiral is the feasibility of including a spiral-in component with a reasonable TE for single-shot high resolution fMRI, which can substantially mitigate susceptibility dropout in brain regions near air-tissue interfaces (20,34). Current single-shot spiral fMRI at higher resolutions is often performed using a spiral-out trajectory (35,36) as a fully Archimedean spiral-in trajectory would necessitate a TE at or beyond the limit of optimal BOLD sensitivity, and an acquisition window of approximately twice that value.
With the variable-density trajectory, there is a tradeoff between the duration of the readout and the severity of high spatial-frequency undersampling. The incomplete sampling of high spatial frequencies raises concerns that aliased signal may obscure true activation or generate false positives. Such effects were not apparent in the current experiments (Figs. (Figs.55 and and6),6), perhaps due to the relatively conservative choice of variable-density parameters (r, α) and because the energy in the high spatial frequencies is small. While it is not possible to gauge absolutely the accuracy of an activation map, the robust localization of activation to the sensory cortices, along with the relatively high correspondence between the maps obtained with the VD and AR2 sequences on a single-subject basis despite potential across-scan variability, suggest that the effect is minor. In addition, though we used a fixed-width convolution kernel for the gridding reconstruction, it has been shown that varying the kernel extent across k-space as a function of the sampling density can further reduce aliasing energy in variable-density trajectories (37).
The measured SNR values agreed closely with the theoretical calculation; the SNR of the VD sequence was approximately 15% less than that of the AR2, and the single-shot Archimedean sequence had an SNR reduction of over 30%, primarily from the longer TE necessitated by the associated readout duration of the spiral-in component. Despite the ~15% reduction in uniform-brain SNR of the VD sequence compared to AR2, the fMRI activation results indicate that the contrast-to-noise ratio in fact improved in the SM-event task, likely because of the greater fidelity of the more rapid single-shot acquisition in sampling the hemodynamic response.
In Comparison 2, the VD sequence yielded equal or larger spatial activation volumes in 4 of the 5 subjects, but the advantage did not survive a group-level t-test. Mean t-scores within suprathreshold voxels were, however, significantly greater at the group level with the VD sequence. For some subjects, the two sequences paradoxically yielded comparable activation magnitudes and yet quite different activation volumes (e.g. Subjects 1 and 5; Table 3). In interpreting this result, however, note that larger activation volumes (i.e. number of voxels exceeding a predefined statistical threshold) do not necessarily imply higher mean t-scores, and vice versa. For Subject 5, fewer voxels were detected as activated for VD compared to AR, but among those voxels that did surpass the threshold, the mean t-score value was actually higher than that of the AR scan. Beyond detecting whether a voxel is “activated”, the t-score magnitude is also important as it affects the outcome of group-level voxel-wise statistics.
The fMRI experiments herein used parameter values of r=0.5 and α=3.6 for the VD trajectory. In general, the selection of these parameters can be guided by the desired readout time/TE and minimum FOV, as well as by the PSFs of their associated trajectories. As seen in Fig. Fig.2b2b and and33 and discussed above, varying (r, α) results in different distributions of aliased energy. The intermediate values (0.5, 3.6) were chosen here for demonstration purposes; however, preliminary experiments revealed that robust fMRI activation was achievable across a range of parameter settings. Allowing r to vary from 0.4 to 0.6 in studies of a small number of subjects did not produce results whose differences exceeded natural inter-session variability, though further studies with larger subject populations might indeed reveal an effect. Reducing the TE from 48 ms to 42 ms in one subject was also found to yield comparable results. Optimal values in this parameter space are likely to vary with the particular spatial and temporal resolution requirements, k-space characteristics of the imaged sample, and brain regions of interest in the fMRI study.
Here, the VD sequence was compared against a 2-shot fully-Archimedean spiral-in/out sequence; for the latter, multiple (2) shots were used in order to attain acceptable image quality at the desired spatial resolution. Yet, whether it is the most appropriate sequence against which to compare the VD method is perhaps open to discussion. A single-shot fully Archimedean sequence with a spiral-out readout (with TE=30ms) was included in preliminary experiments but was, as expected (20), found to produce inferior results compared to the two-shot spiral-in/out, and with excessive signal dropout in air/tissue interfaces. An echo-planar imaging sequence with parallel imaging acceleration was also considered as a point of comparison, and preliminary studies were performed. However, this sequence was ultimately not included in the study due to the SNR disadvantages of accelerated methods at high spatial resolution where thermal noise dominates.
While the present study focused on the use of variable-density spiral for efficient high resolution fMRI, there are other potential applications and alternate implementations. The reduction in readout time may prove useful for signal recovery at higher field strengths (e.g. 7 Tesla) where T2* is severely reduced; at conventional field strengths and lower spatial resolutions, it may reduce signal dropout in brain regions that are compromised by susceptibility-induced field gradients. As an illustration of the latter, Figure 7 shows the SFNR of the spiral-out images for 2 slices acquired at 2.5 × 2.5 mm2 in-plane resolution with a single-shot Archimedean spiral-in/out trajectory (minimum TE=39 ms) and a single-shot variable-density spiral-in/out trajectory (minimum TE = 29 ms). Improved SNR recovery in the orbitofrontal region is observed with the variable-density sequence, due primarily to the shorter TE. In addition, the proposed variable-density design may be implemented in interleaved 2D sequences or in 3D fMRI sequences. 3D methods have been shown to be effective for high-resolution fMRI, where thermal noise dominates physiological noise (38).
Variable-density spiral sampling appears to be a promising technique for increasing the efficiency and performance of high-resolution fMRI. In addition to achieving a reduction in readout time and avoiding artifacts associated with multi-shot approaches, the variable-density trajectory permits single-shot spiral-in/out acquisition at high resolution, which is beneficial for signal recovery in regions affected by susceptibility dropout. The parameterization and design of the proposed trajectory can flexibly accommodate the desired readout duration, sampling density, and minimum FOV, and is suitable for implementation in conventional resolution, 3D, and interleaved fMRI sequences in addition to the single-shot high-resolution design presented here.
The authors gratefully acknowledge support from NIH grants F31-AG032168 (to CC) and P41-RR009784 (to GHG). We thank Prachi Pandit for helpful comments on a draft of this manuscript, and Jonathan Winawer for assistance with preliminary retinotopic mapping experiments (not included here).