WO2018137199A1

WO2018137199A1 - Real-time phase-contrast flow mri with low rank modeling and parallel imaging

Info

Publication number: WO2018137199A1
Application number: PCT/CN2017/072670
Authority: WO
Inventors: Aiqi SUN; Bo Zhao; Rui Li
Original assignee: Tsinghua University
Priority date: 2017-01-25
Filing date: 2017-01-25
Publication date: 2018-08-02
Also published as: CN109791617A; CN109791617B

Abstract

A novel model-based imaging method is proposed to enable high-resolution real time phase-contrast magnetic resonance imaging (PC-MRI) with sparse sampling. The real-time phase-contrast flow magnetic resonance imaging method comprises: acquiring real-time PC-MRI data, which includes training data and imaging data; performing a low-rank based image reconstruction, which includes estimation of temporal subspace and spatial subspace based on the acquired training data and imaging data; calculating velocity maps based on the reconstructed real-time flow images; and performing quantitative flow analysis based on the calculated velocity maps. The proposed method enables high-resolution real-time PC-MRI at 2D and for the first time at 3D without electro-cardiogram (ECG) gating and respiration control.

Description

REAL-TIME PHASE-CONTRAST FLOW MRI WITH LOW RANK MODELING AND PARALLEL IMAGING

TECHNICAL FIELD

The present invention relates， in general， to magnetic resonance imaging (MRI) ， and in particular， to high-resolution real time phase-contrast magnetic resonance imaging (PC-MRI) with sparse sampling.

BACKGROUND OF THE INVENTION

Over the past few decades， phase-contrast magnetic resonance imaging (PC-MRI) has been developed into a powerful tool for quantification and visualization of blood flow dynamics in the heart and large vessels [1-5] . It has advanced the understanding and diagnosis of various cardiovascular diseases， such as atherosclerosis [6] ， aneurysms [7] ， and arteriovenous malformation [8] . Conventional PC-MRI [9， 10] employs electro-cardiogram (ECG) synchronized cine acquisitions with respiration control to acquire data from multiple cardiac cycles， from which spatiotemporally averaged velocity maps are reconstructed. Although this approach is widely used in biomedical research and clinical practice， it suffers from a number of well-known limitations. For example， it often requires periodic or quasi-periodic cardiac motion to ensure efficient data acquisition； rejection of data caused by irregular cardiac motion often leads to prolonged acquisition times. Additionally， due to the underlying assumption， this approach only obtains averaged flowinformation over multiple cardiac cycles， failing to resolve beat-by-beat flow variations associated with irregular cardiac motion (e.g.， cardiac arrhythmia) . Capturing physiological and/or pathological flow variabilities has long been an important goal of PC-MRI research [11-14] .

Real-time PC-MRI [15-17] without ECG gating and respiration control is a promising direction to address these limitations； however， it requires a much higher imaging speed， posing significant challenges for both data acquisition and image reconstruction. A number of techniques have been developed to advance real-time PC-MRI in two spatial dimension (2D) with one velocity encoding. For example， advanced acquisition methods， such as echo-planar [18] ， radial [19] ， and spiral [20-23] acquisition schemes， have been employed for real-time PC-MRI. Besides， several acceleration methods by recovering images from undersampled data sets have been proposed. As with the emergence of parallel imaging， sensitivity encoding (SENSE) [24] and generalized autocalibrating partially parallel acquisitions (GRAPPA) [25] have also been applied [26-30] . And more recently， model-based reconstruction methods [31， 32] using regularized nonlinear inversion [33] have been developed， achieving 2D real-time flow imaging with a spatial resolution of 1.5 mm and a temporal resolution of 25.6 ms by jointly reconstructing a proton density map， a phase map and a set of coil sensitivities. Although so many kinds of methods have been developed to enable 2D real-time PC-MRI， 3D real-time flow MRI has not been reported yet. Considering the complex flow patterns and blood vessel geometry， it is highly desirable to perform 3D real-time flow imaging. However， 3D real-time flow imaging generally involves a more challenging trade-off between spatial resolution， temporal resolution， and imaging time， and a significantly more challenging computational problem.

In the present embodiment， we propose a new model-based method for real-time PC-MRI with sparse sampling. It is based on the integration of a novel low-rank model with parallel imaging. The proposed method enables 2D real-time PC-MRI without ECG gating and respiration control， and well resolves the beat-by-beat flow variations that cannot be obtained from the conventional cine method. And we have also extended this imaging technology to enable 3D real-time PC-MRI with three directional flow encodings by further incorporating a sparse modeling. Most importantly， we demonstrate the feasibility of 3D real-time PC flow MRI for the first time. The effectiveness of the proposed method has been systematically evaluated in 2D and 3D real-time PC-MRI using in vivo experiments.

SUMMARY OF THE INVENTION

An aspect of the present invention provides a real-time phase-contrast flow magnetic resonance imaging (MRI) method， comprising： acquiring real-time phase-contrast MRI (PC-MRI) data， which includes training data and imaging data； performing a low-rank based image reconstruction based on the acquired training data and imaging data； calculating velocity maps based on the reconstructed real-time flow images； and performing quantitative flow analysis based on the calculated velocity maps.

In the method according to the present disclosure， performing a low-rank based image reconstruction further comprises： performing a temporal interpolation on the training data； performing a temporal subspace estimation on the interpolated training data； performing an estimation of coil sensitivities on the acquired imaging data； performing a spatial subspace estimation based on the estimated temporal subspace， estimated coil sensitivities and the acquired imaging data.

The method according to the present disclosure could be used for a 2D and 3 D real-time PC-MRI.

BRIEF DESCRIPTION OF THE DRAWINGS

The accompanying figures， where like reference numerals refer to identical or functionally similar elements throughout the separate views， together with the detailed description below， are incorporated in and form part of the specification， and serve to further illustrate embodiments of concepts that include the claimed invention， and explain various principles and advantages of those embodiments.

