 Technical notes
 Open Access
 Published:
Multipoint 5D flow cardiovascular magnetic resonance  accelerated cardiac and respiratorymotion resolved mapping of mean and turbulent velocities
Journal of Cardiovascular Magnetic Resonance volume 21, Article number: 42 (2019)
Abstract
Background
Volumetric quantification of mean and fluctuating velocity components of transient and turbulent flows promises a comprehensive characterization of valvular and aortic flow characteristics. Data acquisition using standard navigatorgated 4D Flow cardiovascular magnetic resonance (CMR) is timeconsuming and actual scan times depend on the breathing pattern of the subject, limiting the applicability of the method in a clinical setting.
We sought to develop a 5D Flow CMR framework which combines undersampled data acquisition including multipoint velocity encoding with lowrank image reconstruction to provide cardiac and respiratorymotion resolved assessment of velocity maps and turbulent kinetic energy in fixed scan times.
Methods
Data acquisition and datadriven motion state detection was performed using an undersampled Cartesian tiny Golden angle approach. Locally lowrank (LLR) reconstruction was implemented to exploit correlations among heart phases and respiratory motion states. To ensure accurate quantification of mean and turbulent velocities, a multipoint encoding scheme with two velocity encodings per direction was incorporated. Velocityvector fields and turbulent kinetic energy (TKE) were obtained using a Bayesian approach maximizing the posterior probability given the measured data. The scan time of 5D Flow CMR was set to 4 min.
5D Flow CMR with acceleration factors of 19 .0 ± 0.21 (mean ± std) and velocity encodings (VENC) of 0.5 m/s and 1.5 m/s per axis was compared to navigatorgated 2x SENSE accelerated 4D Flow CMR with VENC = 1.5 m/s in 9 subjects. Peak velocities and peak flow were compared and magnitude images, velocity and TKE maps were assessed.
Results
While net scan time of 5D Flow CMR was 4 min independent of individual breathing patterns, the scan times of the standard 4D Flow CMR protocol varied depending on the actual navigator gating efficiency and were 17.8 ± 3.9 min on average. Velocity vector fields derived from 5D Flow CMR in the endexpiratory state agreed well with data obtained from the navigated 4D protocol (normalized rootmeansquare error 8.9 ± 2.1%). On average, peak velocities assessed with 5D Flow CMR were higher than for the 4D protocol (3.1 ± 4.4%).
Conclusions
Respiratorymotion resolved multipoint 5D Flow CMR allows mapping of mean and turbulent velocities in the aorta in 4 min.
Background
The hemodynamic assessment of blood flow through heart valves and in the great vessels using 4D Flow cardiovascular magnetic resonance (CMR) has received increasing attention and a number of important applications have emerged [1]. While standard parameters including volume flow and peak velocities have been widely used, a number of additional readouts have been suggested. Among these readouts, probing the intravoxel velocity standard deviation [2] holds considerable promise as it permits estimating turbulent contributions in blood flow [3,4,5].
A number of clinical studies have indicated the added value of measuring turbulent velocities alongside mean velocities using 4D Flow CMR [2, 6, 7], e.g. by the assessment of turbulent kinetic energy (TKE) which indicates how much energy is stored in turbulent flow and which has shown elevated values for patients with bicuspid aortic valves or a dilated aorta relative to healthy subjects [7]. Ongoing studies promise further insights into the clinical relevance of turbulent flow estimation [8].
A key challenge of encoding mean and turbulent velocities relates to achieving appropriate measurement sensitivity to both quantities. For the expected ranges of mean and turbulent velocities invivo, multiple velocity encodings per spatial direction are preferably obtained and combined during image reconstruction [4]. This approach, however, requires additional data and hence data undersampling is necessary to arrive at clinically acceptable scan times.
Significant scan time reductions in 4D Flow CMR have been achieved using combinations of parallel imaging and compressed sensing [9], kt methods [10, 11], radial [12] or spiral [13] acquisitions. Various modifications to the way data are sampled and reconstructed have been described, e.g. [14,15,16,17,18,19]. In general, however, scans remain too long to be performed in breathholds. Moreover, data acquisition without respiratory motion compensation leads to decreased image quality [14]. Accordingly, navigatorbased respiratory gating is typically employed [20]. Unfortunately, respiratory gating makes scan times unpredictable, as gating efficiencies vary among subjects and often even during a single scan. Clinically, this may lead to conflicts as scan slots have typically fixed durations.
Instead of using respiratory gating, respiratorymotion resolved imaging [21,22,23,24] has been proposed to address respiratorymotion related image artifacts. To this end, data are acquired continuously throughout the respiratory motion cycle and sorted afterwards during the image reconstruction task into discrete motion states based on a respiratory motion signal surrogate [23, 25]. Accordingly, data correlation not only along the cardiac phase dimension but also along the respiratory motion dimension can be exploited.
Data correlations lead to lowrank properties of the decoding problem to be solved in image reconstruction. Various studies have demonstrated that these lowrank characteristics are most efficiently exploited using local rather than global lowrankedness using patchbased decomposition of the multidimensional data [26, 27].
The objective of the present work was to implement and test respiratorymotion resolved Bayesian multipoint 5D Flow CMR for mapping both mean and turbulent velocities in the aorta with a fixed scan time.
Methods
Data acquisition and respiratory motion binning
The data acquisition and motion binning process is illustrated in Fig. 1. Data are sampled using a pseudoradial Cartesian sampling pattern which maps radial spokes onto a Cartesian grid [28, 29]. Using the Golden angle principle [30, 31], spatial incoherence of temporally adjacent frames is ensured. To reduce eddy current effects and keep the acoustic noise level during the measurement to a minimum, the tiny Golden angle concept was employed to avoid large jumps in kspace [31].
Datadriven respiratory motion detection [32] was performed as illustrated in Fig. 2a) using repetitive sampling of profiles through the kspace center d(k_{x}, t) = d(k_{x}, k_{y} = 0, k_{z} = 0, t) as part of the pseudoradial Cartesian sampling pattern. Upon inverse Fourier transform of d(k_{x}, t), the projection X(x, t) of the excited volume resolved along the frequency encode direction x and over time t was obtained, showing signal modulation depending on respiratory and cardiac motion. Given TR = 3.3 ms and with one in 21 readouts through the kspace center, a respiratory sampling rate of ca. 14 Hz was obtained. The signal of each coil channel was treated as a separate signal \( {\mathrm{X}}_{\mathrm{i}}\in {\mathbb{C}}^{N_x\times {N}_t};i\in \left\{1,2,\dots, nCoils\right\} \), with N_{x} being the number of samples in readout dimension and N_{t} the number of readouts through the kspace center. Only the magnitude was considered, as the main variation in signal phase was found to be related to the different velocity encodings which were alternated per heart beat (beatinterleaved) and which led to phase variations which were too large to be eliminated by lowpass filtering.
Singular value decomposition was performed and the 2nd singular vector of abs(X^{T}X) was selected, as it provided the predominant motion for each coil channel (Fig. 2b). To isolate respiratory motion, a lowpass filter was applied to limit the spectrum to frequencies below 0.6 Hz (Fig. 2c). Next, coil channel clustering [28] was performed to determine the predominant dynamics among all channels which were assumed to be the respiratory motion signal (Fig. 2d). Finally, the average over all these signals was calculated to provide the final respiratory motion signal surrogate (Fig. 2e).
Motionresolved image reconstruction
Using the respiratory motion surrogate signal, the acquired data were distributed across discrete motion bins. A total of four motion bins was defined and the width of each motion bin was set to achieve a similar undersampling factor for each respiratory motion state and cardiac phase by setting the bin boundaries to the 0.25, 0.5 and 0.75quantile of the respiratory motion signal. For some combinations of respiratory motion state, heart phase and velocity encoding, this approach led to very low sampling rates. Therefore, a data sharing approach [33, 34] was used to fill up frames with acceleration factors higher than 20 with data from neighboring respiratory motion states.
Accordingly, the final data vector reads \( \boldsymbol{d}\in {\mathbb{C}}^{N_x{N}_c{N}_{hp}{N}_{k_v}{N}_{rs}} \), with N_{x} being the number of samples in the spatial frequency domain, N_{c} denoting the number of coils, N_{hp} heart phases, \( {N}_{k_v} \) velocity encodings and N_{rs} respiratory motion states. A locallylowrank (LLR) [26, 27] approach was used to exploit correlations among heart phases (hp) and respiratory motion states (rs). Imagedomain data for each velocity encoding \( {\boldsymbol{i}}_{kv}\in {\mathbb{C}}^{N_x{N}_{hp}{N}_{rs}} \) were reconstructed separately as:
with the undersampling operator Ω, Fourier transform \( \mathcal{F} \), coil sensitivities \( \mathcal{S} \), kspace data d_{kv} and regularization weight λ. As illustrated in Fig. 3a), the operator \( {\mathcal{R}}_b \) selects the bth out of N_{b} blocks of size n_{x} × n_{y} × n_{z} in the image from all N_{hp} heart phases and N_{rs} respiratory motion states and transforms them into a Casorati matrix with dimensions n_{x}n_{y}n_{z} × N_{hp}N_{rs}. By penalizing the nuclear norm of this local Casorati matrix, lowrankedness is enforced locally.
For comparison, respiratorymotion resolved data were also reconstructed by penalizing total variation (TV) along the respiratory motion states and heart phases [22]:
where λ_{hp} and λ_{rs} denote the regularization weights along cardiac phases and respiratory motion states, respectively.
Multipoint velocity encoding and decoding
Assuming a Gaussian distribution of intravoxel velocities, the following signal model of mean velocity \( \overline{\boldsymbol{v}} \) and fluctuating component with standard deviation σ [2] was used:
where k_{v} is related to the first gradient moment of a bipolar velocity encoding gradient by \( {\boldsymbol{k}}_{\boldsymbol{v}}=\gamma {\int}_0^Tt\boldsymbol{G}(t) dt \), with T being the time of application of gradient waveform G(t) and determines the encoding velocity as \( {v}_{enc}=\frac{\pi }{k_v} \).
The signal model implies a tradeoff between a sufficiently high encoding velocity v_{enc} to avoid phase wraps [1], and a sufficiently high encoding efficiency of the fluctuating velocity components [5] such that v_{enc}~σ. It has been demonstrated that this tradeoff may lead to insufficient sensitivity to fluctuating velocity components and hence a multipoint encoding scheme using two instead of one k_{v} encoding point per axis was implemented [4]. Using the probability map \( P\left(\overline{\boldsymbol{v}},\boldsymbol{\sigma} D,S\right) \), pixelwise estimates of mean velocity \( \overline{\boldsymbol{v}} \) and intravoxel standard deviation σ given the model S and the measured image data D were obtained [4]. Data acquired with different k_{v}′ s were combined by choosing the values for \( \overline{\boldsymbol{v}} \) and σ maximizing the posterior probability of the measured data:
An illustration of multipoint acquisition and Bayesian decoding is provided in Fig. 3b).
Invivo study
All invivo work was performed upon written informed consent of the subjects and according to local ethics regulations. Data of the ascending aorta of nine healthy subjects (30 ± 11 years, 6 male) were acquired on a 3 T CMR system (Ingenia; Philips Healthcare, Best, the Netherlands) using the proposed 5D Flow CMR and a standard, navigatorgated 4D Flow CMR approach [1] with a spatial resolution of 2.5 × 2.5 × 2.5mm^{3}, T_{E} = 3.3 ms, T_{R} = 4.9 ms, flip angle = 8 and 25 cardiac phases. The different velocity encodings were alternated per heart beat. The image acquisition matrix was N_{readout} = 100 ± 9, N_{phase} = 101 ± 8 and N_{slice} = 20 ± 1 (mean ± std). Using the proposed method, velocities were encoded with v_{enc} = 50 cm/s and 150 cm/s per axis and an additional v_{enc} = 0 cm/s measurement. Scan time was fixed to 4 min. After sorting the data into 4 respiratory motion states, acceleration factors were 19.0 ± 0.21 (mean ± std). On average the viewsharing approach led to 61% of the total number of acquired samples being shared among motion bins in the study cohort.
Standard 4D Flow CMR was recorded with a single venc (v_{enc} = 150 cm/s) per axis (4D Flow Ref). Twofold acceleration using parallel imaging [35] and standard pencilbeam navigator gating on the diaphragm (5 mm gating window) was employed with the 4D Flow reference protocol. Accordingly, the effective scan time depended on the breathing pattern of the subjects. All measurements were performed with retrospective cardiac gating and a 28channel receiver coil.
Prior to image reconstruction, the 28channel data were compressed to 8 virtual channels [36]. The locally lowrank (LLR) reconstruction (Equation 1) was solved using an implementation of the fast iterative shrinkagethresholding algorithm (FISTA) [37] provided with the Berkeley advanced reconstruction toolbox (BART) [38]. A patch size of n_{x} = n_{y} = n_{z} = 8 was used for the locally lowrank constraint. The TV reconstruction (Equation 2) was performed with an adapted version of the code provided with [22]. The standard 4D Flow reference data were reconstructed with MRecon (GyroTools LLC, Zurich, Switzerland). Sensitivity maps were estimated from a separate scan using the ESPIRiT method [39].
The regularization hyperparameters λ in Equation 1 and λ_{hp} and λ_{rs} in Equation 2 were optimized for best agreement of velocity components in systole with data of the 4D Flow reference measurement:
Equation 4 was minimized using the bayesopt function in MATLAB R2017b (The MathWorks, Natick, Massachusetts, USA) [40]. Resulting hyperparameters were λ = 0.01 for Equation 1 and λ_{hp} = 0.04 and λ_{rs} = 0.05 for Equation 2. For both optimizations, the number of iterations was set to 80.
Data analysis
Upon image reconstruction, concomitant gradient terms were corrected [41] and thirdorder background phase correction was applied [42, 43]. Segmentation of the aorta was performed using ITKSNAP [44]. The same masks were used for 4D Flow reference data and data in expiration obtained with 5D Flow LLR and 5D Flow TV. Maximum intensity projections (MIP) of velocity magnitude and TKE maps in systole were calculated as:
In order to avoid noise at the vessel boundaries the segmentation masks were eroded by one pixel for MIP projections.
The normalized rootmeansquare error (nRMSE) of the velocity magnitude in the segmented aorta was calculated according to:
where ROI corresponds to the number of voxels in the segmentation mask.
TKE was calculated voxelwise as:
where \( {\sigma}_j^2 \) denotes the variance of velocity deviations per direction j and ρ is the fluid density (here we assumed 1060 kg/m^{3} for blood).
Moreover, multiplanar reslicing along the aorta was conducted using VMTK (www.vmtk.org) and peak flow and peak throughplane velocities (i.e. the maximum velocity component in the corresponding crosssection of the aorta) were assessed for each plane. Peak throughplane velocities were determined upon application of a 3x3x3 median filter to reduce contributions from noise.
The effectiveness of respiratorymotion resolved 5D Flow CMR was assessed by comparing it to reconstructions using only data from the endexpiratory bin (EXP) and to reconstructions using no respiratory binning (NG).
To assess the performance of navigator gating in the standard 4D Flow CMR protocol, the positions measured by the navigator on the diaphragm were interpolated and each readout was attributed its corresponding navigator position. Then, mean and standard deviations of the acceptance window per heart phase were calculated. Due to the coarse temporal resolution of the navigator signal (one navigator per heart beat), linear interpolation was chosen in order to provide a lower bound on the width of the acceptance window for different cardiac phases.
Statistical analysis
BlandAltman analysis [45] of peak flow and peak throughplane velocities was performed to assess the agreement of velocity fields obtained with 5D Flow LLR and the 4D Flow reference.
Results
Scan durations for 4D Flow CMR were 17.8 ± 3.7 min compared to 4 min for the 5D Flow protocol.
Figure 4 provides an analysis of the acceptance window for respiratory gating with one pencil beam navigator per cardiac cycle as being used in the 4D Flow reference measurements. If the navigator signal was within the acceptance window, data of all cardiac phases were accepted. This results in a wider acceptance window for later cardiac phases leading to corresponding uncertainty about the actual motion states and increasing image artifacts for later cardiac phases.
Reconstructions exploiting respiratory motion states (5D Flow LLR) are compared to using only data in endexpiration (EXP LLR) and without gating (NG LLR) in Fig. 5. Magnitude images and TKE maps from 5D Flow LLR reconstructions show less artifacts compared to EXP LLR and NG LLR reconstructions.
Magnitude and velocity magnitude images for the endexpiratory motion state obtained with 5D Flow LLR and 5D Flow TV are compared relative to the 4D Flow reference in Fig. 6. It can be seen that maximum intensity projections of velocity magnitudes are similar for all methods while TKE maps and magnitude images provided by 5D Flow TV show more residual aliasing when compared to 5D Flow LLR. The root mean square estimate (RMSE) of the reconstructed images was lower for 5D Flow LLR compared to 5D Flow TV for both magnitude and velocity magnitude (8.7 vs. 13.5% for and 7.9 vs 9.8%, respectively). Maps of TKE for all 9 subjects obtained with 5D Flow LLR and the 4D Flow reference are compared in Fig. 7. In general, TKE maps derived from the 4D Flow reference show a higher noise level compared to 5D Flow data.
Exemplary slices in the ascending and descending aorta show throughplane velocity profiles to be in agreement among the 4D Flow reference and 5D Flow LLR whereas 5D Flow TV shows artifacts in the descending aorta.
Exemplary magnitude and velocity magnitude images for different respiratory bins reconstructed with 5D LLR are shown in Fig. 8 along with an exemplary respiratory curve. Bins in expiration are narrower as the duration of expiration was longer compared to inspiration. The displacement of the heart and aorta is clearly seen in both magnitude and velocity maps. While image quality is acceptable, there is degradation present from end expiration to end inspiration.
Figure 9 compares velocity fields, peak velocities and peak flow obtained for all 9 subjects. The nRMSE between the velocity magnitudes obtained with 5D Flow LLR and the 4D Flow reference was 8.9 ± 2.1%. On average, 5D Flow LLR shows good agreement of peak velocities and peak flow relative to the 4D Flow reference (peak velocities: 3.1 ± 4.4%, peak flow: −2.4 ± 6.9%).
Discussion
In this work, a respiratorymotion resolved Bayesian multipoint 5D Flow CMR approach has been implemented based on pseudoradial tiny Golden angle Cartesian sampling in conjunction with locally lowrank image reconstruction to map mean and turbulent velocities in the aorta in a fixed scan time of 4 min.
By exploiting data from all respiratory motion states, the duration of 5D Flow CMR becomes independent of the individual respiratory motion patterns of the subjects. In comparison to standard 4D Flow CMR protocols [1], 5D Flow CMR is about 4.5 times shorter on average. Whereas pencilbeam based navigator gating in standard 4D Flow CMR cannot ensure motionfree data across the entire cardiac cycle, as demonstrated in Fig. 4, repetitive datadriven respiratory motion detection and continuous data acquisition of 5D Flow CMR provide motion estimates throughout the cardiac cycle.
It has been demonstrated that the LLR model allows to improve reconstruction accuracy relative to using data from end expiration only, as shown in Fig. 5. Of note, lowrank models have been used to exploit correlations in dynamic CMR [46,47,48] and have already been proposed for reconstructing respiratorymotion resolved data. However, existing approaches confine the reconstruction to a subspace defined by basis functions along the respiratory motion dimension either by combining cardiac and respiratory motion in a Casorati matrix [21] or by means of lowrank tensor factorization [49, 50]. These approaches either require the acquisition of a separate scan to estimate a subspace or a subspace has to be estimated from undersampled kspace data, e.g. via lowrank tensor factorization [49]. The first option increases scan time and can lead to problems when there is a mismatch between both scans (e.g. bulk motion) while the second option adds complexity to the reconstruction as another hyperparameter needs to be optimized for subspace estimation prior to reconstruction of the actual image. In comparison, the LLR model proposed in the present work can be applied straightforwardly by penalizing the nuclear norm of the data rearranged in a Casorati matrix.
When comparing the 5D Flow LLR reconstruction with the 5D Flow TV approach, both methods showed good agreement of velocity data relative to the 4D Flow reference. However, TKE maps and magnitude images exhibit less noise when using the 5D Flow LLR approach. The inferior performance of the TV approach is associated with the piecewise constant solutions favored by TV which in turn leads to smoothing of peak values. To reduce underestimation of peak velocities, hyperparameters may be determined for best agreement with reference velocities. However, this led to nonoptimal results for magnitude images which therefore still showed a high degree of aliasing noise.
Invivo, good agreement of peak velocities and peak flow values of 5D Flow LLR with data derived from standard 4D Flow CMR was found. A small bias towards higher velocities was observed for 5D Flow LLR compared to the 4D Flow reference. This can be related to previous findings that respiratory motion leads to blurring of the image [20, 51] which corresponds to spatial lowpass filtering and is therefore likely to reduce peak velocities.
A multipoint encoding scheme [4] was incorporated into the 5D Flow approach. Thereby the velocity vector field was encoded with different sensitivities to mean and fluctuating velocities to provide an accurate assessment of mean and turbulent velocities over a large dynamic range. The 7point 5D flow LLR implementation yielded TKE maps as they can be expected for healthy, young volunteers, with moderate values of up to 300 J/m^{2} and turbulence mainly occurring in the region of the flow jet in the proximal aorta. In comparison, TKE maps derived from the 4D Flow reference with 4point encoding appeared noisy, due to only using a single encoding velocity of 150 cm/s. Moreover, the 4D Flow reference showed TKE values in the descending aorta which were as high as in the ascending aorta. This is considered unrealistic as transient/turbulent flow typically occurs in the region downstream of the aortic valve and not in the descending aorta in healthy, young volunteers.
The acquisition time for 7point encoded 5D Flow was fixed to 4 min. For 10point encoding, as it has been used when assessing stenotic valves [7], this would correspond to an acquisition time of ca. 6 min for 5D Flow, compared to 17.2 ± 4.7 min for accelerated 10point 4D Flow acquisitions [7].
Image quality of data in inspiration as shown in Fig. 8 was lower compared to data reconstructed in the endexpiratory motion state with 5D Flow LLR. This can be explained with inspiration taking a smaller fraction of the duration of the respiratory cycle than expiration. Moreover, the respiratory motion curve does not reach the maximum value in inspiration for every respiratory cycle, whereas it reaches an endexpiratory plateau in most cases. As respiratory motion states were defined to obtain similar acceleration factors in each motion state, this led to a wider range of motion in inspiratory states than in expiratory states.
The present work assessed mean and fluctuating velocity vector fields in healthy subjects only. When applying the proposed 5D Flow method in patients with aortic stenosis, an extension to 10point encoding [7] needs to be implemented for sufficient dynamic range. Adding 3 more points would lead to an increase in scan time by ca. 42% and can thus still easily be integrated into clinical workflows.
While the present study was focused on turbulence encoding, the proposed technique can be readily applied for the acquisition of standard 4D flow CMR data in the aorta and in the heart. For different applications the scan time will change as a function of the size of the acquisition matrix, but respiratory motion is still expected to be accurately modelled as a lowrank problem. Therefore, fixing scan time independently of respiratory motion should still be feasible. Moreover, the suggested multipoint encoding approach provides improved accuracy for lower velocities, similar to other multiVENC approaches [52, 53]. In addition, the protocol can be used to measure flow in different respiratory motion states to assess respiration dependent flowpatterns, e.g. the Fontan circulation [54, 55]. To obtain better results in inspiration, one might define bins with equal range of motion, thus leading to varying acceleration factors for different bins. In order to still obtain sufficient sampling rates in each motion state, one would therefore need to increase scan time beyond the current 4 min limit.
A limitation of the present study is the lack of a ground truth data. The 4D Flow reference data were acquired in a separate scan and consistency of flow conditions cannot be expected over the entire scan session. In addition, the 4D Flow reference data showed considerable respiratorymotion related artifacts in later cardiac phases since the pencilbeam navigator was played out right after Rwave detection and only once per cardiac cycle (acceptance rates varied between 42 and 72%). Accordingly, the quantitative comparisons with the 4D Flow data were limited to systole in order to avoid respiratorymotion corrupted cardiac phases as much as possible.
Finally, a drawback of both 5D Flow LLR and 5D Flow TV are considerably longer reconstruction times when compared to standard 4D Flow CMR (ca. 35 min for 5D LLR and ca. 60 min for 5D TV vs. ca. 3 min for conventional 4D Flow CMR on a workstation with two 14 Core Intel Xeon E5–2680 CPUs and 256 GB RAM). In the future, this relative disadvantage may be addressed by employing variational neural network based reconstructions and their deployment on dedicated hardware [56, 57].
Conclusion
Respiratory motion resolved multipoint 5D Flow CMR allows for breathingpattern independent mapping of mean and turbulent velocities in 4 min. The reduction in scan time allows for integration of the sequence into standard clinical workflows. Further invivo studies are now warranted to assess the performance of the method in relevant patient populations.
Availability of data and materials
The data that support the findings of this study are available from the corresponding author upon reasonable request subject to restriction on use by the Ethics Committee of the Canton of Zurich. A demo script with exemplary data will be provided online upon acceptance of this manuscript.
Abbreviations
 CMR:

Cardiovascular magnetic resonance
 CS:

Compressed Sensing
 EXP:

endexpiratory
 FISTA:

fast iterative shrinkage thresholding
 LLR:

Locallylowrank
 MIP:

Maximum Intensity Projection
 MP:

Multipoint
 NG:

nongated
 nRMSE:

normalized rootmeansquare error
 RMSE:

root mean square error
 TKE:

Turbulent kinetic energy
 TV:

Total Variation
 VENC:

Velocity encoding
References
 1.
Dyverfeldt P, Bissell M, Barker AJ, Bolger AF, Carlhäll CJ, Ebbers T, et al. 4D flow cardiovascular magnetic resonance consensus statement. J Cardiovasc Magn Reson. 2015;17(1).
 2.
Dyverfeldt P, Sigfridsson A, Kvitting JPE, Ebbers T. Quantification of intravoxel velocity standard deviation and turbulence intensity by generalizing phasecontrast MRI. Magn Reson Med. 2006;56:850–8.
 3.
Dyverfeldt P, Gårdhagen R, Sigfridsson A, Karlsson M, Ebbers T. On MRI turbulence quantification. Magn Reson Imaging. 2009;27:913–22. https://doi.org/10.1016/j.mri.2009.05.004.
 4.
Binter C, Knobloch V, Manka R, Sigfridsson A, Kozerke S. Bayesian multipoint velocity encoding for concurrent flow and turbulence mapping. Magn Reson Med. 2013;69:1337–45.
 5.
