Investigation of prior image constrained compressed sensing-based spectral X-ray CT image reconstruction

Zhou Zhengdong Yu Zili Zhang Wenwen Guan Shaolin

(Department of Nuclear Science and Engineering, Nanjing University of Aeronautics and Astronautics, Nanjing 210016, China)

Abstract:To improve spectral X-ray CT reconstructed image quality, the energy-weighted reconstructed image and the separable paraboloidal surrogates (SPS) algorithm are proposed for the prior image constrained compressed sensing (PICCS)-based spectral X-ray CT image reconstruction. The PICCS-based image reconstruction takes advantage of the compressed sensing theory, a prior image and an optimization algorithm to improve the image quality of CT reconstructions. To evaluate the performance of the proposed method, three optimization algorithms and three prior images are employed and compared in terms of reconstruction accuracy and noise characteristics of the reconstructed images in each energy bin. The experimental simulation results show that the image is the best as the prior image in general with respect to the three optimization algorithms; and the SPS algorithm offers the best performance for the simulated phantom with respect to the three prior images. Compared with filtered back-projection (FBP), the PICCS via the SPS algorithm and as the prior image can offer the noise reduction in the reconstructed images up to 80.46%, 82.51%, 88.08% in each energy bin, respectively. Meanwhile, the root-mean-squared error in each energy bin is decreased by 15.02%, 18.15%, 34.11% and the correlation coefficient is increased by 9.98%, 11.38%, 15.94%, respectively.

Key words:spectral X-ray CT; prior image; compressed sensing; optimization algorithm; image reconstruction

With the development of detector technology, spectral X-ray CT has become an emerging technology recently. Spectral X-ray CT equipped with photon counting detector provides energy-resolved information about objects, which has become a research focus in recent years[1-3]. One challenge in spectral CT is the inherent trade-off between image noise and energy resolution. The better energy resolution that is required, the more energy bins required, and each energy bin becomes smaller. Therefore, the number of photons per energy bin decreases, and image noise becomes higher. The generally used filtered back-projection (FBP) reconstruction algorithm meets the growing clinical requirements of image quality with difficulty.

The introduction of the prior image to the constrained compressed sensing for image reconstruction was proposed recently[4-7], which has been shown to have the capability of suppressing image noise effectively. The usage of the prior image constrained compressed sensing (PICCS) method has been investigated for spectral X-ray CT image reconstruction recently[8], where the prior image was generated by the FBP algorithm, and the performance of both steepest descent (SD) and conjugate gradient (CG) algorithms is evaluated. The experimental results show a good improvement in terms of noise level and error of reconstructed images.

In this paper, an alternative optimization method, i.e. the separable paraboloidal surrogates (SPS) algorithm[9-12] and the image derived by the energy weighted image reconstruction method[13-14]are proposed to quantitatively investigate the performance of PICCS-based spectral X-ray CT image reconstruction. The energy-weighted reconstructed image can provide a perfect prior image and the SPS optimization algorithm is a fully parallelizable algorithm, which is suitable for PICCS-based image reconstruction.

1 Methods and Materials

1.1 Objective function of PICCS algorithm

1.1.1 Prior image

The FBP algorithm was employed for the image reconstruction of each energy bin. The image was derived by the energy weighted image reconstruction method, i.e. a linear combination of the FBP images of each energy bin:

(1)

where m is the index of energy bin; and wm are the FBP image and weight of the m-th energy bin, respectively; wm is proportional to the contrast-to-noise variance ratio in the FBP image of each energy bin, which can be expressed as[13]

(2)

where σm is the noise standard deviation of ; is the average attenuation coefficient difference of contrast material xc,m and background material xb,m in .

1.1.2 Unconstrained objective function

According to the compressed sensing (CS) theory, the signal can be recovered from a limited number of samples by exploiting the sparse nature of the signals in a transform domain. The recovery of an image with sparsity can be described as

s.t. Ax=y x∈RMN×1

(3)

where is the reconstructed image; ψ(x) is a sparsifying transform of image x; A is the system matrix that relates projection measurements vector y=log(f/n) to image vector x; f is the incident fluence per detector and n is the measured projections. l can be expressed by the total variation v of x. Given the image discretization x∈RMN×1 and Λ={Mn:n∈N}, v(x) is defined as the l1-norm of an image bidimensional spatial gradient l2-norm:

(4)

The PICCS can be formulated as a constrained or an unconstrained minimization problem. The reconstruction of an image using the constrained approach can be formally expressed as[5-6]

s.t. Ax=y x∈RMN×1

(5)

(6)

where xp is the prior image; α∈[0,1] is a scalar that controls the relative weight of the prior image.

The above minimization problem can be converted to an unconstrained minimization problem, and the final objective function can be expressed as[5-6]

(7)

(8)

where D=diag{n}; λ is a data consistency parameter.

1.2 Minimization of objective function

The partial derivative of the objective function at the i-th pixel can be expressed as[5-6]

λ[AT(Ax-y)]

(9)

where the gradient of v(x) is

(10)

To remove the singularities in formula (10), a modified total variation υε(x) can be employed as

(11)

where ε should be large enough to remove the singularities in the denominator and be small enough to preserve the shape of the function.

Since there is the gradient of the objective function fuc(x), the gradient-based optimization methods such as the steepest descent (SD) algorithm, conjugate gradient (CG) algorithm and separable paraboloidal surrogates (SPS) algorithm can be used to solve the minimization problem[8-9], among which the SD and CG algorithms include two steps: determination of the search direction d and step size η. For the SD algorithm, d is set to be the negative gradient of the object function fuc(x). For the CG algorithm, d is set to be a set of conjugate gradient of the object function fuc(x). Besides, both step sizes can be determined by a line search method, which is a simple algorithm that uses the criterion of the Wolfe condition. In this paper, all the three algorithms were performed with the same stop criterion, i.e. the relative error e of the adjacent iteration results:

(12)

where k is the iteration number; xk is the k-th iteration result; and the ideal error should be 0 for all algorithms.

1.3 Metrics of image quality evaluation

In order to evaluate the quality of the reconstructed images, several image evaluation metrics including standard deviation σ, root-mean-squared error erms and the correlation coefficient cc are used:

(13)

(14)

(15)

where j=1,2,…,M, and M is the total number of pixels; and represent the pixel values of the reconstructed image and reference image, respectively; and their corresponding mean pixel values are rec and ref, respectively. σ and erms are used to evaluate noise level and the accuracy of images in energy bins, respectively. The less σ or erms, the better the reconstructed image; the higher cc stands for better correlation and greater similarity of the texture structure of images.

Besides, the Laplacian gradient was used to reflect the change of the pixel gray level of the reconstructed image. The greater LS means the more obvious profile, and the more abundant image detail[15].Given the image discretization x∈RMN×1, Λ={MnM(n-1)+1:n∈N}, the Laplacian gradient Ls can be defined as follows:

(16)

where

xi1=3xi-xi-M-1-xi-1-xi+M-1

xi2=2xi-xi-M-xi+M

xi3=3xi-xi-M+1-xi+1-xi+M+1

1.4 Experimental dataset

Genat4 simulation Toolkit[16-17] was used to simulate spectral X-ray CT to generate projections in three energy bins. The simulations were based on single-slice fan-beam scanning. The fan angle was 20.2°, the scanned object was a cylindrical phantom; and the related cross-section parameters of the phantom are given in Tab.1. The equidistant linear detector was set to be 380 mm including 190 detector elements; the distances from the source to the center of phantom and detector were set to be 506 and 1 066 mm, respectively; the rotation angle was in the range [0,2π] and the scanning was performed every 1° with 106 entrance photons emitted by a 140 KeV X-ray spectrum. The detector pixels allowed a binning of the detected pulses into three energy bins and the energy thresholds were 58.5 and 86.5 KeV, where bin1, bin2 and bin3 are used to represent the first, second and third energy bins. The calculated proportion of incident photon quantity of each energy bin was 58.3%, 26.1% and 15.6%, respectively, according to the emission spectrum. All images were generated on a 190×190 pixel matrix, and the reference image was the noiseless FBP reconstruction image from the simulated projection data.

Tab.1 Geometric parameters and components of cylindrical phantom in cross-section

Centercoordinates/mmRadius/mmMaterialcomposites(0.0,0.0)100C5H8O2(0.0,50.0)2514.8%C+11%H+74.2%O(-30.0,-30.0)1522.5%Ca+20.4%C+9.4%H+47.7%O(40.0,-40.0)151.2%I+98.8%H2O