FIG. 1 illustrates the proposed (k， t) -space sampling scheme according to the embodiment of the present disclosure；

FIG. 2 illustrates the pipeline for the proposed real-time PC-MRI method according to the embodiment of the present disclosure；

FIG. 3 shows comparisons of 2D real-time flow imaging with 2D cine flow imaging for two healthy subjects according to the embodiment of the present disclosure；

FIG. 4 shows reconstructed velocity waveforms from 2D real-time flow imaging for a healthy subject according to the embodiment of the present disclosure；

FIG. 5 shows Bland-Altman analysis between 2D real-time PC-MRI and 2D cine PC-MRI；

FIG. 6 shows 2D real-time PC-MRI for a 23-year-old arrhythmic patient according to the embodiment of the present disclosure；

FIG. 7 shows 2D real-time PC-MRI for a 72-year-old arrhythmic patient according to the embodiment of the present disclosure；

FIG. 8 shows velocity maps derived from the conventional 3D cine flow imaging and the proposed 3D real-time flow imaging for a healthy subject according to the embodiment of the present disclosure；

FIG. 9 shows the reconstructed flow waveforms from the proposed 3D real-time flow for a healthy subject according to the embodiment of the present disclosure；

FIG. 10 shows Bland-Altman analysis of peak velocities and stroke volumes between 3D real-time PC-MRI and 3D cine PC-MRI； and

FIG. 11 shows the reconstructed flow waveforms and pathline visualization derived from 3D Real-time PC-MRI for an arrhythmic patient according to the embodiment of the present disclosure.

Skilled artisans will appreciate that elements in the figures are illustrated for simplicity and clarity and have not necessarily been drawn to scale. For example， the dimensions of some of the elements in the figures may be exaggerated relative to other elements to help to improve understanding of embodiments of the present invention.

The method and device components have been represented where appropriate by conventional symbols in the drawings， showing only those specific details that are pertinent to understanding the embodiments of the present invention so as not to obscure the disclosure with details that will be readily apparent to those of ordinary skill in the art having the benefit of the description herein.

DETAILED DESCRIPTION OF THE INVENTION

Abbreviations

PC-MRI： phase-contrast magnetic resonance imaging；

ECG： electro-cardiogram；

SENSE， sensitivity encoding；

GRAPPA： generalized autocalibrating partially parallel acquisitions；

VENC： encoding velocity；

AAo： ascending aorta；

DAo： descending aorta.

At the end of the detailed description of the invention， a list of references is provided. Please note that all references as listed and mentioned in this disclosure are incorporated herein by reference in their entirely.

Embodiment： Real-time PC-MRI

Theory

Ignoring flow during readout time， the imaging equation for real-time PC-MRI can be modeled as follows：

d_v， i (k， t) ＝ ∫ S_i (r) ρ_v (r， t) e ^-j2πk·r dr +η_v， i (k， t) ， (1)

where ρ_v (r， t) denotes the dynamic image associated with either the flow-compensated (i.e.， v ＝ 1) or flow-encoded image sequence (i.e.， v ＝ 2， ...， N_v) ， S_i (r) the sensitivity map for the i th receiver coil (i ＝ 1， 2， ...， N_c) ， and d_v， i (k， t) and η_v， i (k， t) respectively the (k， t) -space measured data and measurement noise. In addition， v denotes different flow encodings (v ＝ 1： the flow-compensated data； v ＝ 2， ...， N_v： the flow-encoded data for each flow encoding direction) ， k the k-space sampling location， r the corresponding spatial location， and t the time instant. Here， the goal is to reconstruct ρ_v (r， t) from the undersampled data {d_v， i (k， t) } ， and then calculate the velocity maps as

where Δφ (r， t) ＝ ∠ρ_v (r， t) -∠ρ₁ (r， t) denotes the phase difference between the flow-encoded and flow-compensated image sequences， and VENC the pre-specified encoding velocity. Since in real-time PC-MRI， there is no data sharing with ECG gating， (k， t) -space data is often highly undersampled. Direct inversion of {d_v， i (k， t) } can incur significant aliasing artifacts and lead to inaccurate velocity measurements.

Here we introduce a low-rank model-based reconstruction method with parallel imaging to address the problem. For convenience， we consider a discrete image model， in which each flow image sequence can be represented as a spatiotemporal Casorati matrix [34] ， i.e.，

Similar to cardiac imaging applications [35-37] ， each C_v admits a low-rank approximation due to strong spatiotemporal correlation of time-series images. Moreover， due to the nature of flow encoding， there is also strong spatial and temporal correlation among different flow image sequences. To exploit such correlation， the following joint Casorati matrix is introduced：

on which we enforce the low-rank structure， i.e.， rank (C) ≤ L， where L refers to the rank of the matrix C， and satisfies L ≤ min (M， N) . There are a number of ways of imposing low-rank constraints [34， 36， 38， 39] . Here， we use an explicit rank constraint via matrix factorization， i.e.， C ＝ UV， where

and

In this low-rank representation， the columns of U and rows of V respectively span the spatial subspace and temporal subspace of C.

Next， we formulate the low-rank constrained reconstruction problem. First， note that with matrix-vector notation， Eq. (1) can be written as：

d_i ＝ Ω (F_sS_iC) + n_i， (4)

where d_i denotes the measured data， Ω the sparse sampling operator， F_s the spatial Fourier transform matrix， and S_i and n_i respectively the sensitivity map and measurement noise. Imposing the low-rank constraint， the image reconstruction problem can be formulated as

This problem is a non-convex optimization problem， for which a number of algorithms can be applied (e.g.， [40， 41 ] ) .