Knobloch V, Binter C, Gülan U, Sigfridsson A, Holzner M, Lüthi B, et al. Mapping mean and fluctuating velocities by Bayesian multipoint MR velocity encodingvalidation against 3D particle tracking velocimetry. Magn Reson Med. 2014;71:1405–15.
 6.
Dyverfeldt P, Hope MD, Tseng EE, Saloner D. Magnetic resonance measurement of turbulent kinetic energy for the estimation of irreversible pressure loss in aortic stenosis. JACC Cardiovasc Imaging. 2013;6:64–71.
 7.
Binter C, Gotschy A, Manka R, Kozerke S. Turbulent kinetic energy assessed by multipoint 4D flow MRI provides additional information relative to echocardiography for the determination of aortic stenosis severity. Circ Cardiovasc Imaging. 2017.
 8.
Ha H. Kvitting JPE, Dyverfeldt P, Ebbers T. 4D flow MRI quantification of blood flow patterns, turbulence and pressure drop in normal and stenotic prosthetic heart valves. Magn Reson Imaging. 55 July 2018;2019:118–27. https://doi.org/10.1016/j.mri.2018.09.024.
 9.
Tariq U, Hsiao A, Alley M, Zhang T, Lustig M, Vasanawala SS. Venous and arterial flow quantification are equally accurate and precise with parallel imaging compressed sensing 4D phase contrast MRI. J Magn Reson Imaging. 2013;37:1419–26.
 10.
