Rank-One and Transformed Sparse Decomposition for Dynamic Cardiac MRI

It is challenging and inspiring for us to achieve high spatiotemporal resolutions in dynamic cardiac magnetic resonance imaging (MRI). In this paper, we introduce two novel models and algorithms to reconstruct dynamic cardiac MRI data from under-sampled k − t space data. In contrast to classical low-rank and sparse model, we use rank-one and transformed sparse model to exploit the correlations in the dataset. In addition, we propose projected alternative direction method (PADM) and alternative hard thresholding method (AHTM) to solve our proposed models. Numerical experiments of cardiac perfusion and cardiac cine MRI data demonstrate improvement in performance.


Introduction
Dynamic magnetic resonance imaging (MRI) reconstructs a temporal series of images to resolve the motion or the variation of the imaged object. It is often used in cardiac, perfusion, functional, and gastrointestinal imaging. The dynamic contrast variations in the myocardium are typically imaged 40-60 s in cardiac perfusion imaging; however, most patients cannot maintain a breath-hold for such long durations. The data in MRI, samples in − space of the spatial Fourier transform of the object, are acquired sequentially in time. Therefore, achieving high spatiotemporal resolutions is challenging and inspiring in dynamic MRI due to the hardware limitations and the risk of peripheral nerve stimulation.
Compressed sensing (CS) [1,2] has been successfully applied to accelerate MRI, for example, [3][4][5][6][7] and references therein. Lustig et al. [3] systematically exploited the sparsity which is implicit in MRI and developed a framework for using CS in MRI. Liang [4] utilized low-rank matrix completion to dynamic MRI by considering each temporal frame as a column of a space-time matrix. Lingala et al. [5] and Zhao et al. [6] combined CS and low-rank matrix completion which finds a solution that is both low-rank and sparse. Gao et al. [8] proposed a different method which decomposes dynamic MRI matrix as a superposition of a low-rank component and a sparse component. It is worth noting that the low-rank component can model the temporally correlated background and the sparse component can model the dynamic information that lies on the top of the background. Otazo et al. [7] presented low-rank and sparse reconstructions for dynamic MRI using joint multicoil reconstruction for Cartesian and non-Cartesian -space sampling.
Most of the above work is concerned about convex relaxation; then it comes to an interesting question whether we can deal with dynamic MRI using ℓ 0 -minimization. Many authors have made great effort; see [9][10][11][12][13], for example. Blumensath and Davies [9] were the first to propose the iterative hard thresholding (IHT) to solve a type of ℓ 0regularized problems and showed that the IHT method converges to a local minimizer. Beck and Eldar [10] introduced and analyzed three kinds of necessary optimality conditions for ℓ 0 -constrained problems. Bahmani et al. [14] studied ℓ 0constrained optimization in cases where nonlinear models are involved or the cost function is not quadratic. Recently, Lu and Zhang [11] presented a penalty decomposition method for solving a more general class of ℓ 0 -constrained and ℓ 0regularized minimization problems. Lu [12] considered ℓ 0regularized convex cone problems and then gave an IHT method and its convergence. More recently, Pan et al. [13] provided the first-and second-order necessary and sufficient optimality conditions for ℓ 0 -constrained optimization.
Motivated by the recent advances in the use of CS in MRI and ℓ 0 -minimization, we propose rank-one and transformed sparse models and algorithms to significantly accelerate dynamic MRI. In contrast to classical low-rank and transformed sparse decomposition, this paper has the following merits: (i) We introduce the use of rank-one matrix in dynamic cardiac MRI data which avoids singular value decomposition (SVD).
(ii) We assume the dynamic component has a sparse representation in some known sparsifying transformation, not itself.
(iii) We develop new ℓ 0 -constrained model and ℓ 0regularized model, rather than their convex relaxation.
(iv) We design two ADM-based algorithms, projected alternative direction method (PADM) and alternative hard thresholding method (AHTM). Numerical experiments of dynamic cardiac MRI illustrate the efficiency of our proposed algorithms.
The remainder of this paper is organised as follows. Section 2 reviews some existing related work and proposes our new models. Section 3 describes PADM and AHTM to solve our proposed models. Section 4 reports experimental results. Section 5 concludes this paper with some final remarks.

Model Analysis
We begin this section by introducing the low-rank and transformed sparse decomposition in [7]; that is, where ‖ ⋅ ‖ * is defined as the sum of all singular values and ‖ ⋅ ‖ 1 is defined as the sum of absolute values of all entries. ∈ R is a real positive weighting parameter to balance the weights of rank and sparsity, and is the undersampled − data. ∈ C × is a low-rank matrix, which can model the temporally correlated background. ∈ C × is a sparse matrix, which can model the dynamic information that lies on top of the background. Generally, the dynamic component has a sparse representation in some known sparsifying transformation (e.g., temporal frequency domain), hence the idea of minimizing and not itself. For acquisition with multiple receiver coils, is given by the frame-by-frame multicoil encoding operator, which performs a multiplication by coil sensitivities followed by a Fourier transform according to the sampling pattern.
Compared with a low-rank or sparse constraint, (1) can significantly increase compressibility of dynamic MRI data and thus enable high acceleration factors. However, (1) is the convex relaxation problem. The original optimization problem associated with (1) is formulated as follows: and its low-rank and sparse constrained optimization problem is as follows: where denotes the rank of , denotes the sparsity of , and ‖ ⋅ ‖ 0 is ℓ 0 -norm, which counts the number of nonzero entries.
Note that if and reduce to identity matrix, (3) becomes a low-rank and sparse constrained optimization problem: where = * is given data. Zhou and Tao [15] developed GoDec method to estimate the low-rank part and the sparse part in (4). GoDec alternatively updated by SVD and by hard thresholding. To overcome the difficulty of SVD, they proposed bilateral random projections. Even so, it still costs much time.
We find sometimes that the background is exactly a rank-one matrix, instead of a general low-rank matrix, so Li et al. [16] proposed to replace low-rank matrix with rank-one matrix which avoids any SVD completely. Xiu and Kong [17] extended it to the case of tensor decomposition and showed that it performs well in surveillance video. Specifically speaking, we set = 1 , where ∈ C and 1 denotes the vector in R whose entries are all 1, which leads to the following rank-one and transformed sparse matrix decomposition problem: and its unconstrained version where > 0 is a regularized parameter.

