A CT Reconstruction Algorithm Based on L1/2 Regularization

Computed tomography (CT) reconstruction with low radiation dose is a significant research point in current medical CT field. Compressed sensing has shown great potential reconstruct high-quality CT images from few-view or sparse-view data. In this paper, we use the sparser L1/2 regularization operator to replace the traditional L1 regularization and combine the Split Bregman method to reconstruct CT images, which has good unbiasedness and can accelerate iterative convergence. In the reconstruction experiments with simulation and real projection data, we analyze the quality of reconstructed images using different reconstruction methods in different projection angles and iteration numbers. Compared with algebraic reconstruction technique (ART) and total variance (TV) based approaches, the proposed reconstruction algorithm can not only get better images with higher quality from few-view data but also need less iteration numbers.


Introduction
Since computed tomography (CT) technique was born in 1973, CT has been widely applied in medical diagnose, industrial nondestructive detection, and so forth [1]. In medical CT field, CT reconstruction with low radiation dose is a significant research problem, which needs the reconstruction of high-quality CT images from few-view or sparse-view data. Recently, compressed sensing (CS) [2] theory has been applied in CT image reconstruction; it is possible to reconstruct high-quality images from few-view data under the frame of CS. In CS theory, sparse signal is often regularized with L (0 ≤ ≤ 1) regularization; L 0 regularization is the sparest and most ideal regularization norm (L 0 regularization denotes the number of nonzero signal elements). However, L 0 regularization is susceptible to noise interference and it is difficult to solve equations, and L 1 regularization is usually used as the regularization operator. Theoretically, if regularization is closer to L 0 regularization could get higher-quality CT images in CT reconstruction.
Recently, Xu et al. proposed L 1/2 regularization [3] and a soft thresholding algorithm [4], and Xu and Wang proposed a hybrid soft thresholding algorithm [5]. L 1/2 regularization has some theoretical properties, such as unbiasedness, sparse regularization, and Oracle. It is easier to produce sparser solution compared with the current L 1 regular operator used in CT reconstruction. Theoretically, CT reconstruction method based on L 1/2 regularization will reconstruct higher-quality CT images with few-view data.
On solving L 1 regularization problem, Goldstein and Osher proposed Split Bregman method [6] which is derived from Bregman iterative algorithm [7]. Split Bregman method uses an intermediate variable to split L 1 regularization and L 2 regularization into two equations; L 2 regularization equation can be solved by gradient descent method, and L 1 regularization equation can be solved by thresholding algorithm. Split Bregman method can accelerate iterative convergence and produce better results. Based on Split Bregman method, Vandeghinste et al. proposed Split Bregmanbased sparse-view CT reconstruction [8] and iterative CT reconstruction using shearlet-based regularization [9]. Chu et al. proposed multienergy CT reconstruction based on low rank and sparsity with the Split Bregman method (MLRSS) [10]. Chang et al. proposed a few-view reweighted sparsity hunting (FRESH) method for CT image reconstruction [11].
In this paper, we propose a CT reconstruction algorithm based on L 1/2 regularization and Split Bregman method. In the following section, L 1/2 regularization, Split Bregman method, and the proposed algorithm will be introduced. In the third section, we will use the proposed algorithm to analyze numerical phantom and real projection data. In the last section, we will discuss relevant issues and conclude the paper.