Giese D, Wong J, Greil GF, Buehrer M, Schaeffter T, Kozerke S. Towards highly accelerated Cartesian timeresolved 3D flow cardiovascular magnetic resonance in the clinical setting. J Cardiovasc Magn Reson. 2014;16:42.
 11.
Schnell S, Markl M, Entezari P, Mahadewia RJ, Semaan E, Stankovic Z, et al. Kt GRAPPA accelerated fourdimensional flow MRI in the aorta: effect on scan time, image quality, and quantification of flow and wall shear stress. Magn Reson Med. 2014;72:522–33.
 12.
Gu T, Korosec FR, Block WF, Fain SB, Turk Q, Lum D, et al. PC VIPR: a highspeed 3D phasecontrast method for flow quantification and highresolution angiography. Am J Neuroradiol. 2005;26:743 LP–749 http://www.ajnr.org/content/26/4/743.abstract.
 13.
Sigfridsson A, Petersson S, Carlhäll CJ, Ebbers T. Fourdimensional flow MRI using spiral acquisition. Magn Reson Med. 2012;68:1065–73.
 14.
Bollache E, Barker AJ, Dolan RS, Carr JC, van Ooij P, Ahmadian R, et al. Kt accelerated aortic 4D flow MRI in under two minutes: feasibility and impact of resolution, kspace sampling patterns, and respiratory navigator gating on hemodynamic measurements. Magn Reson Med. 2018;79:195–207.
 15.
