Investigation of Different Sparsity Transforms for the PICCS Algorithm in Small-Animal Respiratory Gated CT
Investigation of Different Sparsity Transforms for the PICCS Algorithm in Small- Animal Respiratory Gated CT
Juan F. P. J. Abascal 0 1 2 3
Monica Abella 0 1 2 3
Alejandro Sisniega 0 1 2 3
Juan Jose Vaquero 0 1 2 3
Manuel Desco 0 1 2 3
0 1 Departamento de Bioingenieria e Ingenieria Aeroespacial, Universidad Carlos III de Madrid , Madrid , Spain , 2 Instituto de Investigacion Sanitaria Gregorio Maranon (IiSGM) , Madrid , Spain , 3 Centro de Investigacion en Red de Salud Mental (CIBERSAM) , Madrid , Spain
1 Funding: This work was partially funded by the RIC- RETIC network (RD12/0042/0057) from the Ministerio de Economia y Competitividad (www.mineco.gob.es/) and projects TEC2010-21619-C04-01 and PI11/ 00616 from Ministerio de Ciencia e Innovacion (www. micinn.es/). The research leading to these results was supported by funding from the Innovative Medicines Initiative (www.imi.europa.eu) Joint Undertaking under grant agreement n115337, the resources of which comprise financial contributions
2 Academic Editor: Qinghui Zhang, University of Nebraska Medical Center , UNITED STATES
3 Current address: Department of Biomedical Engineering, Johns Hopkins University , Baltimore, Maryland , United States of America
Respiratory gating helps to overcome the problem of breathing motion in cardiothoracic small-animal imaging by acquiring multiple images for each projection angle and then assigning projections to different phases. When this approach is used with a dose similar to that of a static acquisition, a low number of noisy projections are available for the reconstruction of each respiratory phase, thus leading to streak artifacts in the reconstructed images. This problem can be alleviated using a prior image constrained compressed sensing (PICCS) algorithm, which enables accurate reconstruction of highly undersampled data when a prior image is available. We compared variants of the PICCS algorithm with different transforms in the prior penalty function: gradient, unitary, and wavelet transform. In all cases the problem was solved using the Split Bregman approach, which is efficient for convex constrained optimization. The algorithms were evaluated using simulations generated from data previously acquired on a micro-CT scanner following a high-dose protocol (four times the dose of a standard static protocol). The resulting data were used to simulate scenarios with different dose levels and numbers of projections. All compressed sensing methods performed very similarly in terms of noise, spatiotemporal resolution, and streak reduction, and filtered back-projection was greatly improved. Nevertheless, the wavelet domain was found to be less prone to patchy cartoon-like artifacts than the commonly used gradient domain.
Respiratory gating helps to overcome the problem of breathing motion in cardiothoracic
small-animal imaging. CT imaging is the gold standard in several lung diseases, such as
tuberculosis. Blurring caused by breathing motion can hinder quantification in imaging studies,
from the European Union's Seventh Framework
Programme (FP7/2007-2013) and EFPIA companies
("in kind contribution"). The funders had no role in
study design, data collection and analysis, decision to
publish, or preparation of the manuscript.
Competing Interests: The authors have declared
that no competing interests exist.
which are useful for assessing the degree of infection based on the density and extension of the
lesions. One option for improving image quality is to correct movement blurring using
retrospective gating. If we generate complete data sets for a number of respiratory phases by
acquiring multiple images for each projection angle , the radiation dose delivered to the subject
increases proportionally to the number of respiratory phases. Fig 1(A) shows an example of 4
respiratory phases obtained with 32 images per projection angle using a high-dose protocol. If
we use only 8 frames per projection angle (corresponding to a dose similar to that of a static
imaging protocol), few noisy and irregularly distributed projections are available for the
reconstruction of each respiratory phase, thus leading to streak artifacts in the FBP-reconstructed
images (Fig 1(B)). In previous approaches [2, 3], this problem was solved within an analytical
framework using a variation of the McKinnon-Bates method , which is based on correction
of an initial estimate obtained from the whole data set (combining all respiratory phases) with
the undersampled data from each respiratory phase. Although this approach reduces the
Fig 1. Comparison of static imaging and reference high-dose protocols. FDK reconstructions of gated data obtained with the reference high-dose
protocol (A) and with the static imaging protocol (B) comprising 32 and 8 frames per projection angle, respectively. C: Absolute image difference between two
high-dose respiratory phases. D: Prior image obtained from the addition of the four low-dose respiratory phases and processed with a Gaussian filter to
presence of noise and artifacts, correcting the artifacts present in the initial estimate remains
Correction of respiratory motion has been addressed for image-guided radiotherapy using
non-rigid image registration based on a motion model [5, 6]. In order to reduce the effect of
the streak artifacts, the registration was restricted to a volume of interest, defined by a
boundary set 2 mm outside the body of the subject. Finally, the registered images were combined and
residual streak artifacts further reduced using principal component analysis. The drawback of
this approach is that it requires good image quality to guide non-rigid registration.
From a different perspective, in the compressed sensing (CS) framework, an image can be
accurately reconstructed from few projections using convex optimization, provided that the
image is sparse in a transformed domain . The most commonly used transformed
domain is the gradient that leads to total variation (TV) , which efficiently removes noise
and artifacts caused by undersampling, but leads to patchy images for high undersampling
A combination of both strategies, i.e., prior image and sparsity, is found in the so-called
prior image constrained compressed sensing (PICCS) algorithm , which enables
accurate reconstruction of highly undersampled data. PICCS combines TV, which removes noise,
with a prior image that helps to maintain a natural image texture. This prior image is usually
obtained as the average of all respiratory phases, similar to the procedure followed in the
McKinnon-Bates method . PICCS has been applied to contrast cardiac CT data  and
respiratory gated phantom data . It has also been applied to characterize breathing motion
and requires a lower acquisition time than filtered back projection with McKinnon-Bates
correction . Several works have proposed variations of PICCS: an adaptive PICCS for
longitudinal CT studies , an extension to include a log-likelihoodbased fidelity term , and a
nonconvex approach .
The PICCS algorithm is a constrained optimization based on L1-penalty functions that can
be solved using classic constrained optimization methods. However, since these methods can
be computationally expensive, most algorithms solve the constrained TV problem using
methods that alternate steepest descent for minimization of an unconstrained version of TV with
iterative methods such as the simultaneous algebraic reconstruction technique (SART), which
imposes fidelity on the acquired data. Algorithms such as adaptive steepest descent projection
onto convex sets (ASD-POCS) are based on this approach . Nevertheless, solving an
unconstrained approximated version of the constrained problem requires optimal selection of the
regularization parameter and estimation of the step size. The Split Bregman algorithm was
applied to MRI [25, 26] and proved to be optimal and computationally efficient for the solution
of constrained problems with L1-penalty functions. In addition, this approximation facilitates
the enforcement of constraints and it circumvents the requirement of an optimal selection of
the regularization parameter. Similar alternating methods were also applied to CT .
With regard to sparsity transforms, the preferred choice for CT is the gradient domain.
Although other transforms may be sparser, depending on the application, few studies have
actually used a different choice, and even fewer have offered a comparison. In , the authors
propose the shearlet transform for static CT, and in , wavelet frames were tested on
phantom data. In the case of the PICCS method, to our knowledge, no studies have evaluated
different sparsity transforms.
In this study, we compare three versions of the PICCS algorithm using different transforms
in the prior penalty function term (unitary, gradient, and wavelet). In all three cases, the
problem was solved using the Split Bregman formulation. In addition, positive and support
constraints were added to the standard PICCS method. The evaluation was performed on
smallanimal CT data in terms of contrast in bone and lung tissue, mean square error (MSE) with
respect to the target, image noise, contrast-to-noise ratio (CNR), degree of compensation of the
respiratory movement, and image texture quality. We also analyzed the performance of the
algorithm for different weights of the prior penalty term and studied different settings of X-ray
flux (related to delivered dose) and number of projections. Preliminary results were presented
earlier for a fixed flux and number of projections .
PICCS. The PICCS method can be used to reduce streak artifacts and noise when
reconstructing highly undersampled gated-CT data. With ui as the i-th phase image, PICCS assumes
that ui is sparse in a transformed domain T1 and that there must be a prior image up to ensure
that uiup is sparse in a transformed domain T2. If fi represents the data corresponding to the
ith phase image and F is the forward operator, PICCS is the convex constrained optimization
where I is the total number of respiratory phase bins, accounts for noise in the data and
weights the prior penalty function. The common choice for T1 is the spatial discrete gradient
that leads to TV, ||rui||1, which filters out noise while preserving edges in the image. TV is also
a common choice for T2.
We remark that, for = 0, the problem in Equation (1) corresponds to the minimization of
TV subject to a data constraint. TV assumes that the image is piecewise smooth and has been
shown to yield cartoonish images for a low number of projections . When 0 < 1, the
addition of the prior image prevents the excessive smoothing produced by TV and helps to
maintain the texture of the prior image.
PICCS with positivity and support constraints solved using the Split Bregman
formulation. We extended the PICCS method by adding a positivity constraint [31, 32] and a support
constraint that restrict the reconstruction to a circular field of view, O, defined by the Radon
transform. Thus, the reconstruction problem in Equation (1) becomes
where we use isotropic TV, kruik1 qffiffirffiffiffiffixffiffiuffiffiffiiffiffi2ffiffiffiffiffiffiffiffiffiffirffiffiffiffiyffiuffiffiffiiffiffiffi2ffi.
To solve the problem in Equation (2), we use the Split Bregman formulation, which
efficiently handles L1-based constrained problems [25, 33]. The Split Bregman formulation makes
it possible to split L1-norm terms and L2-norm terms in such a way that they can both be
solved analytically in two separate steps. The part including the L2-norm functionals results in
a linear system that can be solved using linear iterative methods, and the part with L1-norm
functionals is solved using shrinkage formulas, as shown below.
To allow for splitting, we include new variables, dxi, dyi, wi and vi, and formulate a new
problem that is equivalent to Equation (2)
Equation (3) is easily managed using an equivalent unconstrained optimization approach
with constraints imposed by adding a Bregman iteration bi for each constraint. That is,
where F(vi 0, vi O) represents the non-negativity and support constraints, k is the iteration
number and the Bregman iterations are updated as
The Bregman iteration imposes the constraints iteratively by adding the error back into the
constraints. Thus, introducing the Bregman iteration into the unconstrained formulation
[Equation (4)] forces its solution to converge to the solution of the constrained problem
[Equation (2)] for sufficiently small values of the parameters , , and . The data constraint in
Equation (5) leads to a sequence of solutions for which both the solution error norm and the data
fidelity term decrease monotonically. This formulation is more robust than equivalent
approximated unconstrained problems or continuation methods that impose the constraint iteratively
by slowly increasing the regularization parameters .
Note that, as ui and the auxiliary variables are independent of each other, Equation (4) can
now be split into several equations (one for each variable) that are solved sequentially, as
Since solution of ui only involves L2-norm functionals, it can be determined exactly by
differentiating the cost function and equating it to zero. The result is a linear system that
corresponds to a Gauss-Newton step
Note that Equation (7) is an analytical function, so an estimation of the step-size is not
required. This linear system constitutes a very large-scale problem, where K = NxN, being N the
number of pixels, yet it can be solved efficiently using a Krylov solver that involves only
mFT Fui lDxT Dxui lDyT Dyui lT2T T2ui gui ri
Here, we used the biconjugate gradient stabilized method with a threshold in the range of
102 to 104, where = 102 accelerated convergence (four to six iterations).
dxi, dyi, vi and wi are solved analytically using shrinkage formulas, which are thresholding
operations [Goldstein 2009, Wang 2008]
j u bikxj jDyuk1 bikyj ; j x; y
vi jvi2=O 0:
Test data: Simulation of different scenarios. Algorithms were evaluated using
simulations generated from data acquired from a 10-week old adult female Wistar rat weighing 300 g,
anaesthetized with isoflurane. Animals were handled according to the European Communities
Council Directive (86/609/EEC) and with the approval of the Animal Experimentation Ethics
Committee of Hospital General Universitario Gregorio Maran (ES 280790000087).
The CT subsystem used for data acquisition was ARGUS PET/CT (SEDECAL), a
conebeam micro-CT scanner based on a flat panel detector . The acquisition comprised 360
views of 512512 pixels (0.2 mm2 pixel size) covering 360 degrees with 32 images per
projection angle, at a source voltage of 45 kV.
Gated CT data were acquired using a high-dose protocol (four times the dose of a static
imaging protocol) . These high-dose projection data were arranged into four phases using
software-based retrospective gating  and reconstructed with an FDK-based algorithm . The
resulting images were selected as our target.
These data were used to simulate scenarios with different X-ray flux levels and number of
projections. Low-dose acquisitions were simulated by selecting a smaller field-of-view of
350x350 pixels (in order to reduce the computational cost), randomly taking 120 projections
or less for each phase, and adding Poisson noise by modelling the measurements fi as
independently distributed Poisson random variables:
where u(x,y,z) is the high-dose reconstruction, I0 is the number of photons emitted by the
xray source and M is the number of measured projections.
We simulated five different scenarios varying the number of photons emitted by the x-ray
source (I0 = 4.5 104, I0/2, and I0/4 photons) and the number of projections per phase (120, 80
and 60 projections). I0 was chosen so as to obtain a noise figure for the prior image similar to
that of the target (high-dose data). The number of projections per phase will depend on the
time resolution of the respiratory cycle, i.e. the number of respiratory phases. Simulations were
computed using the IRT code (J.A. Fessler, Image reconstruction toolbox (IRT), 2011, retrieved
Fig 1(A) shows four respiratory phases for the high-dose protocol and Fig 1(C) shows the
image difference between two phases. Most differences between phases are within the lung area
due to the respiratory movement and result in blurring of the prior image. Reconstruction of
the low-dose respiratory phases with FDK led to images with noise and streak artifacts (Fig 1
(B)). The prior was obtained by adding data from all phases and applying a Gaussian filter with
= 5 to reduce noise (Fig 1(D)).
Comparison of methods and analysis of images. The low-dose data were reconstructed
with three algorithms based on the PICCS approachthe unitary transform (L1-PICCS), the
gradient transform (TV-PICCS), and the wavelet transform (WT-PICCS)by varying
transform T2. The three variations of PICCS based on Equation (2) were solved using the Split
Bregman formulation [Equations (49)]. For the wavelet transform, we used the symmlet-8 base
It is necessary to select the reconstruction parameters in Equation (6), namely, k, , , ,
and . The iteration number k was chosen as the number of iterations that yielded the
minimum mean-square error with respect to the reference high-dose image. The regularization
parameter, , was selected following suggestions from previous studies [31, 33], which showed
that for sufficiently small values of ( 10 in our case), the problem converges to the same
solution, albeit at a higher iteration number. The regularization parameter has the opposite
effect. Low values ( 0.1) resulted in noisy images, as lower weight was given to the TV. For
large values ( 1), the problem converges to similar results, although at a different iteration
number. With regard to the regularization parameter , we found that low values ( 1) were
preferable because higher values impaired convergence. Given these considerations, we
empirically selected = 10, = 1, and = 0.1.
Finally, to select the value of the prior weight, , we analyzed the effect by computing the
two methods for = 0.2, = 0.5, and = 0.8. In order to use the same range of parameters for
all data sets, the algorithm normalizes the data item f as f/f/n, where n is the square root of
the number of pixels in the image, following the suggestions from [Tom Goldstein. Split
Bregman. Retrieved in 2009 from http://www.ece.rice.edu/~tag7/Tom_Goldstein/Split_Bregman.
Images were compared in terms of several quantitative parameters: 1) Contrast in lung and
bone areas measured as peak-to-valley (PV) on image profiles. 2) Noise, measured as the
coefficient of variation in three different oval 312-pixel ROIs inside soft tissue. 3) Contrast-to-noise
ratio, measured as the absolute difference between the value within a vessel and the value in a
lung-tissue ROI divided by the noise (Fig 2 left). 4) Reconstruction error, assessed as
meanFig 2. . Masks for quantitative analysis. Left: Mask used to measure contrast-to-noise ratio as the absolute difference between the yellow and green ROIs
divided by the noise measured in the blue ROI. Middle: Mask to compute MSE in the lung area. Right: Masks to compute MSE in the bone area.
Fig 3. Zoomed-in images of the low-dose protocol data for respiratory phase one. Reconstructions were performed with L1-PICCS, TV-PICCS, and
WT-PICCS for equal to 0.2, 0.5, and 0.8, for 120 projections and X-ray flux corresponding to a number of photons I0 = 4.5 104. Arrows point at locations
where differences in texture are more noticeable.
square error with respect to the target image (high-dose protocol) for both the lung area and
bone, measured in the masks shown in Fig 2. 5) Movement compensation was assessed by
drawing two profiles on moving areas, one across the vertebrae and another crossing a vessel
inside the lung, and comparing them with the same profiles in the target image. In addition,
image texture was evaluated by visual inspection.
Fig 3 shows zoomed-in images of the best result for each CS method, for different values of the
prior weight . For small values of , the prior term has a small influence and all CS methods
converge to a solution similar to that provided by spatial TV, with a noticeable patchy pattern
( = 0.2, first column of Fig 3). Larger values increase weight to the prior and leads to
differences in the image texture depending on the sparsity transform used for this term. For = 0.5
or 0.8, some differences in image texture are visible: L1-PICCS shows salt-and-pepper artifacts
Fig 4. Images of several low-dose protocol data for respiratory phase one reconstructed with FDK, TV-PICCS, and WT-PICCS ( equal to 0.8). Each
column represents a different scenario: 120 projections and flux corresponding to a maximum number of photons I0 (I0 = 4.5 104), 60 projections and number
of photons I0, and 120 projections and number of photons I0/4 (from left to right).
and TV-PICCS shows a patchy-like pattern (arrows in column 3 of Fig 3), while WT-PICCS
shows a more natural texture.
Analysis of the influence of X-ray flux and number of projections
Fig 4 shows the differences between FDK, TV-PICCS, and WT-PICCS reconstructions when
decreasing the X-ray flux and the number of projections.
Fig 5. Zoomed-in images of low-dose protocol data for respiratory phase one, reconstructed with FDK, TV-PICCS, and WT-PICCS. Columns
represent the different scenarios: 120 projections and dose corresponding to a maximum number of photons I0 (I0 = 4.5 104), 60 projections and number of
photons I0, and 120 projections and number of photons I0/4 (from left to right). TV-PICCS and WT-PICCS were obtained with = 0.8.
Fig 5 shows zooms of Fig 4 to better depict differences between TV-PICCS and WT-PICCS
when varying the source flux and number of projections. When decreasing the X-ray flux
(third column in Fig 5), both algorithms converged in fewer iterations towards more blurred
images (the loss of resolution can also be seen in the profile plotted in Fig 6). As the iteration
number increases, the sparsity is imposed in the different domains. In the case of TV-PICCS,
the algorithm leads to patchy artifacts for all noise levels tested, although with high noise the
patchy artifact is less evident. WT- PICCS is more robust against noise level and maintains a
more natural texture for all three scenarios. Varying the number of projections has a similar
effect on texture: TV-PICCS presents patchy artifacts while WT-PICCS maintains a more
natural texture. As the number of projections decreases, the missing data produce streaks that are
not removed by any of the algorithms due to the coherent nature of this artifact, which cannot
be removed by the TV term common to both algorithms.
Fig 6 illustrates the effect of the different number of projections and X-ray flux on an image
profile drawn over lung and bone areas. Lowering either the number of projections or the
Xray flux reduces the recovered contrast in lung and bone. Decreasing the number of projections
led to a 40% reduction in PV ratio in the lung profile with respect to the target, with no
differences between the methods. For bone, decreasing the number of projections led to a 40%
reduction in PV ratio using WT-PICCS and a 27% reduction using TV-PICCS.
Fig 7 shows the image noise for all scenarios. Increasing the X-ray flux leads to less image
noise for both TV-PICCS and WT-PICCS. Decreasing the number of projections does not
have a significant influence on the image noise.
Fig 8 shows MSE for bone and soft tissue in different scenarios. For bone tissue, TV-PICCS
and WT-PICCS showed a MSE reduction of 83% for the lowest the X-ray flux and of 67% for
Fig 6. Influence of the number of projections and X-ray flux on image profiles. Top: Profile (yellow line) drawn over heart and lung areas (left figure) and
profile drawn over a bone area (right figure), overimposed on the high-dose protocol image. Middle and bottom: Normalized profiles for reference high-dose
FDK (target) and for the low-dose protocol reconstructed with TV-PICCS (TV) and WT-PICCS (WT) for different number of projections (middle) and different
X-ray flux values (bottom).
60 projections, with respect to FDK reconstruction. There are no large differences between
TV-PICCS and WT-PICCS, although WT-PICCS presented a slightly lower MSE. For lung
tissue both MSE and PICCS showed an MSE sixty times lower than that of FDK.
Fig 9 shows the contrast-to-noise ratio of nodules in the lung. Both TV-PICCS and
WT-PICCS lead to a ten-fold increase in CNR with respect to FDK, where TV-PICCS has
slightly higher CNR than WT-PICCS.
Analysis of the respiratory movement compensation
Fig 10 shows profiles along a line containing lung tissue and vessels (Fig 10, left) and along
bone tissue (Fig 10, right) for the reference high-dose FDK and WT-PICCS reconstructions of
Fig 7. Coefficient of variation measured in images reconstructed with TV-PICCS and WT-PICCS for
the different scenarios. Plots show mean and standard deviation of the coefficient of variation measured in
three different 312-pixel ROI inside soft-tissues. Left panel shows different X-ray flux values for 120
projections; right panel represents different number of projections for a flux corresponding to a number of
photons I0 = 4.5 104.
respiratory phases 1 and 3 using 120 projections and I0 = 4.5 104. The profiles reveal the
existence of respiratory movement for the two respiratory phases; measuring the separation
between profiles for frames 1 and 3 provides an estimate of motion of 0.5 mm. However, profiles
for WT-PICCS fit the reference case well, with an error of less than 6 m for most points in
Fig 8. MSE with respect to the reference high-dose image for FDK, TV-PICCS and WT-PICCS for the
different scenarios. Plots show MSE in bone tissue (top) and lung tissue (bottom) in the ROIs defined in Fig
2. Left panel shows different X-ray flux values for 120 projections; right panel represents different number of
projections for an X-ray flux corresponding to a number of photons I0 = 4.5 104.
Fig 9. . Contrast-to-noise ratio measured in lung tissue in images reconstructed with FDK, TV-PICCS
and WT-PICCS for the different scenarios. Left panel shows in the results for different X-ray flux values for
120 projections; right panel represents different number of projections for an X-ray flux corresponding to a
number of photons I0 = 4.5 104.
Fig 10. Respiratory artifact analysis. Profiles along the yellow lines in soft tissue (left) and bone tissue (right) for reference high-dose FDK (target) and for
the low-dose protocol reconstructed with WT-PICCS corresponding to respiratory phases 1 and 3. The analysis shows that the reconstruction can follow the
movement of the lung and vessels in the two respiratory phases.
Fig 11. Images of respiratory phases one and three corresponding to the reference protocol for static studies reconstructed with FDK and
Fig 11 shows phases one and three, reconstructed with FDK and WT-PICCS. While FDK is
highly affected by incomplete projections and noise, which hinders details and differences
between phases, WT-PICCS is able to remove streak and motion artifacts, recovering the
differences between phases.
We evaluated the suitability of different sparsity transforms (unitary, gradient, and wavelet)
within the PICCS formulation for reducing dose in CT respiratory gating. Performance was
assessed in different scenarios, corresponding to different X-ray flux levels and number of
projections, and for different weights of the prior penalty term. Overall, our results show that the
selection of the sparsity transform for the prior term does not affect spatial resolution, temporal
resolution or noise performance, but has an influence on the final image texture: the wavelet
transform showed a more natural pattern than the gradient and unitary transforms. While
decreasing the X-ray flux leads to higher image noise, decreasing the number of projections does
not have a large influence on the image noise but leads to more streak artifacts, although we
did not detect significant differences between TV-PICCS and WT-PICCS. Decreasing the
number of projections or the X-ray flux lowered the recovered contrast in bone and lung.
While TV-PICCS and WT-PICCS presented similar profiles across soft tissue and lung,
WT-PICCS led to more blurred edges for small structures than TV-PICCS probably because of
the patchy-like pattern of TV regularization. Both TV-PICCS and WT-PICCS greatly
For the three reconstruction methods the results depended on the prior weight, . For low
, all the methods converged to a solution very similar to that of the spatial TV method, which
presented a patchy pattern. Increasing the prior weight reduced this pattern to some extent for
all the methods. The best results for TV-PICCS and WT-PICCS were obtained using a large
prior weight ( = 0.8), while for L1-PICCS the best results appeared at an intermediate weight
( = 0.5), as a large weight was prone to pepper and salt artifacts. Overall, the WT-PICCS
proved to be more robust against the value, regarding the production of artifacts. These
results are consistent with prior reports on the performance of the TV-PICCS method, which
found optimal values in the range of 0.5 to 0.8, with no decrease in spatial resolution and a
texture similar to that of the prior [17, 18, 20, 23, 37].
Few studies compared sparsity transforms. In this work, we compared the gradient
transform, generally used for CT, with the pixel domain and the symmlet-8 wavelet transform.
Other transforms such as overcomplete wavelet transform, shearlets, curvelets, and dictionary
learning based sparse representation have been proposed [28, 38, 39] but were not included in
There are some limitations in this study. Although we have evaluated the methods under
several scenarios, varying the x-ray source intensity and the number of projections per phase,
the effect of some variables such as image pixel size has not been assessed. In theory, smaller
pixel size for lower order basis functions would result in similar results as larger pixels with
higher order basis functions. Thus, the smoother results obtained by using WT or another
transform instead of TV may be less relevant for smaller pixel size. However, decreasing the
pixel size demands unnecessary higher computational cost and in  it was found to lead to
poor results because the system becomes more underdetermined.
With regard to the reconstruction algorithm, we studied the influence of the sparsity
parameter and verified that varying the rest of parameters (, , and ) within a certain range did
not noticeably change the results. The most important parameter is ,which weights the data
constraint and thus affects the convergence speed. Higher values of speed up convergence,
although they must remain sufficiently small to guarantee convergence . The values are
usually chosen empirically . In fact, it has been shown that results for the Split Bregman
method are independent of the actual value of , provided that it is sufficiently small and that
the number of iterations is large . Thus, there is no need to carefully optimize the weighting
parameters, as opposed to unconstrained optimization methods, where regularization
parameters have to be cautiously selected (for example with the L-curve or U-curve method [40, 41].
This is an additional advantage of the Split Bregman formulation. However, one still has to
choose the number of iterations. In this work, since the high-dose data are known, we selected
the iteration number that minimized the solution error taking the high-dose images as a
reference, but further work is required to select a suitable number of iterations when the target
image is not available.
Further improvements could be made to the proposed method. With regard to the Split
Bregman formulation, in  the authors included statistical data modeling in a similar
alternating approach, thus improving convergence. As for the prior image, although we used one
prior image based on the average of data for all phase bins, other priors, such as a running
average, could be used . Furthermore, in this study, compressed sensing was modeled as a
convex L1-norm problem. An alternative could be to propose an equivalent non-convex L0-norm
problem and solve it using memetic algorithms, an improved version of evolutionary
algorithms, which exploit the available information on the problem [42, 43].
In conclusion, we compared different methodologies for the reconstruction of low-dose CT
data with respiratory gating based on the PICCS approach using different transforms for the
prior term. Our results show that, although the gradient transform is widely used, the wavelet
transform could reduce the formation of patchy cartoon-like artifacts.
Conceived and designed the experiments: MA JA. Performed the experiments: MA AS.
Analyzed the data: JA MA. Contributed reagents/materials/analysis tools: JA MA JJV MD. Wrote
the paper: JA MA MD JJV AS.
18. Tang J, Hsieh J, Chen GH. Temporal resolution improvement in cardiac CT using PICCS (TRI-PICCS):
performance studies. Med Phys., 37(8):437788. Med Phys. 2010; 37(8): 437788. PMID: 20879597
1. Chavarras C , Vaquero JJ , Sisniega A , Rodrguez-Ruano A , Soto-Montenegro ML , Garca-Barreno P , et al. Extraction of the respiratory signal from small-animal CT projections for a retrospective gating method . Phys Med Biol . 2008 ; 53 ( 17 ): 4683 - 4695 . doi: 10.1088/ 0031 - 9155 /53/17/015 PMID: 18695300
2. Sawall S , Bergner F , Lapp R , Mronz M , Karolczak M , Hess A , et al. Low-dose cardio-respiratory phasecorrelated cone-beam micro-CT of small animals . Med Phys . 2011 ; 38 ( 3 ): 1416 - 24 . PMID: 21520853
3. Leng S , Zambelli J , Tolakanahalli R , Nett B , Munro P , Star-Lack J , et al. Streaking artifacts reduction in four-dimensional cone-beam computed tomography . Med Phys . 2008 ; 35 ( 10 ): 4649 - 4659 . PMID: 18975711
4. McKinnon GC , Bates RHT . Towards imaging the beating heart usefully with a conventional CT scanner . IEEE Trans Biomed Eng . 1981 ; 28 ( 2 ): 123 - 127 . PMID: 7287019
5. Zhang Q , Pevsner A , Hertanto A , Hu Y , Rosenzweig K , Ling C , et al. A patient-specific respiratory model of anatomical motion for radiation treatment planning . Med Phys . 2007 ; 34 ( 12 ): 4772 - 81 . PMID: 18196805
6. Zhang Q , Hu Y , Liu F , Goodman K , Rosenzweig K , Mageras G. Correction of motion artifacts in conebeam CT using a patient-specific respiratory motion model . Med Phys . 2010 ; 37 ( 6 ): 2901 - 9 . PMID: 20632601
7. Cands EJ , Romberg J. Practical signal recovery from random projections . Wavelet Applications in Signal and Image Processing XI, Proc. SPIE Conf . 2005 ; 5914 .
8. Candes EJ , Romberg J , Tao T. Robust uncertainty principles: exact signal reconstruction from highly incomplete frequency information . IEEE Trans Inf Theory . 2006 ; 52 : 489 - 509 .
9. Cands E , Romberg J. Sparsity and incoherence in compressive sampling . Inverse Problems . 2007 ; 23 ( 3 ): 969 .
10. Bruckstein AM , Donoho DL , Elad M. From Sparse Solutions of Systems of Equations to Sparse Modeling of Signals and Images SIAM J Appl Math . 2009 ; 51 ( 1 ): 34 - 81 .
11. Pan X , Sidky EY , Vannier M. Why do commercial CT scanners still employ traditional, filtered back-projection for image reconstruction . Inverse Probl . 2009 ; 25 ( 12 ): 123009 . PMID: 20376330
12. Song J , Liu QH , Johnson GA , Badea C. Sparseness prior based iterative image reconstruction for retrospectively gated cardiac micro-CT . Med Phys . 2007 ; 34 ( 11 ): 4476 - 83 . PMID: 18072512
13. Yu H , Wang G . Compressed sensing based interior tomography . Phys Med Biol . 2009 ; 54 ( 9 ): 2791 - 2805 . doi: 10.1088/ 0031 - 9155 /54/9/014 PMID: 19369711
14. Ritschl L , Bergner F , Fleischmann C , Kachelriess M. Improved total variation-based CT image reconstruction applied to clinical data . Phys Med Biol . 2011 ; 56 ( 6 ): 1545 - 61 . doi: 10.1088/ 0031 - 9155 /56/6/ 003 PMID: 21325707
15. Tang J , Nett BE , Chen GH . Performance comparison between total variation (TV)-based compressed sensing and statistical iterative reconstruction algorithms . Phys Med Biol . 2009 ; 54 ( 19 ): 5781 - 804 . doi: 10.1088/ 0031 - 9155 /54/19/008 PMID: 19741274
16. Chen GH , Tang J , Leng S. Prior image constrained compressed sensing (PICCS): A method to accurately reconstruct dynamic CT images from highly undersampled projection data sets . Med Phys . 2008 ; 35 ( 2 ): 660 - 3 . PMID: 18383687
17. Nett BE , Brauweiler R , Kalender W , Rowley H , Chen GH . Perfusion measurements by micro-CT using prior image constrained compressed sensing (PICCS): initial phantom results . Phys Med Biol . 2010 ; 55 ( 8 ): 2333 - 50 . doi: 10.1088/ 0031 - 9155 /55/8/014 PMID: 20360635
19. Leng S , Tang J , Zambelli J , Nett B , Tolakanahalli R , Chen GH . High temporal resolution and streak-free four-dimensional cone-beam computed tomography . Phys Med Biol . 2008 ; 53 ( 20 ): 5653 - 73 . doi: 10. 1088/ 0031 - 9155 /53/20/006 PMID: 18812650
20. Qi Z , Chen GH . Performance studies of four-dimensional cone beam computed tomography . Phys Med Biol . 2011 ; 56 : 6709 - 6721 . doi: 10.1088/ 0031 - 9155 /56/20/013 PMID: 21965275
21. Lee H , Xing L , Davidi R , Li R , Qian J , Lee R. Improved compressed sensing-based cone-beam CT reconstruction using adaptive prior image constraints . Phys Med Biol . 2012 ; 57 : 2287 - 2307 . doi: 10. 1088/ 0031 - 9155 /57/8/2287 PMID: 22460008
22. Stayman JW , Zbijewski W , Otake Y , Uneri A , Schafer S , Lee J , et al. Penalized-likelihood reconstruction for sparse data acquisitions with unregistered prior images and compressed sensing penalties . Proc Physics of Medical Imaging, SPIE Medical Imaging . 2011 ; 7961 : 79611L.
23. Ramirez-Giraldo JC , Trzasko J , Leng S , Yu L , Manduca A , McCollough CH . Nonconvex prior image constrained compressed sensing (NCPICCS): Theory and simulations on perfusion CT . 2157 ( 2011 ). Med Phys . 2011 ; 38 ( 4 ): 2157 - 67 . PMID: 21626949
24. Sidky EY , Pan X. Image reconstruction in circular cone-beam computed tomography by constrained, total-variation minimization . Phys Med Biol . 2008 ; 53 ( 17 ): 4777 - 4807 . doi: 10.1088/ 0031 - 9155 /53/17/ 021 PMID: 18701771
25. Goldstein T , Osher S. The Split Bregman Method for L1 Regularized Problems . SIAM Journal on Imaging Sciences . 2009 ; 2 ( 2 ): 323 - 343 .
26. Montesinos P , Abascal JFPJ , Cuss L , Vaquero JJ , Desco M. Application of the compressed sensing technique to self-gated cardiac cine sequences in small animals . Magnetic Resonance in Medicine . 2013 .
27. Ramani S , Fessler JA. A Splitting-Based Iterative Algorithm for Accelerated Statistical X-Ray CT Reconstruction. IEEE Trans Med Imaging . 2012 ; 31 ( 3 ): 677 - 688 . doi: 10.1109/TMI. 2011 .2175233 PMID: 22084046
28. Vandeghinste B , Goossens B , Van Holen R , Vanhove C , Piurica A , Vandenberghe S , et al. Iterative CT reconstruction using shearlet-based regularization . SPIE Medical Imaging . 2012 .
29. Dong B , Li J , Shen Z. X-Ray CT Image Reconstruction via Wavelet Frame Based Regularization and Radon Domain Inpainting . Journal of Scientific Computing . 2012 : 1 - 17 .
30. Abascal J , Sisniega A , Chavarras C , Vaquero J , Desco M , Abella M. Investigation of different Compressed Sensing Approaches for Respiratory Gating in Small Animal CT . IEEE Nuclear Science Symposium and Medical Imaging Conference Record . 2012 : 3344 - 3346 .
31. Abascal J , Chamorro-Servent J , Aguirre J , Arridge S , Correia T , Ripoll J , et al. Fluorescence diffuse optical tomography using the split Bregman method . Med Phys . 2011 ; 38 ( 11 ): 6275 - 6284 . doi: 10.1118/ 1.3656063 PMID: 22047393
32. Setzer S , Steidl G , Teuber T. Deblurring Poissonian images by split Bregman techniques . J. Vis. Comun. Image Represent . 2010 ; 21 ( 3 ): 193 - 199 .
33. Osher S , Burger M , Goldfarb D , Xu J , Yin W. An iterative regularization method for total variation-based image restoration . SIAM J Multiscale Model Simul . 2005 ; 4 ( 2 ): 460 - 489 .
34. Vaquero JJ , Redondo S , Lage E , Abella M , Sisniega A , Tapias G , et al. Assessment of a New High-Performance Small- Animal X-ray Tomograph. IEEE Trans Nucl Sci . 2008 ; 55 ( 3 ): 898 - 905 .
35. Abella M , Vaquero JJ , Sisniega A , Pascau J , Udas A , Garca V , et al. Software Architecture for MultiBed FDK-based Reconstruction in X-ray CT Scanners . Comput Methods Programs Biomed . 2012 ; 107 ( 2 ): 218 - 32 . doi: 10.1016/j.cmpb. 2011 . 06.008 PMID: 21908068
36. Buckheit JB , Chen S , Donoho DL , Johnstone IM , Scargle JD . WaveLab . Reference Manual . ftp://playfair.stanford.edu/pub/wavelab/WaveLabRef.ps. 1995 .
37. Lauzier PT , Tang J , Chen GH . Time-resolved cardiac interventional cone-beam CT reconstruction from fully truncated projections using the prior image constrained compressed sensing (PICCS) algorithm . Phys Med Biol . 2012 ; 57 ( 9 ): 2461 - 76 . doi: 10.1088/ 0031 - 9155 /57/9/2461 PMID: 22481501
38. Starck JL , Elad M , Donoho DL . Image decomposition via the combination of sparse representations and a variational approach . IEEE Trans Image Process 2005 ; 14 ( 10 ): 1570 - 82 . PMID: 16238062
39. Xu Q , Yu H , Mou X , Zhang L , Hsieh J , Wang G. Low-dose X-ray CT reconstruction via dictionary learning . IEEE Trans Med Imaging . 2012 ; 31 ( 9 ): 1682 - 97 . doi: 10.1109/TMI. 2012 .2195669 PMID: 22542666
40. Hansen PC , O'Leary DP . The use of the L-curve in the regularization of discrete ill-posed problems . SIAM Journal on Scientific Computing . 1993 ; 14 ( 6 ): 1487 - 1503 .
41. Chamorro-Servent J , Aguirre J , Ripoll J , Vaquero JJ , Desco M. Feasibility of U-curve method to select the regularization parameter for fluorescence diffuse optical tomography in phantom and small animal studies . Optics Express . 2011 ; 19 ( 12 ): 11490 - 11506 . doi: 10.1364/OE.19.011490 PMID: 21716381
42. Moscato P , Cotta C. A Gentle Introduction to Memetic Algorithms , in Handbook of Metaheuristics. International Series in Operations Research & Management Science , Glover F. and Kochenberger G. , Editors., Kluwer Academic Publishers: 2003
43. Yanga S , Chenga K , Wangb M , Xiea D , Jiaoa L. High resolution range-reflectivity estimation of radar targets via compressive sampling and Memetic Algorithm . Information Sciences . 2013 ; 252 ( 10 ): 144 - 156 .