The image reconstruction problem can be further simplified. Extending the early work in cardiac imaging [34， 36， 37] ， we can pre-estimate the temporal subspace V by acquiring training data with a specialized data acquisition scheme. Specifically， as shown in FIG. 1 (a) ， we design an interleaved sampling pattern， in which both training data and imaging data are collected. Here， the training data are sampled from the central k-space， while the imaging data are acquired from the remaining (k， t) -space region with a random sampling scheme. With this sampling scheme， the two sets of data provide the complementary information for the low-rank model： the training data have high temporal resolution， while the imaging data have high spatial resolution. From the training data， we estimate the temporal subspace using the principal component analysis [34， 42] . With the imaging data， we estimate the spatial subspace U. To match the timing between the two sets of data， proper temporal interpolation scheme as shown in FIG. 1 (b) is applied， which interpolates the training data into those at the same time instants as the imaging data. Note that with such a scheme， the temporal resolution for the proposed method is 2 × N_v × TR. Moreover， note that the coil sensitivities S_i can be estimated from temporal averaged (k， t) -space data from the flow-compensated image sequence.

Withwe can determine U by solving the following convex optimization problem：

Due to the temporal subspace estimation， the low-rank matrix recovery problem has been reduced to a simple least-squares problem. By solving

the joint Casorati matrix can be reconstructed as

from which we can obtain each flow image sequence and estimate the flow velocities， where U /V /C denote the underlying true quantities， whereas

denote the corresponding reconstructions.

Eq. (6) could be sufficient for 2D real-time flow imaging. However， for 3D real-time flow imaging， it usually results in a significantly more challenging computational problem associated with highly-undersampled data. Here， we further utilize a sparsity constraint for solving the spatial subspace U considering that the joint Casorati matrix C often admits sparse representation in the spatial-spectral domain [36， 43] . Then the reconstruction problem for spatial subspace can be formulated as

where F_i denotes the temporal Fourier transform matrix， λ the regularization parameter， and vec (·) the operator concatenating the columns of a matrix into a vector.

A diagram summarizing the pipeline of the proposed method is shown in FIG. 2 according to the embodiment of the present disclosure. This pipeline consists of three major components： data acquisition， image reconstruction， and post processing.

Implementation

We performed in vivo studies to evaluate the performance of the proposed method for 2D and 3D real-time PC-MRI. The experiments were conducted on a 3.0 T whole body MR scanner (Achieva， Philips Medical System， Best， The Netherlands) ， equipped with a 32-channel cardiovascular coil. A gradient-echo (GRE) based pulse sequence was adapted to implement the proposed real-time acquisition scheme as shown in FIG. 1 (a) . Here neither ECG gating nor respiration control was used to aid data acquisition. Additionally， we performed conventional cine PC-MRI using a vendor-provided GRE-based pulse sequence， in which retrospective ECG gating was used.

In the 2D real-time PC-MRI experiments， ten healthy volunteers (7 males， age： 22-29 years， median： 25 years) ， who had no symptoms of cardiovascular diseases， were recruited. In addition， we performed a feasibility study of applying the proposed method for arrhythmia detection， and recruited two patients (2 males， age： 23-year old and 72-year old) . This study was approved by the Institutional Review Board at Tsinghua University， and all the subjects gave written informed consent. Both cine and real-time flow experiments were performed on the planes perpendicular to the ascending aorta (AAo) and descending aorta (DAo) during free breathing， and with one directional velocity encoding along the FH direction. For the cine acquisition， the retrospective ECG gating was set according to an estimate of each subject’s heartbeat period， and three averages were performed to mitigate respiratory motion artifacts. For both cine and real-time imaging experiments， the following imaging parameters were used： FOV ＝ 240 mm × 225 mm， matrix size ＝ 132 × 124， spatial resolution ＝ 1.80 mm × 1.80 mm， slice thickness ＝ 5 mm， TR/TE ＝ 4.5/2.8 ms， flip angle ＝ 10°， and VENC ＝ 200 cm/s. For the real-time flow imaging， the temporal resolution is 4 × TR ＝ 18 ms， whereas for the cine imaging， the temporal resolution is around 36 ms (with 28 cardiac phases) . The total acquisition times were around 94 s for both experiments.

Additionally， we also performed 3D real-time flow imaging experiments using the proposed method， in which five healthy volunteers and an arrhythmic patient were recruited. For comparison， we also acquired 3D cine flow imaging data with 2x SENSE [24] . Both cine and real-time experiments were conducted on the whole aorta during free breathing with the following parameters： FOV ＝ 180 × 256 × 43 mm³ (FH/AP/RL) ， spatial resolution ＝ 2.4 × 2.4 × 2.4 mm³， matrix size ＝ 76 × 108 × 18， repetition time/echo time ＝ 4.4/2.6 ms， flip angle ＝ 5°， and encoding velocity ＝ 200/150/150 cm/s (FH/AP/RL) ， and temporal resolution is 35 ms.

For cine flow imaging， the flow-compensated and flow-encoded images were simply reconstructed from the fully-sampled data. For the proposed real-time flow imaging， we followed the procedure in FIG. 2. Specifically， we first performed the temporal interpolation and estimated the temporal subspace V. We then estimated the coil sensitivity maps Si from the temporally averaged (k， t) -space measured data. Followed by forming the time-series images for flow-compensated and flow-encoded images， we further determined the spatial subspace U by solving Eq. (6) for 2D real-time imaging and by solving Eq. (7) for 3D real-time flow imaging. To potentially enhance the computational efficiency， proper coil compression (e.g.， [44] ) can be adopted. After image reconstruction， phase correction [45] was performed to correct the phase offsets caused by eddy currents. The velocity maps were then extracted for quantitative flow analysis.

For the in vivo experiments with healthy subjects， we evaluated the degree of agreement between the flow measurements from the cine approach and those from the proposed method. Specifically， we performed a Bland-Airman analysis， as well as a paired Student’s t-test， on the peak velocity and stroke volume obtained from the two methods. Here the peak velocity is defined as the maximum velocity within one cardiac cycle， and the stroke volume is the integral of the flow velocity over one cardiac cycle within the ascending aorta. For the experiments with arrhythmic patients， we evaluated the flow variabilities captured by the proposed method with reference to an external ECG recording of cardiac motion.

