Abstract
The quality of dynamic magnetic resonance imaging reconstruction has heavy impact on clinical diagnosis. In this paper, we propose a new reconstructive algorithm based on the model. In the algorithm, the norm is substituted by the norm to approximate the norm; thus the accuracy of the solution is improved. We apply an alternate iteration method to solve the resulting problem of the proposed method. Experiments on nine data sets show that the proposed algorithm can effectively reconstruct dynamic magnetic resonance images.
1. Introduction
Image reconstruction is widely applied in the medical field, and most clinical diagnoses depend on computer hardware equipment; thus the improvement on image reconstruction algorithm has great significance. At present, among all kinds of detection methods, magnetic resonance imaging (MRI) and computed tomography (CT) are the most common and important ways for clinical diagnosis. When checking organs by MR, the final image is affected in varying degrees by artifacts, which results in degrading in quality and affects the diagnoses. Therefore, reconstruction algorithms with high quality and fast calculation have become one of the research focuses in this field.
In recent years, in order to ensure image quality and to speed up the pace of reconstruction, there emerge many improved methods, such as multicoil parallel imaging [1], keyhole imaging technology [2, 3], unfold method [4], and SENSE [5]. After the theory of compressed sensing (CS) [6, 7] was proposed, many scholars applied it to the dynamic MR reconstruction [5, 8, 9]. These methods process the original data by downsampling method instead of traditional full-sampling and reduce the computing time.
Since the development of robust principal component analysis (RPCA) theory [10], breakthroughs have emerged in the field of MR image reconstruction. This theory suggests that images with spatiotemporal correlation can be decomposed into low-rank (LR) matrices and sparse matrices. Many researchers introduced this theory to reconstruct dynamic MR images to achieve better results. In [11], a method was proposed based on accelerated diffusion weighted sequences for image reconstruction by the decomposition. In [12], 3D cardiac MRI was reconstructed by combining decomposition with prior knowledge. In [13], a new model was proposed, which decomposed a series of dynamic MR images with temporal and spatial correlation into LR matrices and sparse matrices. As a result, the dynamic MRI could be divided into the foreground and background to help diagnosis. In order to reconstruct dynamic MR images with better quality, based on the traditional model, a rank-1 and sparse model was proposed in [14].
Recently, a lot of researchers focused on the problem of norm approximation [15–17]. In particular, in [15], Candés et al. pointed out that a nonconvex problem can be transformed into a convex optimization problem using the norm instead of norm. However, the norm and norm are still different in numerical values. For example, Sun et al. suggested that if the nuclear norm of the matrix is regarded as a norm and the rank of the matrix is regarded as a norm, the error between them can not be ignored in [18]. The algorithm in [19] achieved good approaching effect by introducing a nonconvex regularization which was like a norm. In [20], Quach et al. proposed a nonconvex online RPCA model which enhanced the sparsity by minimization of the norm and achieved good results in many fields, such as face modeling and online background removal.
Inspired by the above discussion, combining RPCA theory with a nonconvex idea, we propose a new MR image reconstruction model. The contributions of this paper are as follows: in order to better approximate norm, we introduce a nonconvex regularization using the norm instead of the norm. The alternating direction method of multipliers (ADMM) is used to solve the nonconvex problem. In addition, we point out that the solution is convergent. The framework of this proposed method is shown in Figure 1.
The rest of this paper is organized as follows: in Section 2, a review of the model is presented. In Section 3, we firstly introduce the proposed model and then provide the detailed implementation of the proposed model. Finally, we illustrate the convergence of the proposed algorithm. In Section 4, in order to verify the effectiveness of our algorithm, comparisons between our method and the state-of-the-art algorithms are provided. In Section 5, we conclude our works.
2. Related Works
According to the theory of RPCA [10], we focus on the model where is the observed data, and are the LR and the sparse matrix, respectively, , and denotes the regularization parameter.
Since (1) is a nonconvex optimization problem, the effective solution of the problem is not easy to obtain [21]. According to the theory of RPCA, Otazo et al. [13] used the nuclear norm of and the norm of instead of the rank of and the norm of in (1), respectively, and proposed a convex optimization problem ( model): where is the nuclear norm. denotes norm, which is the sum of absolute values of the elements in matrix .
Let the data of an image correspond to the data of the space by introducing an encoding operator . And assume to be the sparse transformation of ; that is, the dynamic component can be expressed by the known sparse basis . So the minimization problem (2) can be written as
Singular Value Thresholding (SVT) algorithm [15] is introduced to solve the minimization problem (3).
3. The Proposed Model
3.1. Modeling
Although the above methods have good theoretical guarantee and have made major breakthroughs in dynamic MR image reconstruction, norm can not well approximate norm [15, 16, 20]. In addition, authors of [20, 22] pointed out that, among common nonconvex functions, such as -norm [23], SCAD [24], log [15], MCP [25], ETP [26], and Geman [27], -norm is the best approximation to norm. Inspired by this, we propose a new model as follows: where is a LR matrix, is a sparse matrix, and denotes norm. is an undersampled observed data in the space and denotes the ideal reconstructed MR image, where denotes the size of the image. is an encoding operator, which makes the data of an image correspond to the data of the space.
Defining the norm as , the minimization problem (4) can be written as where denotes the elements of located in th row and th column, , and . And assume that the nonconvex function is continuous and monotonically increasing on .
Combining the properties of supergradient for nonsmooth points with Taylor expansion [20], the second term of (5) can be approximated as where is the number of iterations.
Therefore, the new model can be further expressed as
3.2. Solving the Model
By introducing a Lagrangian multiplier, (7) turns to be an unconstrained minimization problem as where and are regularization parameters which balance the nuclear norm and the norm. denotes a Lagrangian multiplier, and is a penalty parameter. In order to solve (8), we divide it into two subproblems as When is obtained, we solve (9) by the Singular Value Thresholding (SVT) algorithm [15], and the solution is where , , and is a positive constant.
When is computed, the minimization problem (10) is equivalent to where and is a small positive real number to avoid the case that the gradient is 0. The minimization problem (12) is solved by the soft-thresholding algorithm [20, 28] as follows: where is a soft-thresholding operator which is defined as and is a threshold.
In the iteration process, after and are computed, we update by where represents artifacts, which is subtracted from the updating as a residual, resulting in a clearer image [22].
Finally, the parameters and are updated as follows: where . The complete process of this algorithm is shown in Algorithm 1.
According to [20], the above solving process is convergent.
Theorem 1 (see [20]). Let be the sequences generated by Algorithm 1, and is monotonically decreasing; that is, . Then the sequences are bounded and have at least one gathering point.
4. Numerical Experimental Comparisons and Analysis
In order to evaluate the performance of the proposed model, we compare our algorithm with the SENSE [5], XD-GRASP [9], and the LplusS [13]. All experiments are run on a PC with Intel Core i5-4690 processor, 8 GB of memory, Win 7 64-bit operating system, and MATLAB 2014a. The quality of the reconstructed image is evaluated by the root mean square error (RMSE), which is defined as where is the observation data and and are the LR matrix and sparse matrix for the reconstructed image, respectively. The smaller the RMSE is, the better the MR image is reconstructed.
The iteration stop indicator for all algorithms of this paper is linked to the relative error (Err), which is defined as The iteration stops when the Err reduction of the consecutive iterations is less than a predefined threshold.
4.1. Experimental Data and the Selection of Parameters
There are nine sets of test data in this paper, three of them are open-data which come from [13], and the others are from a research institution, which are shown in Figure 3. The gray value of all the images ranges from 0 to 255. The values of , , and are the choices of the best reconstruction after repeated tests. For all data in this paper, we set , , and .
One of the three sets of open-data is the dynamic cardiac perfusion data, which consists of 40 images with size . The second is the dynamic cardiac cine data which is made up of 25 images with size . And the third is dynamic abdominal data which is composed of 28 images with size . The parameter combinations for data set one to three are , , and , respectively. When the Err is less than , we stop the iteration for the dynamic cardiac perfusion and cardiac cine images. And, for dynamic abdomen images, the iteration stopping threshold is . The other six data sets from the research institution are all composed of 15 images with size for each set. For image 1, we set and , for image 2 we set and , and, for the remaining four images, we let and , and the Err is set as .
For norm , the choice of heavily influences the convergence rate. Figure 2 shows the convergence rate with , 0.4, and 0.8. The greater results in the faster convergence rate. A large number of experiments also prove that the smaller leads to the smaller Err and the more precise solution. Through repeated trials, we set in all experiments for our algorithm.
(a)
(b)
(c)
(d)
(e)
(f)
4.2. Comparisons and Analysis
In this section, we compare our method with the SENSE, XD-GRASP, and LplusS algorithm in subjective and objective aspects. Figures 4–9 present the visual comparisons, and the differences between the reconstructions by four algorithms are depicted with white arrows. The comparisons on RMSE, Err, and time consumption are listed in Tables 1 and 2, where the boldfaced numbers are the best results of all.
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
(a)
(b)
(c)
(d)
Figures 4 and 5 are the reconstructed comparisons of the first frame and the fifth frame of dynamic cardiac perfusion, respectively. Figure 4(d) is the best reconstructed image of Figure 4. In Figures 4(a) and 4(b), the edge contours of organs at the top right corner are blurred, and the details pointed out by the white arrows are not clear. Only Figure 4(c) is comparable to Figure 4(d), but, at the places of the bottom left and upper right arrows, local details contain more noise, which causes the characteristics of organ unobtrusiveness and is not conducive to disease diagnosis. Figure 5(d) is obviously clearer than Figures 5(a) and 5(b) and superior to Figure 5(c) in the place that the white arrows indicate.
Figures 6 and 7 represent the reconstructions about the dynamic abdomen image with the second frame and the fifth frame, respectively. By comparing, the proposed algorithm is better than the other three algorithms in brightness and resolution.
Figure 8 shows the reconstructed images of the 20th frame of the cardiac cine. The reconstruction effects of all four algorithms have no obvious visual differences, but, through careful observation on the enlarged regions to which the white arrows point, the proposed algorithm performs better.
Figure 9 shows the reconstructed effects of images 1, 2, 4, and 6, respectively. In images 1 and 2, we find that the features of our algorithm are the most obvious among four compared methods, especially in the regions where the white arrows point. In images 4 and 6, the outline and details of the organ of the proposed algorithm are clearer than those of the other reference algorithms, and the result is brighter.
Tables 1 and 2 list the numerical results of the four algorithms. The proposed algorithm ranks the first on the Err and RMSE among the four comparison methods, which denotes that the reconstruction quality is the best. Of all experiments, the time consumption of our algorithm is the second best except for cardiac cine in Table 1, on which the proposed algorithm ranks the first.
5. Conclusions
In this paper, we propose a nonconvex model for reconstructing high-quality dynamic MR images. In the new model, based on the RPCA theory, the norm is substituted by the norm to approximate the norm; thus the accuracy of the solution is improved. An alternate iteration method is applied to solve the proposed model. Experimental results show that the proposed algorithm has advantages in terms of visual effect, RSME, and Err over three reference algorithms. But, in terms of time consumption, the proposed algorithm is slightly inferior to LplusS, which requires further improvement.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This project is partially supported by the National Natural Science Foundation of China (61362021, 61661017, 61272216, and 61572147), Guangxi Natural Science Foundation (2013GXNSFDA019030, 2014GXNSFAA118003, and 2016GXNSFAA380043), Guangxi Colleges and Universities Key Laboratory of Intelligent Processing of Computer Images and Graphics (GIIP201408, GIIP201503), Guangxi Colleges and Universities Key Laboratory of Data Analysis and Computation (LDAC201704), Basic Research Capacity Promotion Project for Youth and Middle-Aged Teachers in Colleges and Universities of Guangxi (ky2016YB162), and Project of Education Innovation Project of Guilin University of Electronic Science and Technology (2017YJCX84).