Total variation with modified group sparsity for CT reconstruction under low SNR

Abstract

BACKGROUND AND OBJECTIVE:

Since the stair artifacts may affect non-destructive testing (NDT) and diagnosis in the later stage, an applicable model is desperately needed, which can deal with the stair artifacts and preserve the edges. However, the classical total variation (TV) algorithm only considers the sparsity of the gradient transformed image. The objective of this study is to introduce and test a new method based on group sparsity to address the low signal-to-noise ratio (SNR) problem.

METHODS:

This study proposes a weighted total variation with overlapping group sparsity model. This model combines the Gaussian kernel and overlapping group sparsity into TV model denoted as GOGS-TV, which considers the structure sparsity of the image to be reconstructed to deal with the stair artifacts. On one hand, TV is the accepted commercial algorithm, and it can work well in many situations. On the other hand, the Gaussian kernel can associate the points around each pixel. Quantitative assessments are implemented to verify this merit.

RESULTS:

Numerical simulations are performed to validate the presented method, compared with the classical simultaneous algebraic reconstruction technique (SART) and the state-of-the-art TV algorithm. It confirms the significantly improved SNR of the reconstruction images both in suppressing the noise and preserving the edges using new GOGS-TV model.

CONCLUSIONS:

The proposed GOGS-TV model demonstrates its advantages to reduce stair artifacts especially in low SNR reconstruction because this new model considers both the sparsity of the gradient image and the structured sparsity. Meanwhile, the Gaussian kernel is utilized as a weighted factor that can be adapted to the global distribution.

Keywords

Computed tomography (CT)image reconstruction total variation (TV)group sparsity signal-to-noise ratio (SNR).

1 Introduction

Computed tomography (CT) is referred to as an important version of the problem which is that of obtaining the attenuated coefficients within the object from multiple X-ray projection data [1]. Owing to the emerging of CT, it has been extensively utilized in medicine and industry [2, 3]. For instance, for new medical application, the dedicated CT systems have been developed in breast imaging [4]. Under the circumstance, radiation exposure to the chest is a key consideration [5], since the chest issues are sensitive to the radiation. In view of the X-ray dose, lowering X-ray tube current intensity is an effective way to lower the exposure dose [6], but it is unfortunate that the obtained projection data will be of low Signal-to-Noise ratio (SNR). Moreover, in industry, reducing exposure time for the fast scanning speed may yet lead to the low SNR projection data. Reconstructing high-quality image from projection data with low SNR is a challenging and ill-posed problem [7].

Several efforts have been made to improve the quality and preserve edges of the image reconstructed from noisy projection data. The corresponding iterative algorithms has been reported. For instance, the classical simultaneous algebraic reconstruction technique (SART) as a major refinement of the algebraic reconstruction technique (ART) [8], has the more potential to acquire decent image quality using projections with low-level noise, compared to the commercial analytical algorithm such as filtered back-projection [9]. However, SART will not work well in condition of the projection data with the high-level noise as shown in Fig. 1 (b).

Fig. 1

(a) is the reference image. (b) represents the image reconstructed from 512 projection views with added Gaussian noise whose average value is zero and standard deviation is 2% maximum projection in current view by using SART.

It is noteworthy that the regularization is an effective technique for addressing severely ill-posed problem using a certain prior information of the image to be reconstructed [10]. Total variation (TV) regularization as a representative is state-of-the-art in handling low-dose and few-view CT [11, 12], and the other sparsity regularization such as wavelet frames [13], Haar transform [14], S-transform [15] and so on. Currently, the noise suppression-guided image filtering reconstruction algorithm was proposed to deal with the low SNR problem [16].

TV based regularization is already commercialized in reconstructing image from noisy and incomplete projection data [11, 17]. But the assumption of TV regularization is that the image to be reconstructed consists of the areas, which are piece-wise constant [18]. TV performs well for the projection data with relatively low-level noise. Whereas, when the SNR of the obtained projection data is low, TV-based methods sometimes yield the undesirable stair artifacts and results in stair artifact [18]. Figure 2 (b) exhibits the performance of TV regularization for disposing noisy projection data. It can be observed that although TV algorithm has ability to deal with the high-level noise, it leads to the stair artifacts and distorts some details.

Fig. 2

Most of the regularization techniques usually put emphasis on the global characteristics of the transformed image, and its constraint according to the image is posed on the whole. Latterly, the regional structure information of the reconstructed image is getting more and more attention. The way of group sparsity (GS) captures this key point. Soon afterwards, a great deal of the researches on GS has been reported in concern with denoising and preserving the edges. The approach that introduces group sparsity into TV was proposed for signal denoising [19]. TV with overlapping group sparsity (OGS) was presented for image restoration to alleviate the stair artifacts [20]. Subsequently, TV with OGS was proposed for deblurring Poisson noisy images, and it has ability to avoid stair effect and preserve edges in the restored images [21]. Recently, a novel image reconstruction for electrical impedance tomography was proposed, which utilizes enhanced adaptive group sparsity with TV constraint [22]. A new approach based on enhanced group sparsity and non-convex regularization was presented to perform compressed sensing MRI (magnetic resonance imaging) [23]. More recently, the concept of group sparsity is widely used in many applications, and they all worked very well.