Results for 2D real-time PC-MRI

FIG. 3 shows comparisons of 2D real-time flow imaging with 2D cine flow imaging for two healthy subjects according to the embodiment of the present disclosure. Here， we show the reconstructed magnitude images and velocity maps corresponding to a systolic cardiac phase and a diastolic cardiac phase. As can be seen， the proposed method provides at least comparable reconstruction quality to the cine method. Although both methods can resolve the vessel structure， the real-time imaging method is more motion-robust than the cine method.

In addition， we analyzed the mean flow velocities associated with two regions of interest in AAo and DAo obtained from the two methods. FIG. 4 shows the reconstructed velocity waveforms according to the embodiment of the present disclosure for a healthy subject. Specifically， the velocity waveforms associated with the ascending aorta (AAo) and descending aorta (DAo) over 10 cardiac cycles are shown in FIGs. 4 (a) and (b) . Clearly， the proposed method well resolves beat-by-beat variations. We further evaluated how the velocity waveforms from the real-time imaging are related to the ones from the conventional cine method. We averaged the velocity waveforms over 30 consecutive cardiac cycles from the proposed method into one velocity waveform associated with a synthetic cardiac cycle， and then compared with that from the cine method. The comparisons for AAo are shown in FIG. 4 (c) and for DAo are shown in FIG. 4 (d) . It is evident that the averaged velocity vs.time profiles for AAo and DAo correlate well with the conventional cine method. In particular， both methods yield very similar peak velocities for the AAo and DAo.

We also performed statistical analysis of the results from the two methods for all ten healthy subjects. FIG. 5 shows the Bland-Altman plots of peak velocities and stroke volumes between the two methods. Specifically， Bland-Altman analysis of peak velocities as shown in FIG. 5 (a) and stroke volumes as shown in FIG. 5 (b) is used to compare the proposed real-time imaging method with the conventional cine method. The peak velocities and stroke volumes from real-time imaging are the mean values over 30 consecutive cardiac cycles. As can be seen， the results from the proposed method are in excellent agreement with those from the conventional cine method. In addition， we performed the paired Student’s t-test analysis on the two methods， and the correlation coefficients for peak velocities and stroke volumes are 0.94 (P ＜ 0.0001 ) and 0.90 (P ＝ 0.0002) ， respectively. This further confirms strong correlation between the two methods.

FIG. 6 shows 2D real-time PC-MRI for the 23-year-old arrhythmic patient (with mild cardiac arrhythmia) according to the embodiment of the present disclosure. FIG. 6 (a) shows the ECG recordings and the velocity waveforms of AAo and DAo. FIG. 6 (b) shows the magnitude images and velocity maps for the three representative time frames within an arrhythmic period. As expected， the proposed method is able to reconstruct flow variations over different cardiac cycles. In particular， it nicely captures a sudden flow velocity drop occurring in an arrhythmic period (see FIG. 6 (b) ) . Note that this type of flow dynamics cannot be obtained from the conventional cine method. Further， it is worth noting that the flow velocity variations correlate well with the ECG signal recorded during the acquisition. Besides， we show three snapshot images from the proposed method. Clearly， the velocity maps confirm the dramatic flow variations within the arrhythmic period.

FIG. 7 shows 2D real-time PC-MRI for the 72-year-old arrhythmic patient according to the embodiment of the present disclosure. The velocity waveforms associated with the AAo and DAo from the proposed method are shown in FIG. 7 (a) . Again， the proposed method well captures irregular flow variations， which are more significant than the ones from the previous patient. Moreover， we show the reconstructed magnitude images and velocity maps in FIG. 7 (b) .

Results for 3D real-time PC-MRI

FIG. 8 shows the reconstructed velocity maps from 3D cine flow imaging and the proposed real-time imaging method for a healthy subject. FIG. 8 (a) shows a systolic frame and FIG. 8 (b) shows a diastolic frame. As can be seen， the proposed method provides comparable quality to the cine flowimaging method.

FIG. 9 shows the reconstructed flow waveforms from the proposed 3D real-time flow for a healthy subject. The flow waveforms associated with ascending aorta (AAo) and descending aorta (DAo) marked in FIG. 9 (a) are respectively presented in FIGs. 9 (b) and 9 (c) . As can be seen， the proposed method is able to resolve flow variations. We further averaged the flow waveforms from the proposed 3D real-time imaging method over 10 cardiac cycles into one synthetic flow waveform， and compared it with the one from the cine acquisition. As is shown in FIGs. 9 (d) and 9 (e) ， the averaged flow waveform from the proposed method well correlates with the one obtained from the cine flow imaging method， which is highly desirable.

Moreover， the Bland-Altman analysis of the peak velocities and stroke volumes are shown in FIG. 10 for the healthy subjects. Specifically， the Bland-Altman analysis of peak velocities as shown in FIG. 10 (a) and stroke volumes as shown in FIG. 10 (b) is used to compare the proposed 3D real-time flow imaging with the conventional 3D cine method for five healthy subjects (3 males， mean age： 21-year old) . The peak velocities and stroke volumes from the proposed method are the averaged ones over 10 consecutive cardiac cycles. These results further confirm the excellent agreement between the two methods.

FIG. 11 shows 3D Real-time PC-MRI for an arrhythmic patient. FIGs. 11 (a) and 11 (b) respectively show the reconstructed flow waveforms of AAo and DAo. FIG. 11 (c) shows the 3D pathline visualization for the four representative time frames from a normal period (A and B) and an arrhythmic period (C and D) . As can be seen， the proposed method well resolves the beat-by-beat pathological variations， and specifically， nicely captures a dramatic change of flow occurring during an arrhythmia period. This cannot be obtained from the conventional cine method.

Discussion