Materials and Methods
2.1. L 1/2 Regularization. Generally, CT reconstruction algorithm can be divided into analytic reconstruction algorithm and iterative reconstruction algorithm; the current typical analytic reconstruction algorithm is filter back-projection (FBP), and iterative reconstruction algorithm contains algebraic reconstruction technique (ART) [12], simultaneous algebraic reconstruction technique (SART) [13], expectationmaximization (EM) [14], and so forth. Mathematical model [15] of CT image reconstruction can be expressed as where A = ( ) is the projection matrix, b = ( 1 , 2 , . . . , ) ∈ is the projection data, and u = ( 1 , 2 , . . . , ) ∈ is the reconstruction image. For sparse-view data, it is difficult to reconstruct highquality images using the conventional CT image reconstruction algorithms, especially for analytic reconstruction algorithms which require high completeness of data. Meanwhile, there are also some artifacts in the reconstruction images using the conventional iterative reconstruction algorithms. In 2006, Donoho put forward the compressed sensing (CS) theory [2]; its main idea is that most of signals are sparse in the proper orthogonal transform domain, which means most of signal transformation coefficients are close to zero or equal to zero in orthogonal transformation, such as gradient transformation [16] and shearlet transformation [9]. In CS theory, an image can be reconstructed from a rather limited amount of data as long as an underlying image can be sparsely represented in an appropriate domain and determined from these data. CT sampling signal is a typical sparse signal; the x-ray attenuation coefficients of some regions of tested objects (i.e., human body) are similar or equal. Thus, CT reconstruction approaches based on compressed sensing can reconstruct high-quality CT images from sparse-view data.
CT reconstruction problem can be converted to a constrained optimization problem where ( ) is the regularization function, usually denoted by L 1 norm of wavelet, gradient, and so forth. In order to simplify (2), we can use the penalty function method in the optimization method to convert it into an unconstraint optimization problem where is the weight coefficient, ( ) can be denoted by total variation (TV) which is a research hotspot in current CT research field. TV method has been widely used in sparseview and limited angle CT reconstruction [7].
In compressed sensing theory, L 0 norm is the most ideal regularization norm, but it is difficult to solve equations with L 0 norm, and L 0 regularization is easily interfered by noise in CT reconstruction, so L 0 norm is commonly replaced by L 1 norm. Theoretically, using a regularization norm which is closer to L 0 norm will reconstruct higher-quality CT images. The definitions of L 0 , L 1 , and L 1/2 regularization norm are where = ( 1 , 2 , . . . , ) and denotes the nonzero elements number in matrix .
As shown in Figure 1, L 1/2 regularization norm is closer to L 0 regularization norm than L 1 regularization norm.
Xu et al. proposed a L 1/2 regularization fast algorithm [4]; it can be expressed as where is the regularization coefficient.
If is the optimal solution of (5), then the optimal condition of (5) will be denoted by where ∇(‖ ‖ 1/2 1/2 ) is the gradient of ‖ ‖ 1/2 1/2 ; multiplying coefficient and adding at the both sides of (6), then we have The definition of operator is Then the optimal solution can be represented as where the operator is Please see [4] for more proof details.

Split Bregman Method.
In order to solve (3), Goldstein and Osher proposed Split Bregman method [6], using an intermediate variable to split L 1 regularization and L 2 regularization into two equations; L 2 regularization equation can be solved by gradient descent method and L 1 regularization equation can be solved by thresholding algorithm. Then the unconstrained optimal problem of (3) can be converted into where Φ is the sparse transform and the common used sparse transform contains gradient, wavelet, shearlet, and so forth.
Using an intermediate variable = Φ( ), (11) can be converted into where is the coefficient. Then (12) can be converted into two unconstrained optimal Bregman problems; it can be expressed as Equation (13) can split into two equations: There are several advantages of Split Bregman method. Firstly, Split Bregman method can accelerate iterative convergence and calculate better results. Secondly, Split Bregman method can be widely used in CT reconstruction; it can not only solve L 1 regularization problem but also solve other regularization problems.

CT Reconstruction Algorithm
Based on L 1/2 Regularization. According to aforementioned methods, we propose a CT reconstruction algorithm based on L 1/2 regularization, where L 1/2 norm is used as the regularization norm and gradient as the sparse conversion; then (11) can be expressed as Combine with Split Bregman method to solve (17) as follows.
Step 1. One has Step 2. One has Step 3. One has Step 2 can be solved by the method in Section 2.1 and Step 3 can be solved directly. To solve Step 1, we use gradient descent method Equation (21) is derivable and derivation of L 2 regularization as follows: To derivate (21), we have where represents iteration numbers and parameter can be acquired by the following equation:

Numerical Simulation.
In this section, we study the ART algorithm, TV based ART algorithm (ART-TV), and L 1/2 regularization based Split Bregman method (SpBr-L 1/2 ) and analyze the reconstructed images. In this paper, we test Shepp-Logan phantom as shown in Figure 2, and the size of phantom image is 256 × 256. We assume that the CT system was viewed as a typical parallel-beam geometry, and the scanning range was from 0 ∘ to 180 ∘ with a angular increment; projection angles can be indicated as where view is the number of projection angles. We will compare the reconstruction results from noisefree and noise data and projection numbers view = 60. In the simulation, we add 0.01% Gaussian noise to noise-free projection data, and iteration number for every reconstruction algorithm is 50; the reconstruction results are shown in Figure 3.
From Figures 3 and 4, we can see that the reconstructed images using ART and ART-TV methods contain a lot of noise and artifacts from noise-free and noise data, while the where and * are the reconstructed image and original image, respectively, and the image size is × . As shown in Table 1, the RMSE of reconstructed image using SpBr-L 1/2 method is much smaller than reconstructed images using ART and ART-TV methods, which means SpBr-L 1/2 can reconstruct higher-quality images. As shown in Figure 5, the SpBr-L 1/2 algorithm can reconstruct high-quality images at less iteration numbers, which indicates that SpBr-L 1/2 algorithm can accelerate iterative convergence. From Figure 6, we can see that SpBr-L 1/2 algorithm can reconstruct high-quality images at less projection numbers.

Real Data Study.
In this section, we reconstruct oral images using three different algorithms with real projection data, where projection numbers are 90 and iteration numbers are 50 while the original projection numbers are 360. And as shown in Figure 7(a), we evaluate three images with reconstructed image using ART with original projection data as a standard image. The reconstructed images and RMSE of three reconstructed images are shown in Figure 7 and Table 2. As shown in Figure 7, the reconstructed images using ART and ART-TV method have more noise and artifacts, while the reconstructed image using SpBr-L 1/2 method has less noise and artifacts. And from Table 2, the RMSE of reconstructed image using SpBr-L 1/2 method is smaller than that of reconstructed images using ART and ART-TV methods, which means SpBr-L 1/2 can reconstruct higherquality images with clearer details. From the comparison of  reconstructed images using different algorithms, we can see that SpBr-L 1/2 can reconstruct high-quality images at less projection angles and less iteration numbers.

Discussions and Conclusion
There are several issues worth further discussion in the reconstruction study. Firstly, the thresholding algorithm was not applied to solve the L 1/2 regularization problem. There are two reasons. First, if the thresholding algorithm is applied to solve the L 1 regularization problem, the edge and details of reconstructed images will not be clear, while the SpBr-L 1/2 without thresholding algorithm reconstructs images with clear edges. Then, if the threshold value is settled to solve L 1/2 regularization problem, the reconstruction speed will be reduced a lot. Secondly, there are some artifacts in the reconstructed image using SpBr-L 1/2 method from noise data. The reason is that L 1/2 regularization is closer to L 0 regularization, and the ability of denoising of L 0 regularization is bad; therefore the ability of denoising of L 1/2 regularization is not good enough. Thirdly, in the reconstruction from real projection data, the reconstructed image using ART-TV method lost some details, while details of reconstructed image using SpBr-L 1/2 method are much clearer.
In the further research, we will try to use SpBr-L 1/2 algorithm in interior CT and study the region of interest (ROI) reconstruction, which will reduce radiation dose as much as possible.
In conclusion, we proposed a CT reconstruction algorithm based on L 1/2 regularization; the reconstructed results