In this paper, motivated by the above works, a novel model is proposed that introduces the Gaussian kernel and overlapping group sparsity into TV model for CT image reconstruction from the projection data with low SNR. On one hand, this model extends total variation that aims to include group sparsity characteristics of the derivatives of the reconstructed image. It supposes that the derivatives of the reconstructed image are not only sparse, but also exhibit a simple form of structured sparsity. On the other hand, this model utilizes the Gaussian kernel to differentiate the importance of each point in every group. The Gaussian kernel is a kind of smooth functions, which is easy to deal with and implement numerically. The closer to the center of the group, the greater the weight value is. Using this kernel, not only the case of a single point is considered, but also the local characteristics are considered. Numerical simulated experiments showed that GOGS-TV has more ability to alleviate the stair artifacts and preserve the edges of the reconstructed image, compared with SART and TV.

The rest of the paper is organized as follows. In section 2, it will briefly review the fundamental model of image reconstruction, and then propose the novel model for low SNR CT. Consequently, the corresponding algorithm is given to address the presented model. In section 3, it demonstrates several numerical simulated experiments for image reconstruction and verifies the effectiveness of our algorithm, as compared to the classical FBP, SART and the state-of-the-art TV algorithm. Finally, conclusions are made in section 4.

2 Methods

2.1 Model

Generally, from a mathematical point of view, the image reconstruction problem can be formulated as: $Ax + e = b$ (1)

where x denotes the discrete attenuation coefficients of the object to be reconstructed; A is the system matrix; e represents the noise; b is the measurement of the calibrated and log-transformed projection data. Generally, the problem (1) can obtain a solution by solving the following least square problem: $\min_{x} {∥ Ax - b ∥}_{2}^{2}$ (2)

It can perform well by iterative algorithms, such as ART and SART. SART can suppress the noise when the projection data are complete and contains low-level noise, whereas it will lead to obvious artifacts in reconstructed images when projection data are incomplete or contains the high-level noise. In this situation, the solution generated from SART will not converge to the correct solution.

The regularization methods offer an effective way to surmount the instability of the solution of an ill-posed problem. Nowadays, TV regularization [10] as a state-of-the-art representative can be expressed as: $\min_{x \geq 0} {∥ Tx ∥}_{1}, s . t . {∥ Ax - b ∥}_{2} \leq ɛ$ (3)

where ∥Tx ∥ ₁ denotes the total variation of the image x, x_i,j is the intensity of the image x at the position (i,j), λ is a parameter balancing the fidelity term and regularization term,T denotes the gradient transform and (Tx) _i,j = (x_i+1,j - x_i,j+1 - x_i,j). The model (3) puts emphasis on the global characteristics of the transformed image, and its constraint according to the image is posed on the whole. TV has ability to deal with several reconstruction problems. However, if the projection data contain the large noise and the object to be reconstructed is relatively complex, TV may lead to stair artifacts and obscure the edges.

For the sake of addressing low SNR reconstruction, this paper proposes a novel model (4) that aims to include group sparsity characteristics of the gradient of the reconstructed image. The presented model combines the Gaussian kernel and overlapping group sparsity into TV model (3), which can be noted as GOGS-TV. The presented model is expressed as follows: $\min_{x \geq 0} λ \sum_{k = 1}^{r} \sum_{(i, j) \in G_{k}^{I}} ω_{i, j} (pk) \cdot | {(Tx)}_{i, j} |, s . t . {∥ Ax - b ∥}_{2} ⩽ ɛ$ (4)

where the symbols of the presented model (4) can be explained as follows: I denotes the group size and Ω denotes the whole coordinate positions of the image x; $G_{k}^{I}, k = 1, 2, . . ., r$ are the I-point groups such that $⋃_{k = 1}^{r} G_{k}^{I} = Ω$ ; p_k denotes the central point of the group $G_{k}^{I}$ ; p_i,j is (i,j) coordinate position of the image x; ω_i,j (p_k) = e^{-d²(p_k,p_i,j)/-h²} is the weighted factor, i.s. Gaussian kernel; d (p_k, p_i,j) is the distance between the point p_k and p_i,j; h is the smoothing parameter.

2.2 Method

In this section, it will briefly introduce the methods to address the models (2) and (3), respectively, then dilate on the algorithm to deal with the proposed model (4). To address the model (2), SART [8] is utilized. In this situation, let the image x be arranged into one vector. The system matrix A belongs to the real space R^M^×N, i.e. A is M×N matrix and a_i,j denotes the ith row and jth column element of the matrix A. Likewise, the projection data b belongs to R^M and the image x belongs to R^N. The (n + 1)th iteration is expressed as follows: $x_{j}^{(n + 1)} = x_{j}^{n} + \frac{ϖ}{A_{+, j}} \sum_{i = 1}^{M} \frac{a_{i, j}}{A_{i, +}} (b_{i} - A_{j, +} x^{n}), j = 1, 2, \dots, N$ (5)