In the present embodiment， we introduced a new real-time flow imaging method and systematically demonstrated its effectiveness with in vivo experiments. Here， it is worth reiterating the key characteristics of the proposed method. First， it can be used as a viable alternative to the conventional cine flow imaging method in that it provides comparable (if not superior) image quality and flow information for healthy subjects. Second， it is able to resolve beat-by-beat physiological and/or pathological flow variations， which cannot be obtained from the conventional cine method. Such information is often clinically important (e.g.， for assessing cardiac arrhythmia) .

As with other model-based methods， the proposed method involves model selection (i.e.， selection of the rank L) . Generally， the selection of L needs to balance the model representational power， the number of measurements (i.e.， acquisition time) ， and signal-to-noise ratio [36] . In the present embodiment， we manually selected L to trade off the above factors， and it consistently yielded good reconstruction performance， although principled model selection methods (e.g.， [46， 47] ) can be investigated in future research.

The proposed method for 2D real-time flow imaging is computationally efficient. The algorithm runtime for reconstructing an in vivo dataset (from 94 s real-time acquisition) takes around 10 min on a workstation with 64 GB RAM and 3.47 GHz CPU. However， for 3D real-time flow imaging， the runtime for reconstructing an in vivo dataset (from 20 mins real-time acquisition) would take more than one hour. To enhance the practical utility， computational efficiency may be improved by an implementation on graphical processing units. Such an investigation is worthwhile to explore for future research.

This work is mainly focused on the development of a novel real-time flow imaging technique， which should serve as a foundation for our subsequent clinical studies. Given that the proposed method well resolves beat-by-beat flow variations， it can provide more information on hemodynamics for patients with significant irregular heartbeats. In the future work， we plan to conduct systematic study of the proposed method for various potential clinical applications (e.g.， atrial fibrillation， premature atrial contraction or congenital heart disease) .

Conclusions

The present embodiment presents a new model-based method for high-resolution real-time PC-MRI without ECG gating and respiration control， and for the first time to achieve 3D real-time PC-MRI. It features the novel low-rank model and the integration with parallel imaging， which together enable high-quality reconstruction from highly undersampled (k， t) -space data for real-time PC-MRI. The effectiveness and utilities of the proposed method have been demonstrated for 2D and 3D real-time PC-MRI with in vivo experiments. We expect that the proposed method will enhance the practical utility of real-time PC-MRI for various clinical applications.

In the foregoing specification， specific embodiments have been described. However， one of ordinary skill in the art appreciates that various modifications and changes can be made without departing from the scope of the invention as set forth in the claims below. Accordingly， the specification and figures are to be regarded in an illustrative rather than a restrictive sense， and all such modifications are intended to be included within the scope of present teachings.

The benefits， advantages， solutions to problems， and any element (s) that may cause any benefit， advantage， or solution to occur or become more pronounced are not to be construed as a critical， required， or essential features or elements of any or all the claims. The invention is defined solely by the appended claims including any amendments made during the pendency of this application and all equivalents of those claims as issued.

Moreover in this document， relational terms such as first and second， top and bottom， and the like may be used solely to distinguish one entity or action from another entity or action without necessarily requiring or implying any actual such relationship or order between such entities or actions. The terms ″comprises， ″ ″comprising， ″ “has” ， “having， ” “includes” ， “including， ” “contains” ， “containing” or any other variation thereof， are intended to cover a non-exclusive inclusion， such that a process， method， article， or apparatus that comprises， has， includes， contains a list of elements does not include only those elements but may include other elements not expressly listed or inherent to such process， method， article， or apparatus. An element proceeded by “comprises... a” ， “has... a” ， “includes... a” ， “contains... a” does not， without more constraints， preclude the existence of additional identical elements in the process， method， article， or apparatus that comprises， has， includes， contains the element. The terms “a” and “an” are defined as one or more unless explicitly stated otherwise herein. The terms “substantially” ， “essentially” ， “approximately” ， “about” or any other version thereof， are defined as being close to as understood by one of ordinary skill in the art， and in one non-limiting embodiment the term is defined to be within 10％， in another embodiment within 5％， in another embodiment within 1％and in another embodiment within 0.5％. The term “coupled” as used herein is defined as connected， although not necessarily directly and not necessarily mechanically. A device or structure that is “configured” in a certain way is configured in at least that way， but may also be configured in ways that are not listed.

It will be appreciated that some embodiments may be comprised of one or more generic or specialized processors (or “processing devices” ) such as microprocessors， digital signal processors， customized processors and field programmable gate arrays (FPGAs) and unique stored program instructions (including both software and firmware) that control the one or more processors to implement， in conjunction with certain non-processor circuits， some， most， or all of the functions of the method and/or apparatus described herein. Alternatively， some or all functions could be implemented by a state machine that has no stored program instructions， or in one or more application specific integrated circuits (ASICs) ， in which each function or some combinations of certain of the functions are implemented as custom logic. Of course， a combination of the two approaches could be used.

Moreover， an embodiment can be implemented as a computer-readable storage medium having computer readable code stored thereon for programming a computer (e.g.， comprising a processor) to perform a method as described and claimed herein. Examples of such computer-readable storage mediums include， but are not limited to， a hard disk， a CD-ROM， an optical storage device， a magnetic storage device， a ROM (Read Only Memory) ， a PROM (Programmable Read Only Memory) ， an EPROM (Erasable Programmable Read Only Memory) ， an EEPROM (Electrically Erasable Programmable Read Only Memory) and a Flash memory. Further， it is expected that one of ordinary skill， notwithstanding possibly significant effort and many design choices motivated by， for example， available time， current technology， and economic considerations， when guided by the concepts and principles disclosed herein will be readily capable of generating such software instructions and programs and ICs with minimal experimentation.