Valvano G, Martini N, Huber A, Santelli C, Binter C, Chiappino D, et al. Accelerating 4D flow MRI by exploiting lowrank matrix structure and hadamard sparsity. Magn Reson Med. 2017;78:1330–41. https://doi.org/10.1002/mrm.26508.
 16.
Kim D, Dyvorne HA, Otazo R, Feng L, Sodickson DK, Lee VS. Accelerated phasecontrast cine MRI using kt SPARSESENSE. Magn Reson Med. 2012;67:1054–64. https://doi.org/10.1002/mrm.23088.
 17.
Santelli C, Loecher M, Busch J, Wieben O, Schaeffter T, Kozerke S. Accelerating 4D flow MRI by exploiting vector field divergence regularization. Magn Reson Med. 2016;75:115–25.
 18.
Rich A, Potter LC, Jin N, Ash J, Simonetti OP, Ahmad R. A Bayesian model for highly accelerated phasecontrast MRI. Magn Reson Med. 2016;76:689–701.
 19.
Kwak Y, Nam S, Akçakaya M, Basha TA, Goddu B, Manning WJ, et al. Accelerated aortic flow assessment with compressed sensing with and without use of the sparsity of the complex difference image. Magn Reson Med. 2013;70:851–8.
 20.
Dyverfeldt P, Ebbers T. Comparison of respiratory motion suppression techniques for 4D flow MRI. Magn Reson Med. 2017;78:1877–82.
 21.