where ϖ denotes a relaxation parameter in (0,2), $A_{i, +} = \sum_{j = 1}^{N} a_{i, j}, i = 1, 2, \dots, M$ and $A_{+, j} = \sum_{i = 1}^{M} a_{i, j}, j = 1, 2, \dots, N$ .

For convenience, we let x⁽ⁿ⁺¹⁾ = SART (xⁿ) instead of the equation (5).

For the sake of solving the model (3), SART is first utilized to deal with the fidelity term ∥Ax - b ∥ ₂ ⩽ ɛ, and that is ${\tilde{x}}^{(n)} = SART (x^{n})$ . Then the gradient of the regularization term ∥Tx ∥ ₁ is computed. In order to prevent the singular value of it, a small positive parameter τ is brought in, which is set to 10^–8 in this paper. Then the gradient of ∥Tx ∥ ₁ [24] according to variable x at x_i,j can be expressed as follows:

${v_{i, j}^{n} = \frac{\partial {∥ Tx ∥}_{1}}{\partial x_{i, j}} |}_{x = x^{n}} \approx \frac{(x_{i, j}^{n} - x_{i - 1, j}^{n}) + (x_{i, j}^{n} - x_{i, j - 1}^{n})}{\sqrt{{(x_{i, j}^{n} - x_{i - 1, j}^{n})}^{2} + {(x_{i, j}^{n} - x_{i, j - 1}^{n})}^{2} + τ}} - \frac{x_{i + 1, j}^{n} - x_{i, j}^{n}}{\sqrt{{(x_{i + 1, j}^{n} - x_{i, j}^{n})}^{2} + {(x_{i + 1, j}^{n} - x_{i + 1, j - 1}^{n})}^{2} + τ}} - \frac{x_{i, j + 1}^{n} - x_{i, j}^{n}}{\sqrt{{(x_{i, j + 1}^{n} - x_{i, j}^{n})}^{2} + {(x_{i, j + 1}^{n} - x_{i - 1, j + 1}^{n})}^{2} + τ}}$ (6)

Finally, the gradient descent method is used to obtain the reconstructed result. The (n + 1)th iterative result can be formulated as: $x_{i, j}^{n + 1} = x_{i, j}^{n} - λ ν_{i, j}^{n}$ (7)

Before addressing the model (4), let’s analyze the various situations of (4). Its three main forms can be represented as follows:

If $⋂_{k = 1}^{r} G_{k}^{I} = Φ$ and ω_i,j (p_k) = 1 for all i and j, then the model (4) is identical to the model (3);

If $⋂_{k = 1}^{r} G_{k}^{I} \sim = Φ$ and ω_i,j (p_k) = e^{-d²(p_k,p_i,j)/-h²} for each i and j, it exhibits the weighted and non-overlapping group sparsity modality;

If there exists k₁ and k₂ such that G_{k
₁} ∩ G_{k
₂} ≠ Φ and ω_i,j (p_k) = e^{-d²(p_k,p_i,j)/-h²}, then it denotes the weighted and overlapping group sparsity modality.

In this study, it mainly deals with the third case. Just like TV algorithm, SART is used to address the fidelity term, i.e. ${\tilde{x}}^{n} = SART (x^{n})$ . Then the gradient of the regularization term $\sum_{k = 1}^{r} \sum_{(i, j) \in G_{k}^{I}} ω_{i, j} (p_{k}) \cdot | {(Tx)}_{i, j} |$ is computed according to variable x. For simplicity, let $G (x) = \sum_{k = 1}^{r} \sum_{(i, j) \in G_{k}^{I}} ω_{i, j} (p_{k}) \cdot | {(Tx)}_{i, j} |$ and ${d_{i, j}^{n} = \frac{\partial G (x)}{\partial x_{i, j}} |}_{x = {\tilde{x}}^{n}}$ . Then the (n + 1)th iterative result can be expressed as: $x_{i, j}^{n + 1} = x_{i, j}^{n} - λ d_{i, j}^{n}$ (8)

To summarize above the procedures of our method, let N_iter be the maximum iteration number and N_GOGS - TV be the maximum iteration number of GOGS-TV minimization step, and the pseudo code of the GOGS-TV method is exhibited in the following:

Input: b, x₀, λ, ϖ, h, τ, ɛ₀, n ≤ N_iter, N_GOGS - TV,

Initialization: x⁰ = x₀, n ← 0.

While n ≤ N_iter and ɛ ≤ ɛ₀ do

step 1. SART Updating:

1) Compute ${\tilde{x}}^{n}$ in Eq. (5).

2) Non-negative constraint:

${\tilde{x}}_{i, j}^{n} = {\begin{matrix} {\tilde{x}}_{i, j}^{n}, {\tilde{x}}_{i, j}^{n} ⩾ 0 \\ 0, {\tilde{x}}_{i, j}^{n} < 0 \end{matrix}, \forall (i, j) \in Ω$

step 2. GOGS- TV Minimization:

Initialization: m = 1

While m ≤ N_GOGS - TV and (i, j) ∈ Ω

$y_{i, j}^{m} = {\tilde{x}}_{i, j}^{n},$

${d_{i, j}^{n} = \frac{\partial G (x)}{\partial x_{i, j}} |}_{x = {\tilde{x}}^{n}},$

$y_{i, j}^{m} \leftarrow y_{i, j}^{m} - λ d_{i, j}^{n},$

$m \leftarrow m + 1 .$

End While

$x_{i, j}^{n + 1} = y_{i, j}^{N_{GOGS - TV}} .$

$n \leftarrow n + 1, ɛ = \frac{{∥ x^{n + 1} - x^{n} ∥}_{2}}{{∥ x^{n + 1} ∥}_{2}}$

End While

Output: Final image xⁿ

3 Numerical implementations

In this section, simulated experiments for low SNR reconstruction are implemented to validate and evaluate the effectiveness of GOGS-TV approach. Then the results using the proposed approach are compared with the classical SART and state-of-the-art TV algorithm which has been widely used for several reconstruction problems. All the experiments are implemented on 2.40 GHz intel(R) Core (TM) i5-4210U CPU processor with 4G memory.

3.1 Quantitative assessments

In the experiments, it conducts visual inspection for reconstructed images using the proposed approach, FBP [1], SART [8] and TV [11], respectively. To exhibit the details and edges in the reconstructed images, the results magnify certain selected regions of the reconstructed images. Additionally, quantitative assessments are utilized to characterize the reconstruction quality by considering the root mean squared error (RMSE) [25], peak signal to noise ratio (PSNR) [26] and structural similarity index (SSIM) [27] as follows. $RMSE = \sqrt{\frac{1}{M \times N} \sum_{i = 1}^{N} \sum_{j = 1}^{M} {(x_{i, j} - x_{i, j}^{r})}^{2}}$ (9) $PSNR = 10 lo g_{10} \frac{{∥ x^{r} ∥}_{\infty}^{2}}{RMS E^{2}}$ (10) $SSIM (x, x^{r}) = \frac{2 μ_{x} μ_{x^{r}} (2 Cov {x, x^{r}} + c_{2})}{(μ_{x}^{2} + μ_{x^{r}}^{2} + c_{1}) (σ_{x}^{2} + σ_{x^{r}}^{2} + c_{2})}$ (11)

Where x denotes the reconstructed image; x^ris the reference image; M×N is the total number of the image pixels; μ_x and μ_{x
^r} are the mean values of _x and x^r; σ_x and σ_{x
^r} are the standard deviations of x and x^r; Cov{x,x^r} denotes the covariance between x and x^r; c₁ = (0.01R)² and c₂ = (0.03R)², R is the dynamic range of pixel values. The value of SSIM is closer to 1, then the similarity between the reference image and the reconstructed image is higher.

3.2 Simulated chest phantom results

3.2.1 Simulated scanning parameters and parameter selections for the methods

A chest phantom [27] is used to test the proposed approach. The simulated geometrical scanning parameters for low SNR CT are listed in Table 1. this paper adds the Gaussian noise to the projection data such that it can verify the effectiveness of the proposed approach, where the average value of the Gaussian noise is zero, and its standard deviation is 1% and 2% maximum projection in current view (denoted as 1% max and 2% max). The experiments take the complete range 0°∼360° into consideration, and the iterations are all set to 50, i.s. N_iter = 50. The relative error ε₀ is set to 10^–6.

Table 1
Geometrical scanning parameters for simulated CT imaging system

Parameter Value

The distance between source and rotation center 981 mm

The diameter of field of view 279.5 mm

The angle interval of projection views 0.705⁰

The numbers of detector bin 560

The angle interval of rays 0.0005⁰

Pixel size 0.5×0.5 mm²

Image size 512×512

Parameter	Value
The distance between source and rotation center	981 mm
The diameter of field of view	279.5 mm
The angle interval of projection views	0.705⁰
The numbers of detector bin	560
The angle interval of rays	0.0005⁰
Pixel size	0.5×0.5 mm²
Image size	512×512

In the simulated experiments, the parameters are chosen by trial and error for better image quality. For different degrees of low SNR projection data and the methods, the corresponding parameters are selected differently. However, the parameter selections are in connection with the noise level of the projection data. The specific choice in this paper is shown below.

For SART [8], the relaxation parameter ϖ is both set to 0.05 for the two simulated situations. In general, ϖ is set to 1, but for the low SNR situation, it will lead to high noise in the reconstructed results. For TV algorithm [11], there are three main parameters to be confirmed which are the relaxation parameter ϖ, the regularization parameter λ and the iteration number of TV minimization step (N_TV). For both simulated situations, N_TV = 20 and ϖ = 0.1. Separately, λ is set to 0.1 for the 1% max, and λ is set to 0.2 for the 1% max.

For GOGS-TV algorithm, five parameters are needed to be confirmed, which are the size of the group I, the smoothing parameter h, the relaxation parameter ϖ, the regularization parameter λ and the iteration number of GOGS-TV minimization step (N_GOGS - TV). For both simulated experiments, there are four cases: I = 3×3, h = 0.4; I = 3×3, h = 0.5; I = 5×5, h = 0.4; I = 5×5, h = 0.5. In these situations, ϖ = 0.1 and N_GOGS - TV = 20 for all situations, meanwhile, λ = 0.05 for the 1% max and λ = 0.09 for the 2% max.

3.2.2 Simulated chest results

Figure 3 (a) shows the reference image. The 1st and 2nd red rectangles of Fig. 3 (a) are to be locally zoom-in regions such that it can further confirm the effectiveness of the proposed approach GOGS-TV. The rest of Fig. 3 are the results reconstructed using FBP, SART, TV and GOGS-TV from projection data with added 1% max Gaussian noise. It can be observed that although Fig. 3 (b) and (c) has a good contrast, but it contains a lot of noise. Thus, sometimes it may affect some small details. Figure 3 (d) shows the result using TV algorithm, which inhibits the noise very well, but simultaneously it does not preserve the edges and details very well.

Fig. 3

(a) is the reference image and its red rectangles stand for the locally zoom-in positions. (b), (c) and (d) are the results reconstructed using FBP, SART and TV from the 1% max projection data, respectively. (e1), (e2), (e3) and (e4) represent the results using GOGS-TV with different parameters from the 1% max projection data. (e1): I = 3×3, h = 0.4; (e2): I = 3×3, h = 0.5; (e3): I = 5×5, h = 0.4; (e4): I = 5×5, h = 0.5. Grayscale window is [0,255].

From Fig. 3 (e1) to Fig. 3 (e4), they are reconstructed by GOGS-TV with the different parameters. In vision, they all have ability to suppress the noise and preserve the edges except for some tiny places. To better demonstrate the details, Fig. 4 and Fig. 5 exhibit the local zoom-in images of Fig. 3. It can be found that the results using GOGS-TV have better ability to suppress the noise while preserving the edges, compared with the classical SART and the state-of-the-art TV algorithm. For better distinguishing the reconstruction effects, the quantitative assessments are implemented as shown in Table 2. From Table 2, it is seen that for the 1% max situation GOGS-TV (I = 5×5, h = 0.4) can obtain the best assessments.

Fig. 4

It shows the locally zoom-in images corresponding to the 1st red remark of Fig. 3 (a). Grayscale window is [0,255].

Fig. 5

It shows the locally zoom-in images corresponding to the 2nd red remark of Fig. 3 (a). Grayscale window is [0,255].

Table 2

Quantitative assessments for the 1% max situation

Algorithm	RMSE	PSNR	SSIM
FBP	21.7961	21.3632	0.6409
SART	9.4374	28.6337	0.8590
TV	5.0265	34.1055	0.9543
GOGS-TV (I = 3×3, h = 0.4)	4.0115	36.0646	0.9644
GOGS-TV (I = 3×3, h = 0.5)	4.0521	35.9773	0.9640
GOGS-TV (I = 5×5, h = 0.4)	3.9244	36.2553	0.9647
GOGS-TV (I = 5×5, h = 0.5)	4.0690	35.9411	0.9621

In order to better verify the effectiveness of GOGS-TV for addressing low SNR, this paper added higher level noise to the projection data. Figure 6 (a) shows the reference image. The 1st and 2nd red rectangles of Fig. 6 (a) are to be locally zoom-in regions which can further confirm the superiority of GOGS-TV. The rest of Fig. 6 are the results reconstructed using FBP, SART, TV and GOGS-TV from projection data with added 2% max Gaussian noise. From Fig. 6, it can be observed that Fig. 6 (a) and (b) represents the result using SART and FBP respectively, which contains so much noise that its edges and details are distorted. In comparison, TV algorithm can obtain a better result as shown in Fig. 6 (d). However, due to the projection data with the high-level noise, TV algorithm brings about the stair artifacts and makes the edges fuzzy.

Compared to SART and TV algorithm, GOGS-TV has more advantages to deal with the low SNR projection data. From Fig. 6 (e1) to (e4), it can be observed that they have the better ability to preserve the edges and suppress the noise, especially Fig. 6. (e4). Table 3 shows that the result using GOGS-TV (I = 5×5, h = 0.5) can obtain the best quantitative assessments.

Fig. 6

(a) is the reference image and its red rectangles stand for the locally zoom-in positions. (b), (c) and (d) are the results reconstructed using FBP, SART and TV from the 2% max projection data, respectively. (e1), (e2), (e3) and (e4) represent the results using GOGS-TV with different parameters from the 2% max projection data. (e1): I = 3×3, h = 0.4; (e2): I = 3×3, h = 0.5; (e3): I = 5×5, h = 0.4; (e4): I = 5×5, h = 0.5. Grayscale window is [0,255].

Table 3

Quantitative assessments for the 2% max situation

Algorithm	RMSE	PSNR	SSIM
FBP	25.4449	20.0188	0.4631
SART	19.0365	22.5391	0.7262
TV	6.7102	31.5961	0.9071
GOGS-TV (I = 3×3, h = 0.4)	6.2360	32.2327	0.9204
GOGS-TV (I = 3×3, h = 0.5)	5.6059	33.1579	0.9162
GOGS-TV (I = 5×5, h = 0.4)	5.5743	33.2071	0.9200
GOGS-TV (I = 5×5, h = 0.5)	5.3951	33.4908	0.9204

3.3 Simulated few-view and low-dose chest results

In order to verify the effectiveness of the presented approach for the few-view and low-dose CT reconstruction, this paper also utilizes the Chest phantom. The geometrical scanning parameters for simulated CT imaging system is as shown in Table 1. There are three experiments with 20, 30, 40 projection views to be done. The projection data is added by 1% max Gaussian noise. The total number of iterations N_itet is set to 500.

In these experiments, the parameters are chosen by trial and error for better image quality. For different projection views, the corresponding parameters are selected differently. For SART, the relaxation parameter ϖ is both set to 0.5 for the three simulated situations. For TV algorithm, there are three main parameters to be confirmed which are the relaxation parameter ϖ, the regularization parameter λ and the iteration number of TV minimization step (N_TV). For both simulated situations, N_TV = 20 and ϖ = 0.5. λ is set to 0.8, 1 and 1.1 for 20, 30, 40 projection views separately.

For GOGS-TV algorithm, I = 3×3, h = 0.4; I = 3×3, h = 0.5; I = 5×5, h = 0.4; I = 5×5, h = 0.5;. In these situations, ϖ = 0.5 and N_GOGS - TV = 20. Meanwhile, λ is set to 0.78, 0.96 and 1.08 for 20, 30 and 40 projection views, respectively. The results can be shown in Fig. 9, Fig. 10 and Fig. 11. It can be shown that the results reconstructed by FBP and SART contain a lot of noise and the details and edges are distorted; meanwhile, TV can work better than FBP and SART, but they have stair artifacts and some details has been distorted. The reason is that firstly the Chest phantom has a few details and the grayscale of each field is similar; secondly, the obtained projection data contains a high level of noise. Then, when dealing with the noise, the details are distorted. On the contrary, when preserving the details, the noise will affect the reconstructed results. Compared with previous algorithms, GOGS-TV has the ability of address the low-dose with few-view projection data in terms of quantitative indicators as shown in Table 4. But it’s worth noting that the reconstructed results gradually become worse as the few-view number decreasing. The numerical assessments show that the results by GOGS-TV with I = 3×3, h = 0.4; can get the best RMSE, PSNR and SSIM.

Fig. 7

It shows the locally zoom-in images corresponding to the 1st red remark of Fig. 6 (a), respectively. Grayscale window is [0,255].

Fig. 8

It shows the locally zoom-in images corresponding to the 2nd red remark of Fig. 6 (a), respectively. Grayscale window is [0,255].

Fig. 9

It shows the results reconstructed from 40 projection views. (a) is the reference image. (b), (c) and (d) are the results reconstructed using FBP, SART and TV from 40 projection views with added the 1% max Gaussian noise, respectively. (e1), (e2), (e3) and (e4) represent the results using GOGS-TV with different parameters from the 1% max projection data. (e1): I = 3×3, h = 0.4; (e2): I = 3×3, h = 0.5; (e3): I = 5×5, h = 0.4; (e4): I = 5×5, h = 0.5. Grayscale window is [0,255].

Fig. 10

It shows the results reconstructed from 30 projection views. (a) is the reference image. (b), (c) and (d) are the results reconstructed using FBP, SART and TV from 30 projection views with added the 1% max Gaussian noise, respectively. (e1), (e2), (e3) and (e4) represent the results using GOGS-TV with different parameters from the 1% max projection data. (e1): I = 3×3, h = 0.4; (e2): I = 3×3, h = 0.5; (e3): I = 5×5, h = 0.4; (e4): I = 5×5, h = 0.5. Grayscale window is [0,255].

Table 4

Quantitative assessments for few-view and low-dose Chest phantom

Views	Algorithm	RMSE	PSNR	SSIM
	FBP	35.1236	17.2188	0.2896
	SART	16.1342	23.9758	0.7926
	TV	8.8212	29.2202	0.8742
40	GOGS-TV (I = 3×3, h = 0.4)	8.7598	29.2809	0.8749
	GOGS-TV (I = 3×3, h = 0.5)	8.8247	29.2168	0.8719
	GOGS-TV (I = 5×5, h = 0.4)	8.8083	29.2330	0.8751
	GOGS-TV (I = 5×5, h = 0.5)	8.7715	29.2693	0.8731
	FBP	∼37.4319	16.6660	0.2613
	SART	17.8733	23.0867	0.7706
	TV	9.8844	28.2318	0.8612
30	GOGS-TV (I = 3×3, h = 0.4)	9.8130	28.2948	0.8637
	GOGS-TV (I = 3×3, h = 0.5)	9.9331	28.1891	0.8593
	GOGS-TV (I = 5×5, h = 0.4)	9.8717	28.2430	0.8630
	GOGS-TV (I = 5×5, h = 0.5)	9.8855	28.2308	0.8589
	FBP	41.9822	15.6695	0.2402
	SART	20.5335	21.8815	0.7292
	TV	11.3527	27.0288	0.8429
20	GOGS-TV (I = 3×3, h = 0.4)	11.3463	27.0337	0.8437
	GOGS-TV (I = 3×3, h = 0.5)	11.3585	27.0244	0.8419
	GOGS-TV (I = 5×5, h = 0.4)	11.3838	27.0050	0.8426
	GOGS-TV (I = 5×5, h = 0.5)	11.3345	27.0428	0.8415

3.4 Walnut real results

To further confirm the effectiveness and accuracy of the presented approach for low-dose CT reconstruction. Some details of the Walnut are difficult to reconstruct, which includes a dense, layered shell. The Walnut is acquired from a CT system with tube current 200μA and tube voltage 80 kv. Its projection data and geometrical parameters can be obtained at reference [28]. The scanning range [0,181°] is investigated, and the corresponding projection views are 604. The size of reconstruction image is 656×656.

Similarly, the parameter selections are by trial and error. For SART, the relaxation parameter ϖ is both set to 0.2 for the real Walnut. For TV algorithm, there are three main parameters to be confirmed which are the relaxation parameter ϖ = 1, the regularization parameter λ = 0.2 and the iteration number of TV minimization step (N_TV = 20). For GOGS-TV algorithm, five parameters are needed to be confirmed, which are the size of the group I, the smoothing parameter h, the relaxation parameter ϖ, the regularization parameter λ and the iteration number of GOGS-TV minimization step (N_GOGS - TV). For the Walnut projection data, there are four cases: I = 3×3, h = 0.4; I = 3×3, h = 0.5; I = 5×5, h = 0.4; I = 5×5, h = 0.5;. In these situations, ϖ = 0.2, λ = 0.0001 and N_GOGS - TV = 20.

From the reconstructed results in Fig. 12, FBP can reconstruct the relatively good details but include the noise. It shows that SART, TV and GOGS-TV can obtain the good results. From some small details (such as the center of the Walnut), GOGS-TV can address that better.

Fig. 11

It shows the results reconstructed from 20 projection views. (a) is the reference image. (b), (c) and (d) are the results reconstructed using FBP, SART and TV from 20 projection views with added the 1% max Gaussian noise, respectively. (e1), (e2), (e3) and (e4) represent the results using GOGS-TV with different parameters from the 1% max projection data. (e1): I = 3×3, h = 0.4; (e2): I = 3×3, h = 0.5; (e3): I = 5×5, h = 0.4; (e4): I = 5×5, h = 0.5. Grayscale window is [0,255].

Fig. 12

(a) is the result reconstructed by FBP. (b) and (c) are the results reconstructed using SART and TV from the real Walnut projection data, respectively. (d1), (d2), (d3) and (d4) represent the results using GOGS-TV with different parameters. (d1): I = 3×3, h = 0.4; (d2): I = 3×3, h = 0.5; (d3): I = 5×5, h = 0.4; (d4): I = 5×5, h = 0.5. Grayscale window is [0,1].

4 Discussion and conclusion

In recent decades, regularization based iterative image reconstruction algorithms have been widely developed and demonstrated their capable of addressing many different imaging problems. TV as one of the most powerful tools has been extensively utilized in low SNR situation. However, TV yields the undesirable stair artifacts when the projection data contains the high-level noise, since it puts emphasis on the global characteristics of the gradient image and it constrains the whole image.

In order to capture the shortcomings of the TV algorithm, the idea of group sparsity is utilized in this paper. This not only considers the sparsity of the gradient image, but also include the structured sparsity. Meanwhile, the Gaussian kernel is utilized to be the weighted factor, which can be adapted to the global distribution. The simulated experiments verify the superiority of the proposed approach GOGS-TV, compared to the classical FBP, SART and state-of-the-art TV algorithm. Furthermore, quantitative assessments similarly confirm this point. As adding the high-level noise, GOGS-TV still owns the merits of suppressing the noise and preserving the edges. However, for the different noise levels, the parameters of GOGS-TV need to be redetermined. In the future, the adaptivity of the parameters of the approach will be considered.

Footnotes

Acknowledgments

This work is supported by the Science and Technology Research Program of Chongqing Municipal Education Commission (Grant No. KJQN2019013), the Scientific Research Foundation of Chongqing University of Arts and Sciences (Grant No. R2019FSC17), the Natural Science Foundation of Chongqing Municipal Science and Technology Commission (Grant numbers: cstc2019jcyj-msxmX0012, cstc2020jcyj-msxm2352), and the Open Project of Key Laboratory No.CSSXKFKTQ202004, School of Mathematical Sciences, Chongqing Normal University.

References

Panetta

, Niccoló Camarlinghi, 3D Image Reconstruction for CT and PET: A Practical Guide with Python, CRC Press, 2020.

Siregar

E.S.B.

, Sutapa

G.N.

and Sudarsana

I.W.B.

, Analysis of radiation dose of patients on CT scan examination using Si-INTAN application, BULETIN FISIKA 21(2) (2020), 53–59.

Khan

S.U.

, Khan

I.U.

, et al., A review of airport dual energy X-ray baggage inspection techniques: image enhancement and noise reduction, Journal of X-ray Science and Technology 28(3) (2020), 481–505.

Vedantham

, Tseng

H.W.

, Konate

, et al., Dedicated cone-beam breast CT using laterally-shifted detector geometry: Quantitative analysis of feasibility for clinical translation, Journal of X-ray Science and Technology 28(3) (2020), 405–426.

Hernandez

A.M.

, Abbey

C.K.

, Ghazi

, et al., Effects of kV, filtration, dose, and object size on soft tissue and iodine contrast in dedicated breast CT, Medical Physics 47(7) (2020), 2869–2880.

Hsieh

, Adaptive streak artifact reduction in computed tomography resulting from excessive X-ray photon noise, Medical Physics 25(11) (1998), 2139–2147.

Xia

J.Q.

, Lo

J.Y.

, Yang

, et al., Dedicated breast computed tomography: volume image denoising via a partial-diffusion equation based technique, Medical Physics 35(5) (2008), 1950–1958.

Chee

, O’Connell

Y.M.

, et al., McSART: an iterative model-based, motion-compensated SART algorithm for CBCT reconstruction, Physics in Medicine and Biology 64(9) (2019), 095013.

Park

, Lee

, et al., Strategies for rapid reconstruction in 3D MRI with radial data acquisition: 3D fast Fourier transform vs two-step 2D filtered back-projection, Scientific Reports 10(1) (2020), 13813.

10.

Engl

H.W.

and Ramlau

, Regularization of Inverse Problems, Encyclopedia of Applied and Computational Mathematics, Springer Berlin Heidelberg, 2015

11.

Kudo

, Yamazaki

, Nemoto

, et al., A very fast iterative algorithm for TV-regularized image reconstruction with applications to low-dose and few-view CT, Proc. SPIE, 2016, 996711.

12.

, Sang

, Zhang

, et al., Few-view CT image reconstruction using improved total variation regularization, Journal of X-ray Science and Technology 27(4) (2019), 739–753.

13.

Zhao

, Gao

, Ding

, et al., Tight-frame based iterative image reconstruction for spectral breast CT, Medical Physics 40(3) (2013), 031905.

14.

Garduño

, Herman

G.T.

and Davidi

, Reconstruction from a Few Projections by l1-Minimization of the Haar Transform, Inverse Problems 27(5) (2011), 55006–55018.

15.

Luo

, Liu

, Li

, et al., Image Reconstruction from Sparse Projections Using S-Transform, Journal of Mathematical Imaging and Vision 43(3) (2012), 227–239.

16.

, Zeng

, Yu

, et al., Noise suppression-guided image filtering for low-SNR CT reconstruction, Medical & Biological Engineering & Computing 58 (2020), 2621–2629.

17.

Chen

, Constrained TV-minimization image reconstruction for industrial CT system, AIP Conference Proceedings 1581 (2014), 1856–1860.

18.

Knoll

, Bredies

, Pock

, et al., Second order total generalized variation (TGV) for MRI, Magnetic Resonance in Medicine 65(2) (2011), 480–491.

19.

Selesnick

I.W.

and Chen

P.Y.

, Total variation denoisingwith overlapping group sparsity, IEEE International Conference on Acoustics, Speech and Signal Processing. IEEE (2013), 5696–5700.

20.

Liu

, Huang

T.Z.

, Selesnick

I.W.

, et al., Image restoration using total variation with overlapping group sparsity, Information Sciences 295 (2015), 232–246.

21.

X.G.

, Jiang

and Liu

, Deblurring Poisson noisy images by total variation with overlapping group sparsity, Applied Mathematics & Computation 289 (2016), 132–148.

22.

Yang

, Wu

and Jia

, Image reconstruction for electrical impedance tomography using enhanced adaptive group sparsity with total variation, IEEE Sensors Journal 17 (2017), 5589–5598.

23.

Liu

S.J.

, et al., MRI reconstruction via enhanced group sparsity and nonconvex regularization, Neurocomputing 272 (2018), 108–121.

24.

Shin

J.W.

, Yoon

, Min

B.J.

, et al., A novel copyright protection for digital images using the gradient of image intensity, 2007 International Symposium on Information Technology Convergence, 2007.

25.

Storath

, Weinmann

, et al., Joint image reconstruction and segmentation using the Potts model, Inverse Problems 31(2) (2014), 25003–25031.

26.

Chow

L.S.

and Paramesran

, Review of medical image quality assessment, Biomedical Signal Processing & Control 27 (2016), 145–154.

27.

Cchest and abdomen phantom website. Available: http://lgdv.cs.fau.de/External/vollib/.

28.

Hmlinen

, Harhanen

, Kallonen

, et al., Tomographic X-ray data of a walnut, arXiv: 1502.04064 (2015).