The Abstract of the Disclosure is provided to allow the reader to quickly ascertain the nature of the technical disclosure. It is submitted with the understanding that it will not be used to interpret or limit the scope or meaning of the claims. In addition， in the foregoing Detailed Description， it can be seen that various features are grouped together in various embodiments for the purpose of streamlining the disclosure. This method of disclosure is not to be interpreted as reflecting an intention that the claimed embodiments require more features than are expressly recited in each claim. Rather， as the following claims reflect， inventive subject matter lies in less than all features of a single disclosed embodiment. Thus the following claims are hereby incorporated into the Detailed Description， with each claim standing on its own as a separately claimed subject matter.

Funding

This work was partially supported by the National Key R&D Program during the “13th Five-Year Plan” (2016YFC1301601 ) ， and National Institute of Health (NIH-RO1-EB013695) .

References

[1] Van Dijk， P.： Direct cardiac NMR imaging of heart wall and blood flow velocity. J Comput Assist Tomogr 8 (3) ， 429-436 (1984)

[2] Nayler， G.， Firmin， D.， Longmore， D.： Blood flow imaging by cine magnetic resonance. J Comput Assist Tomogr 10 (5) ， 715-722 (1986)

[3] Gatehouse， P.D.， Keegan， J.， Crowe， L.A.， Masood， S.， Mohiaddin， R.H.， Kreitner， K. -F.， Firmin， D.N.： Applications of phase-contrast flow and velocity imaging in cardiovascular MRI. Eur Radiol 15 (10) ， 2172-2184 (2005)

[4] Markl， M.， Kilner， P.J.， Ebbers， T.： Comprehensive 4D velocity mapping of the heart and great vessels by cardiovascular magnetic resonance. J Cardiovasc Magn Reson 13 (7) ， 1-22 (2011)

[5] Nayak， K.S.， Nielsen， J. -F.， Bernstein， M.A.， Markl， M.， D. Gatehouse， P.， M. Botnar， R.， Saloner， D.， Lorenz， C.， Wen， H.， S. Hu， B.， Epstein， F.H.， N. Oshinski， J.， Raman， S.V.： Cardiovascular magnetic resonance phase contrast imaging. J Cardiovasc Magn Reson 17 (1) ， 1 (2015)

[6] Markl， M.， Frydrychowicz， A.， Kozerke， S.， Hope， M.， Wieben， O.： 4D flow MRI. J Magn Reson Imaging 36 (5) ， 1015-1036 (2012)

[7] Hope， T.A.， Hope， M.D.， Purcell， D.D.， von Morze， C.， Vigneron， D.B.， Alley， M.T.， Dillon， W.P.： Evaluation of intracranial stenoses and aneurysms with accelerated 4D flow. Magn Reson Imaging 28 (1) ， 41-46 (2010)

[8] Ansari， S.， Schnell， S.， Carroll， T.， Vakil， P.， Hurley， M.， Wu， C.， Carr， J.， Bendok， B.， Batjer， H.， Markl， M.： Intracranial 4D flow MRI： Toward individualized assessment of arteriovenous malformation hemodynamics and treatment-induced changes. Am J Neuroradiol 34 (10) ， 1922-1928 (2013)

[9] Lenz， G.W.， Haacke， E.M.， White， R.D.： Retrospective cardiac gating： A review of technical aspects and future directions. Magn Reson Imaging 7 (5) ， 445-455 (1989)

[10] Dyverfeldt， P.， Bissell， M.， Barker， A.J.， Bolger， A.F.， Carlh¨all， C. -J.， Ebbers， T.， Francios， C.J.， Frydrychowicz， A.， Geiger， J.， Giese， D.， Hope， M.D.， Kilner， P.J.， Kozerke， S.， Myerson， S.， Neubauer， S.， Wieben， O.， Markl， M.： 4D flow cardiovascular magnetic resonance consensus statement. J Cardiovasc Magn Reson 17 (1) ， 1-19 (2015)

[11] Finn， J.P.， Nael， K.， Deshpande， V.， Ratib， O.， Laub， G.： Cardiac MR imaging： State of the technology. Radiology 241 (2) ， 338-354 (2006)

[12] Thavendiranathan， P.， Verhaert， D.， Walls， M.C.， Bender， J.A.， Rajagopalan， S.， Chung， Y. -C.， Simonetti， O.P.， Raman， S.V.： Simultaneous right and left heart real-time， free-breathing CMR flow quantification identifies constrictive physiology. J Am Coll Cardiol： Cardiovasc Imaging 5 (1) ， 15-24 (2012)

[13] Kowallick， J.， Joseph， A.， Unterberg-Buchwald， C.， Fasshauer， M.， Van Wijk， K.， Merboldt， K.， Voit， D.， Frahm， J.， Lotz， J.， Sohns， J.： Real-time phase-contrast flow MRI of the ascending aorta and superior vena cava as a function of intrathoracic pressure (Valsalva manoeuvre) . Br J Radiol 87 (1042) ， 20140401 (2014)

[14]

H.， Barth， P.， Gieseke， J.， Müller， K.， Burchert， W.， Esdorn， H.， Kececioglu， D.， Beerbaum， P.， Laser， K.T.： Impact of respiration on stroke volumes in paediatric controls and in patients after Fontan procedure assessed by MR real-time phase-velocity mapping. Eur Heart J Cardiovasc Imaging 16 (2) ， 198-209 (2015)

[15] Riederer， S.J.， Wright， R.C.， Ehman， R.L.， Rossman， P.J.， Holsinger-Bampton， A.E.， Hangiandreou， N.J.， Grimm， R.C.： Real-time interactive color flow MR imaging. Radiology 181 (1) ， 33-39 (1991)

[16] Liu， C. -Y.， Varadarajan， P.， Pohost， G.M.， Nayak， K.S.： Real-time color-flow MRI at 3 T using variable-density spiral phase contrast. Magn Reson Imaging 26 (5) ， 661-666 (2008)