2 Results and Discussion

2.1 Comparison of the performance of PICCS with different prior images

The objective function was optimized by the SD, CG and SPS algorithms, respectively. For each algorithm, three different images, i.e., the FBP reconstruction image of the first energy bin ), the FBP reconstruction image of the second energy bin ) and the image ) derived by the energy weighted method, were served as the prior image for performance evaluation, where the prior image parameter α was set to be 0.5[6]. When the relative error between the adjacent iteration results was less than 10-4, the optimization procedure ended. The relative errors of the objective function values of adjacent iterations of each energy bin with SD, CG and SPS algorithms vs. the iteration number are plotted in Fig.1. For Figs.1(a), (b) and (c), was set as the prior image. For Figs.1(d), (e) and (f), was set as the prior image. For Figs.1(h), (i) and (j), was set as the prior image. The image correlation coefficients cc and the root-mean-squared errors erms relative to the reference image are listed in Tabs.2, 3 and 4, where the reference image refers to the noiseless FBP reconstruction image of the phantom. Fig.2 shows erms and cc of reconstructed images optimized by the above three optimization methods with respect to the three different prior images, i.e., and .

From the above results, the performance comparison of the three optimization methods and three prior images for the phantom developed in this paper can be concluded as follows:

1) The stability of image reconstruction with regard to the prior image. Compared with the SD and CG algorithms, the SPS algorithm has the least sensitive to the prior image, and generally achieves the best stability with regard to erms and cc, while both the CG and SD algorithms are sensitive to the prior image.

2) The image quality of the three energy bins. In general, erms and cc of reconstructed images optimized by the SPS algorithm are better than those by the SD or CG algorithms, particularly for the 3rd energy bin, and the results with the SPS algorithm are much better than those with the SD or CG algorithms. Meanwhile, image is generally the best one as the prior image for image reconstruction.

2.2 Comparison of FBP and PICCS-based image reconstruction

Based on the above experiments, the SPS algorithm and are employed to optimize the PICCS objective functions where α is set to be 0.5. The weights of the three energy bins in formula (1) are set to be 0.420 3, 0.320 2, and 0.259 5, respectively, to generate prior image . Figs.3(a),(b),(c) and Figs.3(e),(f),(g) show the FBP and PICCS reconstructed images of three energy bins, respectively; Fig.3(d) shows the prior image . From Fig.3, we can see that the PICCS algorithm has the ability to reconstruct energy bin images effectively, and the reconstructed images have clear profiles and less noise. Meanwhile, the edge information of the images is well preserved.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

(i)

Fig.1 Relative error e of the adjacent iteration results of each energy bin vs. iteration number k. (a) SD with xFBP1st as the prior image; (b) CG with xFBP1st as the prior image; (c) SPS with xFBP1st as the prior image; (d) SD with xFBP2nd as the prior image; (e) CG with xFBP2nd as the prior image; (f) SPS with xFBP2nd as the prior image;(g) SD with xWbins as the prior image; (h) CG with xWbins as the prior image; (i) SPS with xWbins as the prior image

Tab.2 Performance of three algorithms in each energy bin with as the prior image

MeasuredBin1SDCGSPSBin2SDCGSPSBin2SDCGSPSerms0.28030.28520.27000.27730.27650.24100.23520.25340.1786cc0.92820.92890.92980.93610.93590.97280.94510.91790.9896

Tab.3 Performance of three algorithms in each energy bin with as the prior image

MeasuredBin1SDCGSPSBin2SDCGSPSBin2SDCGSPSerms0.25630.26250.26420.25310.25980.23850.22740.21590.1771cc0.97670.96050.94190.97780.96180.97840.96360.96770.9918

Tab.4 Performance of algorithms in each energy bin with as the prior image

MeasuredBin1SDCGSPSBin2SDCGSPSBin2SDCGSPSerms0.25490.25100.26360.24630.24930.24130.22590.22250.1791cc0.96500.96640.94340.96820.96810.97420.95920.95150.9898