Sigfridsson A, Wigström L, Kvitting JPE, Knutsson H. Kt2 BLAST: exploiting spatiotemporal structure in simultaneously cardiac and respiratory timeresolved volumetric imaging. Magn Reson Med. 2007;58:922–30.
 22.
Feng L, Axel L, Chandarana H, Block KT, Sodickson DK, Otazo R. XDGRASP: Goldenangle radial MRI with reconstruction of extra motionstate dimensions using compressed sensing. Magn Reson Med. 2016;75:775–88.
 23.
Cheng JY, Zhang T, Alley MT, Uecker M, Lustig M, Pauly JM, et al. Comprehensive multidimensional MRI for the simultaneous assessment of cardiopulmonary anatomy and physiology. Sci Rep. 2017;7:6–9.
 24.
Feng L, Coppo S, Piccini D, Yerly J, Lim RP, Masci PG, et al. 5D wholeheart sparse MRI. Magn Reson Med. 2018;79:826–38. https://doi.org/10.1002/mrm.26745.
 25.
Bastkowski R, Weiss K, Maintz D, Giese D. Selfgated goldenangle spiral 4D flow MRI. Magn Reson Med. 2018;80:904–13.
 26.
Trzasko J, Manduca A, Borisch E. Local versus global lowrank promotion in dynamic MRI series reconstruction. In: Proceedings of the 19th annual meeting of ISMRM; 2011. p. 4371.
 27.