[17] Joseph， A.A.， Merboldt， K. -D.， Voit， D.， Zhang， S.， Uecker， M.， Lotz， J.， Frahm， J.： Real-time phase-contrast MRI of cardiovascular blood flow using undersampled radial fast low-angle shot and nonlinear inverse reconstruction. NMR Biomed 25 (7) ， 917-924 (2012)

[18] Eichenberger， A.， Schwitter， J.， McKinnon， G.， Debatin， J.， Von Schulthess， G.： Phase-contrast echo-planar MR imaging： Real-time quantification of flow and velocity patterns in the thoracic vessels induced by Valsalva’s maneuver. J Magn Reson Imaging 5 (6) ， 648-655 (1995)

[19] Shankaranarayanan， A.， Simonetti， O.P.， Laub， G.， Lewin， J.S.， Duerk， J.L.： Segmented k-space and real-time cardiac cine MR imaging with radial tajectories. Radiology 221 (3) ， 827-836 (2001)

[20] Gatehouse， P.D.， Firmin， D.N.， Collins， S.， Longmore， D.B.： Real time blood flow imaging by spiral scan phase velocity mapping. Magn Reson Med 31 (5) ， 504-512 (1994)

[21] Nayak， K.S.， Pauly， J.M.， Kerr， A.B.， Hu， B.S.， Nishimura， D.G.： Real-time color flow MRI. Magn Reson Med 43 (2) ， 251-258 (2000)

[22] Steeden， J.A.， Atkinson， D.， Taylor， A.M.， Muthurangu， V.： Assessing vascular response to exercise using a combination of real-time spiral phase contrast MR and noninvasive blood pressure measurements. J Magn Reson Imaging 31 (4) ， 997-1003 (2010)

[23] Kowalik， G.T.， Steeden， J.A.， Pandya， B.， Odille， F.， Atkinson， D.， Taylor， A.， Muthurangu， V.： Real-time flow with fast GPU reconstruction for continuous assessment of cardiac output. J Magn Reson Imaging 36 (6) ， 1477-1482 (2012)

[24] Pruessmann， K.P.， Weiger， M.， Scheidegger， M.B.， Boesiger， P.： SENSE： Sensitivity encoding for fast MRI. Magn Reson Med 42 (5) ， 952-962 (1999)

[25] Griswold， M.A.， Jakob， P.M.， Heidemann， R.M.， Nittka， M.， Jellus， V.， Wang， J.， Kiefer， B.， Haase， A.： Generalized autocalibrating partially parallel acquisitions (GRAPPA) . Magn Reson Med 47 (6) ， 1202-1210 (2002)

[26] Weiger， M.， Pruessmann， K.P.， Boesiger， P.： Cardiac real-time imaging using SENSE. Magn Reson Med 43 (2) ， 177-184 (2000)

[27] Hoogeveen， R.， Leone， B.， van der Brink， J.： Real-time quantitative flow using EPI and SENSE. In： Proceedings of the 9th Annual Meeting of the International Society for Magnetic Resonance in Medicine： 21-27 April 2001； Glasgow， Scotland， vol. 9， p. 114 (2001)

[28]

H.， Gieseke， J.， Barth， P.， Hoogeveen， R.， Esdorn， H.，

A.， Meyer， H.， Beerbaum， P.： Flow volume and shunt quantification in pediatric congenital heart disease by real-time magnetic resonance velocity mapping. Circulation 109 (16) ， 1987-1993 (2004)

[29] Wintersperger， B.J.， Nikolaou， K.， Dietrich， O.， Rieber， J.， Nittka， M.， Reiser， M.F.， Schoenberg， S.O.： Single breath-hold real-time cine MR imaging： Improved temporal resolution using generalized autocalibrating partially parallel acquisition (GRAPPA) algorithm. Eur Radiol 13 (8) ， 1931-1936 (2003)

[30] Jones， A.， Steeden， J.A.， Pruessner， J.C.， Deanfield， J.E.， Taylor， A.M.， Muthurangu， V.： Detailed assessment of the hemodynamic response to psychosocial stress using real-time MRI. J Magn Reson Imaging 33 (2) ， 448-454 (2011)

[31] Joseph， A.， Kowallick， J.T.， Merboldt， K. -D.， Voit， D.， Schaetz， S.， Zhang， S.， Sohns， J.M.， Lotz， J.， Frahm， J.： Real-time flow MRI of the aorta at a resolution of 40 msec. J Magn Reson Imaging 40 (1) ， 206-213 (2014)

[32] Tan， Z.， Roeloffs， V.， Voit， D.， Joseph， A.A.， Untenberger， M.， Merboldt， K.D.， Frahm， J.： Model-based reconstruction for real-time phase-contrast flow MRI： Improved spatiotemporal accuracy. Magn Reson Med (2016) . doi： 10.1002/mrm. 26192

[33] Uecker， M.， Hohage， T.， Block， K.T.， Frahm， J.： Image reconstruction by regularized nonlinear inversion-joint estimation of coil sensitivities and image content. Magn Reson Med 60 (3) ， 674-682 (2008)

[34] Liang， Z. -P.： Spatiotemporal imaging with partially separable functions. In： Proceedings of IEEE International Symposium on Biomedical Imaging： April 2007； Washington， DC， USA， pp. 988-991 (2007)

[35] Zhao， B.， Haldar， J.P.， Brinegar， C.， Liang， Z. -P.： Low rank matrix recovery for real-time cardiac MRI. In： Proceedings of IEEE International Symposium on Biomedical Imaging： April 2010； Rotterdam， The Netherlands， pp. 996-999 (2010)

[36] Zhao， B.， Haldar， J.P.， Christodoulou， A.G.， Liang， Z. -P.： Image reconstruction from highly undersampled-space data with joint partial separability and sparsity constraints. IEEE Trans Med Imaging 31 (9) ， 1809-1820 (2012)

[37] Christodoulou， A.G.， Zhang， H.， Zhao， B.， Hitchens， T.K.， Ho， C.， Liang， Z. -P.：High-resolution cardiovascular MRI by integrating parallel imaging with low-rank and sparse modeling. IEEE Trans Biomed Eng 60 (11) ， 3083-3092 (2013)