Tab.5 and Fig.4 show the quantitative evaluation of image noise level σ, erms and cc for the reconstructed image of each energy bin. Since the transmitted photons are distributed into three energy bins, the recorded photons are relatively fewer in each energy bin, thus images reconstructed with the FBP algorithm may contain a significant amount of noise. Furthermore, the noise level in different energy bins is affected by the incident photon flux,and image noise increases as the incident photon flux decreases. Compared with the FBP algorithm, the PICCS has the ability to reduce noise up to 80.46%,82.51%,88.08% in each energy bin, respectively, thus image reconstruction with the PICCS can suppress noise effectively. Moreover, erms of the PICCS reconstructed image in three energy bins is decreased by 15.02%, 18.15%, 34.11% and image cc is increased by 9.98%, 11.38%, 15.94%, respectively. Thus, the image accuracy and similarity relative to the reference image are improved obviously, and the image quality increases as the photon mean energy of each energy bin increases.

(a)

(b)

Fig.2 Variation of erms and cc of reconstructed images optimized by three optimization methods with respect to the priori mages. (a) erms; (b) cc

(a)

(b)

(c)

(d)

(e)

(f)

(g)

Fig.3 Image reconstruction with FBP and PICCS algorithm for each energy-bin. (a) FBP images of the 1st energy bin; (b) FBP images of the 2nd energy bin; (c) FBP images of the 3rd energy bin; (d) Energy-weighted reconstructed image as the prior image; (e) PICCS image of the 1st energy bin; (f) PICCS image of the 2nd energy bin; (g) PICCS image of the 3rd energy bin

Tab.5 σ, ermsand cc of images reconstructed by FBP and PICCS algorithm for each energy bin

MethodσermsccFBPbin10.17350.31020.8578bin20.20010.29480.8747bin30.25250.27180.8537PICCSbin10.03390.26360.9434bin20.03500.24130.9742bin30.03010.17910.9898

(a)

(b)

(c)

Fig.4 FBP and PICCS images in the corresponding energy bins. (a) σ; (b)erms; (c) cc

3 Conclusion

The PICCS algorithm and its performance were extensively investigated for the spectral X-ray CT image reconstruction with three unconstrained optimization algorithms and three prior images including the image derived by energy weighted image reconstruction. The experimental results demonstrate that the SPS algorithm is the best choice among the three mentioned methods to optimize PICCS objective functions. Compared with and , is generally the best choice as the prior image for image reconstruction. Compared with the FBP algorithm, the PICCS algorithm has the better capability to reduce noise and error in reconstructed images of all three energy bins, and preserve a good similarity of image quality relative to the FBP reconstructed image.

Although the usage of Geant4 for the projections simulation can achieve good approximation to the real projection data, the detector and electronic noise should also be included to accurately model the imaging process in the future investigation. Moreover, the number of energy bins and their allocations also need to be investigated to fully exploit the advantages of the PICCS algorithm.

References:

[1]Wang X, Meier D, Mikkelsen S, et al. Micro CT with energy-resolved photon-counting detectors[J]. Physics in Medicine & Biology, 2011, 56(9):2791-2816.DOI:10.1088/0031-9155/56/9/011.

[2]Yu H, Xu Q, Peng H, et al. Medipix-based spectral micro-CT[J]. Computerized Tomography Theory & Applications, 2012, 21(4): 583-596.

[3]Lee S, Choi Y N, Kim H J. A simulation study of high-resolution X-ray computed tomography imaging using irregular sampling with a photon-counting detector[J]. Nuclear Instruments & Methods in Physics Research Section A, 2013, 726(13):175-180.DOI:10.1016/j.nima.2013.05.044.

[4]Chen G H, Tang J, Leng S. Prior image constrained compressed sensing (PICCS): A method to accurately reconstruct dynamic CT images from highly under sampled projection datasets [J]. Medical Physics, 2008, 35(2):660-663.

[5]Lauzier P T, Chen G H. Characterization of statistical prior image constrained compressed sensing. Ⅰ. Application to time-resolved contrast-enhanced CT [J]. Medical Physics, 2012, 39(10):5930-5948.DOI:10.1118/1.4748323.

[6]Lauzier P T, Tang J, Chen G H. Prior image constrained compressed sensing: Implementation and performance evaluation [J]. Medical Physics, 2012, 39(1): 66-80.DOI:10.1118/1.3666946.

[7]Lauzier P T, Chen G H. Characterization of statistical prior image constrained compressed sensing (PICCS): Ⅱ. Application to dose reduction [J]. Medical Physics, 2013, 40(2): 021902.DOI:10.1118/1.4773866.