Zhang T, Pauly JM, Levesque IR. Accelerating parameter mapping with a locally low rank constraint. Magn Reson Med. 2015;73:655–61.
 28.
Zhang T, Cheng JY, Chen Y, Nishimura DG, Pauly JM, Vasanawala SS. Robust selfnavigated body MRI using dense coil arrays. Magn Reson Med. 2016;76:197–205.
 29.
Cheng JY, Hanneman K, Zhang T, Alley MT, Lai P, Tamir JI, et al. Comprehensive motioncompensated highly accelerated 4D flow MRI with ferumoxytol enhancement for pediatric congenital heart disease. J Magn Reson Imaging. 2016;43:1355–68.
 30.
Winkelmann S, Schaeffter T, Koehler T, Eggers H, Doessel O. An optimal radial profile order based on the golden ratio for timeresolved MRI. IEEE Trans Med Imaging. 2007;26:68–76.
 31.
Wundrak S, Paul J, Ulrici J, Hell E, Rasche V. A small surrogate for the golden angle in timeresolved radial MRI based on generalized fibonacci sequences. IEEE Trans Med Imaging. 2015;34:1262–9.
 32.
Larson AC, Kellman P, Arai A, Hirsch GA, McVeigh E, Li D, et al. Preliminary investigation of respiratory selfgating for freebreathing segmented cine MRI. Magn Reson Med. 2005;53:159–68.
 33.
Rasche V, De Boer RW, Holz D, Proksa R. Continuous radial data acquisition for dynamic MRI. Magn Reson Med. 1995;34:754–61.
 34.
Zhang S, Block KT, Frahm J. Magnetic resonance imaging in real time: advances using radial FLASH. J Magn Reson Imaging. 2010;31:101–9.
 35.
Pruessmann KP, Weiger M, Scheidegger MB, Boesiger P. SENSE: sensitivity encoding for fast MRI. Magn Reson Med. 1999;42:952–62.
 36.
Buehrer M, Pruessmann KP, Boesiger P, Kozerke S. Array compression for MRI with large coil arrays. Magn Reson Med. 2007;57:1131–9.
 37.
Beck A, Teboulle M. A fast iterative shrinkagethresholding algorithm. Soc Ind Appl Math J Imaging Sci. 2009;2:183–202.
 38.
Tamir JI, Ong F, Cheng JY, Uecker M, Lustig M. Generalized magnetic resonance image reconstruction using the Berkeley advanced reconstruction toolbox. In: Proceedings of the ISMRM 2016 data sampling and image reconstruction workshop. 2016. p. 9660006. http://www.ismrm.org/workshops/Data16/.
 39.
Uecker M, Lai P, Murphy MJ, Virtue P, Elad M, Pauly JM, et al. ESPIRiT—an eigenvalue approach to autocalibrating parallel MRI: where SENSE meets GRAPPA. Magn Reson Med. 2014;71:990–1001.
 40.
Snoek J, Larochelle H, Adams RP. Practical Bayesian optimization of machine learning algorithms. Adv Neural Inf Process Syst 2012;:2951–2959. doi:1206.2944.
 41.
Bernstein MA, Zhou XJ, Polzin JA, King KF, Ganin A, Pelc NJ, et al. Concomitant gradient terms in phase contrast MR: analysis and correction. Magn Reson Med. 1998;39:300–8.
 42.