[38] Haldar， J.P.， Liang， Z. -P.： Spatiotemporal imaging with partially separable functions： A matrix recovery approach. In： Proceedings of IEEE International Symposium on Biomedical Imaging： April 2010； Rotterdam， The Netherlands， pp. 716-719 (2010)

[39] Lingala， S.G.， Hu， Y.， DiBella， E.， Jacob， M.： Accelerated dynamic MRI exploiting sparsity and low-rank structure： k-t SLR. IEEE Trans Med Imaging 30 (5) ， 1042-1054 (2011)

[40] Haldar， J.P.， Hernando， D.： Rank-constrained solutions to linear matrix equations using power factorization. IEEE Signal Proc Let 16 (7) ， 584-587 (2009)

[41] Shen， Y.， Wen， Z.， Zhang， Y.： Augmented lagrangian alternating direction method for matrix separation based on low-rank factorization. Optim Method Softw 29 (2) ， 239-263 (2014)

[42] Gupta， A.S.， Liang， Z. -P.： Dynamic imaging by temporal modeling with principal component analysis. In： Proceedings of the 9th Annual Meeting of the International Society for Magnetic Resonance in Medicine： 21-27 April 2001； Glasgow， Scotland， UK， p. 10 (2001)

[43] Zhao， B.， Haldar， J.P.， Christodoulou， A.G.， Liang， Z.P.： Further development of image reconstruction from highly undersampled (k， t) -space data with joint partial separability and sparsity constraints. In： Proceedings of IEEE International Symposium on Biomedical Imaging： April 2011； Chicago， USA， pp. 1593-1596 (2011)

[44] Zhang， T.， Pauly， J.M.， Vasanawala， S.S.， Lustig， M.： Coil compression for accelerated imaging with Cartesian sampling. Magn Reson Med 69 (2) ， 571-582 (2013)

[45] Walker， P.G.， Cranney， G.B.， Scheidegger， M.B.， Waseleski， G.， Pohost， G.M.， Yoganathan， A.P.： Semiautomated method for noise reduction and background phase error correction in MR phase velocity data. J Magn Reson Imaging 3 (3) ， 521-530 (1993)

[46] Stoica， P.， Selen， Y.： Model-order selection： a review of information criterion rules. IEEE Signal Process Mag 21 (4) ， 36-47 (2004)

[47] Haldar， J.P.： Constrained imaging： denoising and sparse sampling. PhD thesis， University of Illinois at Urbana-Champaign， Electrical &Computer Engineering Department (2011)

Claims

A real-time phase-contrast flow magnetic resonance imaging (MRI) method， comprising：

acquiring real-time phase-contrast MRI (PC-MRI) data， which includes training data and imaging data；

performing a low-rank based image reconstruction based on the acquired training data and imaging data；

calculating velocity maps based on the reconstructed real-time flow images； and

performing quantitative flow analysis based on the calculated velocity maps.
The method of claim 1， wherein performing a low-rank based image reconstruction further comprises：

performing a temporal interpolation on the training data； and

performing a temporal subspace estimation on the interpolated training data.
The method of claim 2， wherein performing a low-rank based image reconstruction further comprises：

performing an estimation of coil sensitivities on the acquired imaging data； and

performing a spatial subspace estimation based on the estimated temporal subspace， estimated coil sensitivities and the acquired imaging data.
The method of claim 1， wherein performing a low-rank based image reconstruction further comprises：

representing each flow image sequence as a spatiotemporal Casorati matrix：

where ρ_v (r， t) denotes the dynamic image associated with either the flow-compensated (v ＝ 1 ) or flow-encoded image sequence (v ＝ 2， …， N_v， ) ，

introducing the following joint Casorati matrix：

on which the low-rank structure， rank (C) ≤ L， is enforced， where L refers to the rank of the matrix C， and satisfies L ≤ min (M， N) ；

using an explicit rank constraint via matrix factorization：

C＝ UV，

where
and
and wherein the columns of U and rows of V respectively span the spatial subspace and temporal subspace of C；

performing a spatial subspace estimation for 2D real-time PC-MRI by solving the equation：

where d_i denotes the measured data， Ω the sparse sampling operator， F_s the spatial Fourier transform matrix， S_i the sensitivity map， i ＝ 1， 2， …， N_c for the i th receiver coil；

after solving
reconstructing the joint Casorati matrix as：

from which each flow image sequence is obtained and the flow velocities are estimated， where U /V /C denote the underlying true quantities， whereas
denote the corresponding reconstructions.
The method of claim 1， wherein performing a low-rank based image reconstruction further comprises：

representing each flow image sequence as a spatiotemporal Casorati matrix：

where ρ_v (r， t) denotes the dynamic image associated with either the flow-compensated (v ＝ 1 ) or flow-encoded image sequence (v ＝ 2， …， N_v ) ，

introducing the following joint Casorati matrix：

on which the low-rank structure， rank (C) ≤ L， is enforced， where L refers to the rank of the matrix C， and satisfies L ≤ min (M， N) ；

using an explicit rank constraint via matrix factorization：

C＝ UV，where
and
and wherein the columns of U and rows of V respectively span the spatial subspace and temporal subspace of C；

performing a spatial subspace estimation for 3D real-time PC-MRI by solving the equation：

where d_i denotes the measured data， Ω the sparse sampling operator， F_s and F_i respectively the spatial and temporal Fourier transform matrix， S_i the sensitivity map， i ＝ 1， 2， …， N_c for the i th receiver coil， λ the regularization parameter， and vec (·) the operator concatenating the columns of a matrix into a vector；

after solving
reconstructing the joint Casorati matrix as：

from which each flow image sequence is obtained and the flow velocities are estimated， where U /V /C denote the underlying true quantities， whereas
denote the corresponding reconstructions.