[8]Brunner S T. Prior image constrained image reconstruction in emerging computed tomography applications [D]. Madison, WI, USA: School of Medicine and Public Health, the University of Wisconsin-Madison, 2013.

[9]Sotthivirat S, Fessler J A. Image recovery using partitioned-separable paraboloidal surrogate coordinate ascent algorithms[J]. IEEE Transactions on Image Processing, 2002, 11(3):306-317.DOI:10.1109/83.988963.

[10]Fletcher R. A limited memory steepest descent method[J]. Mathematical Programming, 2010, 135(1):413-436.DOI:10.1007/s10107-011-0479-6.

[11]Xiao Y, Zhu H. A conjugate gradient method to solve convex constrained monotone equations with applications in compressive sensing[J]. Journal of Mathematical Analysis & Applications, 2013, 405(1):310-319.DOI:10.1016/j.jmaa.2013.04.017.

[12]Jian J, Han L, Jiang X. A hybrid conjugate gradient method with descent property for unconstrained optimization [J].Applied Mathematical Modeling, 2015, 39(3/4):1281-1290.DOI:10.1016/j.apm.2014.08.008.

[13]Schmidt T G. Optimal “image-based” weighting for energy-resolved CT [J].Medical Physics, 2009, 36(7):3018-3027.DOI:10.1118/1.3148535.

[14]Schmidt T G. CT energy weighting in the presence of scatter and limited energy resolution [J]. Medical Physics, 2010, 37(3):1056-1057.DOI:10.1118/1.3301615.

[15]Du F M. Research and application on clarity-evaluation-model for the cross-line gray-scale image [C]//IEEE International Conference on Computer Science & Automation Engineering. Shanghai, China,2011, 4:322-324.

[16]Agostinelli S, Allison J, Amako K, et al. Geant4—A simulation Toolkit[J]. Nuclear Instruments & Methods in Physics Research Section A, 2003, 506(3):250-303.DOI:10.1016/s0168-9002(03)01368-8.

[17]Geant4 user’s documents: Introduction to geant4 [EB/OL].(2015-12-04)[2016-09-28]. http://geant4.web.cern.ch/geant4/UserDocumentation/UsersGuides/IntroductionToGeant4/html/index.html.

基于先验图像约束和压缩感知的能谱X-CT图像重建

周正东 余子丽 张雯雯 管绍林

(南京航空航天大学核科学与工程系, 南京 210016)

摘要:为了提高能谱 X-CT重建图像的质量,提出了利用能量加权重建图像及可分离抛物面替代法进行基于先验图像和约束压缩感知的能谱X-CT图像重建.利用压缩感知理论、先验图像和优化算法来提高CT重建图像的质量.为了评价所提方法的性能,从重建的各能量段图像精度和噪声特性2个方面比较了3种优化算法及3种先验图像.仿真实验结果表明,对于不同的优化算法,能量加权重建图像作为先验图像总体性能最佳;对于不同的先验图像,可分离抛物面替代法算法性能最佳.与滤波反投影算法相比,在基于先验图像约束和压缩感知的能谱X-CT图像重建算法中,采用SPS算法进行优化,采用能量加权重建图像作为先验图像,重建得到的各能量段的图像噪声分别降低了80.46%,82.51%,88.08%,每个能量段图像的均方根误差分别下降了15.02%,18.15%和34.11%,相关系数分别提高了9.98%,11.38%和15.94%.

关键词:能谱X-CT; 先验图像; 压缩传感; 优化算法; 图像重建

中图分类号:TP391; R318

Foundation items:The National Natural Science Foundation of China (No.51575256), the Fundamental Research Funds for the Central Universities (No.NP2015101, XZA16003), the Priority Academic Program Development of Jiangsu Higher Education Institutions (PAPD).

Citation::Zhou Zhengdong, Yu Zili, Zhang Wenwen, et al. Investigation of prior image constrained compressed sensing-based spectral X-ray CT image reconstruction[J].Journal of Southeast University (English Edition),2016,32(4):420-425.

DOI:10.3969/j.issn.1003-7985.2016.04.005.

DOI:10.3969/j.issn.1003-7985.2016.04.005

Received 2016-06-09.

Biography:Zhou Zhengdong(1969—), male, doctor, associate professor, zzd_msc@nuaa.edu.cn.