Busch J, Giese D, Kozerke S. Imagebased background phase error correction in 4D flow MRI revisited. J Magn Reson Imaging. 2017:1–10. https://doi.org/10.1002/jmri.25668.
 43.
Walker PG, Cranney GB, Scheidegger MB, Waseleski G, Pohost GM, Yoganathan AP. Semiautomated method for noise reduction and background phase error correction in MR phase velocity data. J Magn Reson Imaging. 1993;3:521–30.
 44.
Yushkevich, Paul A., Yang Gao, and Guido Gerig. "ITKSNAP: An interactive tool for semiautomatic segmentation of multimodality biomedical images." 2016 38th Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC). IEEE, 2016.
 45.
Altman DG, Bland JM. Measurement in medicine: the analysis of method comparison studies. Stat. 1983:307–17.
 46.
Liang ZP. Spatiotemporal imaging with partially separable functions. IEEE Int Symp Biomed Imaging. 2007;2:988–91.
 47.
Pedersen H, Kozerke S, Ringgaard S, Nehrke K, Kim W. Kt PCA: temporally constrained kt BLAST reconstruction using principal component analysis. Magn Reson Med. 2009;62:706–16.
 48.
Otazo R, Candès E, Sodickson DK. Lowrank plus sparse matrix decomposition for accelerated dynamic MRI with separation of background and dynamic components. Magn Reson Med. 2015;73:1125–36.
 49.
He J, Liu Q, Christodoulou A, Ma C, Lam F, Liang ZP. Accelerated highdimensional MR imaging with sparse sampling using lowrank tensors. IEEE Trans Med Imaging. 2016;0062(c):1–1. https://doi.org/10.1109/TMI.2016.2550204.
 50.
Christodoulou AG, Shaw JL, Nguyen C, Yang Q, Xie Y, Wang N, et al. Magnetic resonance multitasking for motionresolved quantitative cardiovascular imaging. Nat Biomed Eng. 2018;2:215–26. https://doi.org/10.1038/s415510180217y.
 51.
Wang Y, Riederer SJ, Ehman RL. Respiratory motion of the heart: kinematics and the implications for the spatial resolution in coronary imaging. Magn Reson Med. 1995;33:713–9.
 52.
Lee AT, Bruce Pike G, Pelc NJ. Threepoint phasecontrast velocity measurements with increased velocitytonoise ratio. Magn Reson Med. 1995;33:122–6. https://doi.org/10.1002/mrm.1910330119.
 53.
Herment A, Mousseaux E, Jolivet O, DeCesare A, Frouin F, ToddPokropek A, et al. Improved estimation of velocity and flow rate using regularized threepoint phasecontrast velocimetry. Magn Reson Med An Off J Int Soc Magn Reson Med. 2000;44:122–8.
 54.
Gewillig M. The Fontan circulation. Heart. 2005;91:839–46. https://doi.org/10.1136/hrt.2004.051789.
 55.
Jarvis K, Schnell S, Barker AJ, Rose M, Robinson JD, Rigsby CK, et al. Caval to pulmonary 3D flow distribution in patients with Fontan circulation and impact of potential 4D flow MRI error sources. Magn Reson Med. 2018; February 2017:1205–18.
 56.
Vishnevskiy V, Sanabria SJ, Goksel O. Image reconstruction via Variational network for realtime handheld soundspeed imaging. In: International workshop on machine learning for medical image reconstruction; 2018. p. 120–8.
 57.
Hammernik K, Klatzer T, Kobler E, Recht MP, Sodickson DK, Pock T, et al. Learning a variational network for reconstruction of accelerated MRI data. Magn Reson Med. 2018;79:3055–71.
Acknowledgements
The authors thank Giuseppe Valvano for his valuable advice regarding lowrank modelling of 4D flow and Richard Droste for contributing to initial experiments.
Funding
Financial support by the Swiss National Science Foundation SNF grant 320030_149567 and by ETH PASC is gratefully acknowledged.
Author information
Affiliations
Contributions
SK and JW conceived of the study. JW implemented the sequence and reconstructions. JW and HD performed data acquisition. All authors participated in revising the manuscript and read and approved the final manuscript.
Corresponding author
Correspondence to Jonas Walheim.
Ethics declarations
Ethics approval and consent to participate
The study was approved by the Ethics Committee of the Canton of Zurich, Switzerland, and all subjects provided written informed consent. Imaging was performed at the Zurich University Hospital, Zurich, Switzerland. Anonymized data was analyzed at ETH Zurich with approval by the mentioned authority. Written, informed consent was obtained before the experiment according to ethics approval and institutional guidelines.
Consent for publication
Written informed consent was obtained from all subjects for publication of their individual details and accompanying images in this manuscript. The consent forms are held in the subjects’ notes and are available for review by the EditorinChief.
Competing interests
The authors declare that they have no competing interests.
Additional information
Publisher’s Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Received
Accepted
Published
DOI
Keywords
 4D flow MRI
 Velocity mapping
 Turbulent kinetic energy
 Respiratory motion compensation
 Cartesian Golden angle
 Lowrank image reconstruction
 Cardiovascular magnetic resonance