Algorithm and Convergence
In this section, we present the projected alternative direction method (PADM) for (5) and alternative hard thresholding method (AHTM) for (6). Then we discuss their convergent properties.
BioMed Research International 3 (5). Let us look at the optimization problem (5). Inspired by ADM-based methods, it is easy to divide (5) into three subproblems:

Projected Alternative Direction Method for
where is obtained by enforcing data consistency, where the aliasing artifacts corresponding to the residual in − space data * ( ( 1 + )− ) are subtracted from 1 + . Here, we denote these subproblems in (7) as -subproblem, -subproblem, and -subproblem, respectively. The first -subproblem has the solution where, for any ∈ C × , mean(⋅) ∈ C is defined as Before establishing the solution of the second -subproblem, we first give a proposition about sparse projection.
where B = { : ‖ ‖ 0 ≤ } = max( ) and the third equation is satisfied by introducing = . Then, it holds that The proof is completed.
From Proposition 1, we derive the solution ofsubproblem: Hence the scheme is summarized as the following projected alternative direction method (PADM). In Algorithm 1, -subproblem is obtained by enforcing data consistency, so the convergence properties can be analyzed by considering the iterations of -subproblem andsubproblem. Now we will present the convergence result of PADM as follows.

Alternative Hard Thresholding Method for
Compared with the above subsection, we find that thesubproblem and -subproblem are the same; hereafter we mainly explain how to solve the -subproblem. By applying the separability of the objective function and the technique of operator splitting, the -subproblem can be converted into corresponding single variable minimization problems. Therefore, the following proposition will give the solution of the corresponding single variable minimization problem.
Thus, we immediately obtain the proposition.
It is straightforward to establish the following result.
where ∈ R × is some transform matrix and is a given matrix in R × . Then, the solution is described as follows: where H 0.5 is the hard thresholding operator defined as Proof. The conclusion holds directly from Proposition 3 by separating the objective function.
From Proposition 4, the solution of -subproblem is characterized as Based on the above argument, we have derived the following alternative hard thresholding method (AHTM) for (6).
In the end of this subsection, we will establish the convergence result of Algorithm 2.

Numerical Experiments
In this section, we conduct numerical experiments to compare the performance of IST [7], PADM, and AHTM for dynamic cardiac MRI data. We apply all these methods on MR images: cardiac perfusion and cardiac cine, which are available at the website http://cai2r.net/resources/ software/ls-reconstruction-matlab-code. All experiments are implemented in MATLAB (MathWorks, Natick, MA) on a desktop computer with Intel Core I5 2.60 GHz CPU and 8 GB of RAM. We now address the initialization and the termination criteria for these methods. In particular, we set the stopping criterion as RelErr < 2.5 − 3, where We choose = 0.01 in AHTM as suggested in [7]. However, we do not know the sparse parameter beforehand. Hence, in this paper, we initialize as one percent of pixels in each frame of .

Cardiac Perfusion.
We next conduct numerical experiments to test the performance of IST [7], PADM, and AHTM for dynamic cardiac perfusion MRI data. The computational results of the three methods are presented in Figures 1-3. In detail, we use + (the first row) to denote the original dynamic cardiac perfusion MRI frames, (the second row) to denote the reconstructed background parts, and (the third row) to denote the reconstructed moving parts. Comparing Figures 1-3, we can observe that PADM (see Figure 2) and AHTM (see Figure 3) are easier to recognize the pathological parts than IST (see Figure 1), which is due to the use of ℓ 0 rather than its convex relaxation ℓ 1 . This capacity may be useful to identify lesions that are difficult to visualize in the original images. In addition, Table 1 reports the computational results for dynamic cardiac perfusion MRI data. Here, the number of iterations is given in the second column (i.e., Iteration). The third column (i.e., Time) lists the computing time in seconds. The relative error is reported in the last column (i.e., RelErr).  It is evident to see that the "Time" and "Iteration" of PADM and AHTM are lower than these of IST. This is because we employ rank-one to avoid SVD. Table 2 demonstrate the performance of IST [7], PADM, and AHTM for dynamic cardiac cine MRI data. From the third row, we find that Figures 5 and 6 can get clearer moving parts, because we can tune the sparsity of and .

Conclusions
In this paper, we focused on the application of rank-one and transformed sparse decomposition for dynamic cardiac MRI.  We established the projected alternative direction method and alternative hard thresholding method and showed their convergent results. Finally, numerical experiments demonstrated the efficiency and effectiveness of our proposed methods.

Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.