SPECT Reconstruction with a Non-linear Transformed Attenuation Prototype Sven Barendt, Jan Modersitzki Institut of Mathematics and Image Computing, University of Lübeck sven.barendt@mic.uni-luebeck.de Abstract. This work deals with the single photon emission computed tomography (SPECT) reconstruction process. As a SPECT measure- ment also depends on unknown attenuation properties of the tissue, such a process is challenging. Furthermore, the given attenuation may not be a good approximation to the true attenuation field. Reasons are reposi- tioning or movement of the patient such as relaxation during scan time or even breathing. We propose a novel model for an attenuation correction in SPECT reconstruction, which is a natural extension of an idea of Nat- terer in that way, as the linear transformation of a so-called attenuation prototype is enhanced to an arbitrary transformation. We present nu- merical results for a non-linear spline transformation model which clearly indicate the superiority of the proposed reconstruction model compared to the case of no motion correction and the correction with a linear transformation model. 1 Introduction Single photon emission computed tomography (SPECT) is a nuclear medicine imaging technique which can provide in vivo 3D functional information of tis- sue. More precisely, functional information corresponds to the density of an administered radio-pharmaceutical which has to be reconstructed from a set of projections also known as sinogram. The reconstruction process is challenging as the sinogram also depends on unknown attenuation properties of the tissue. A commonly used simplification is to assume that the attenuation field is given and to solve only for the density [1, 2]. Practically, the attenuation informa- tion is supplied by an extra measurement such as a computed tomography (CT) scan. Drawbacks of such an approach are that it requires a costly additional scan which adds stress to the patient. Furthermore, the given attenuation may not be a good approximation to the true attenuation field. Reasons are repositioning or movement of the patient such as relaxation during scan time or even breathing or heart beat. Because of these limitations, recent approaches aim to simultaneously recon- struct the density and attenuation [3, 4, 5, 6, 7, 8, 9, 2]. However, a fundamental problem is that the combined reconstruction is generally ill-posed [4, 9]. In order to circumvent ill-posedness of the combined reconstruction problem, Chang used a simple model for the unknown attenuation field, which is essentially based on a Non-linear SPECT Reconstruction 415 piecewise constant function [10]. However, the attenuation model is inadequate to cover complex patient anatomy. Dicken addressed ill-posedness by introducing non-linear Tikhonov regularization for both, density and attenuation [9]. How- ever, limited success has been reported in the literature [2]. Natterer phrased so-called consistency conditions on the range of the projection operator [4], which then constrain the set of feasible attenuations. This approach was then used in combination with a linearly transformed prototype for attenuation [6] and a piecewise constant attenuation field [2]. This paper pursues the work of Natterer [6] as it removes the limitations to an affine linear transformation model and enables essentially arbitrary transfor- mations. In addition, we present a continuous mathematical framework for the combined reconstruction model and follow the so-called discretize then optimize paradigm to compute a numerical solution. This paper is organized a follows. A novel reconstruction model is presented in Section 2 and two numerical examples which clearly demonstrate superiority of the non-linear approach are shown in Section 3. Finally, in Section 4, results are discussed and an outlook is given. 2 Materials and Methods In this section we introduce a novel model for the simultaneous reconstruction of the tracer density f and the attenuation field µ, where we assume that the unknowns f, µ ∈ L2 (R2 , R) are compactly supported on a domain Ω ⊂ R2 . As common, the forward projection process is modeled by g = P[f, µ] with the projection operator P (attenuated ray transform) and the measurement denoted by g, e.g., [11]. Our reconstruction model is to find a minimizer (f, µ) of the energy functional J [f, µ] = ∥g − P[f, µ]∥L2 + αR[µ] (1) where α denotes a positive regularization parameter and R a regularizer, which is discussed in detail later. Natterer used a similar approach with the constraint µ(x) = ν(ℓ(x)), where ℓ denotes an affine linear spatial transformation and ν denotes a given prototype of the expected attenuation [6]. Thus, the set of feasible µ consists only of affine linear transformations of the prototype. Moreover, regularization is skipped and consistency conditions are used to determine the affine linear transformation. The new idea is to extend the class of feasible transformation and to replace the consistency assumption by regularization. To simplify our presentation, we restrict our model to free-form transformations, i.e. the transformation can be expanded in terms of basis functions Bi (e.g. monomials of degree less equal to one and B-splines of a certain order) and coefficients ci . Moreover, we also limit our discussion to plain Tikhonov regularization of the coefficient vector. However, it is noticeable that our model is far more general and allows for non-parametric transformations and arbitrary regularization. With the above 416 Barendt & Modersitzki specification, we have { ( ∑ ) } R[µ] = ∥c∥22 where µ ∈ M := µ(c; x) = ν x + i ci Bi (x) , c ∈ Rn (2) The mathematical model is to minimize J 1 subject to µ ∈ M 2. As an ana- lytical solution to this problem is unknown, we employ numerical optimization techniques. In particular, we follow the discretize then optimize paradigm and a straightforward cell-centered discretization of the spatial domain. To eliminate the constraint, we use a reduction approach and the unknown of the reduced model are the density values f = (fi )m i=1 on cell-centers and the coefficients c ∈ Rn . The optimization is performed with a generalized Gauss-Newton scheme and an Armijo Linesearch, e.g. [12] for details. 3 Results The scheme has been implemented using MATLAB R 2010a. In order to have a gold standard for comparison and validation, we use synthetic data generated using the XCAT phantom [13]. To be more precise, from a 4D simulation of a human torso 2D slices are taken. That is, an average activity image f0 and two attenuation images µ0 and ν at different time steps are chosen, such that µ0 and ν differ with respect to respiratory motion. The activity and the attenuation images are calculated by the XCAT software based on a 140 KeV photon energy. The two time steps related to the attenuation images assume a respiration period of 5 seconds and are simulated at 0 seconds and 2.5 seconds. This corresponds to a 0 % and 96 % inhale, according to the default respiration of the XCAT phantom. The maximal motion of the chest is parameterizable and choosen to be 1.2 cm. According to the XCAT phantom a movement of 1.2 cm simulates a normal breathing activity. For the sake of testing (a) f0 (b) ν (c) µ0 (d) g Fig. 1. One of the two test cases are exemplarily shown above. The averaged density of radioactivity over a period of 5 seconds is shown on the left (f0 ). The next two images (ν and µ0 ) show the attenuation at two different time steps. As respiration movement is taken into account, these images differ in a non-linear movement. The image on the right (g) is the simulated SPECT projection (sinogram) which is formally the application of the attenuated ray transform on f0 and µ0 . The task is to reconstruct f0 and µ0 , given the sinogram g and ν. Non-linear SPECT Reconstruction 417 Table 1. Relative error between the gold standard (Fig. 1) without any motion cor- rection, linear, and the proposed non-linear motion correction. err(f ) = ∥f∥f−f 0 ∥2 0 ∥2 err(µ) = ∥µ−µ 0 ∥2 ∥µ0 ∥2 chest movement no correction linear non-linear no correction linear non-linear 1.2 cm 48.43% 31.16% 18.84% 29.13% 32.36% 17.12% 2.4 cm 81.48% 36.50% 18.19% 39.65% 47.06% 23.31% of the proposed motion correction in SPECT reconstruction, the comparatively larger non-linear movement of the chest with a 2.4 cm a.-p. movement is choosen as a second test case. In Figure 1 the latter test case is shown exemplarily. In the following the reconstruction results for the test case presented in Figure 1 are illustrated in Figure 2 in comparison with the gold standard provided by the XCAT phantom. In Table 1 a relative error between the gold standard and different reconstructions is shown. (a) no correction (b) linear (c) non-linear Fig. 2. The absolute value of the difference between the gold standard and different reconstructions. Rows (a) and (b) are related to the density of radioactivity f and the attenuation µ, respectively. For instance, the image in the lower left corner illustrates the absolute difference of the attenuation prototype ν and the gold standard µ0 . 418 Barendt & Modersitzki 4 Discussion The proposed model for a motion correction of SPECT reconstruction is a nat- ural extension of the idea of Natterer as the linear transformation model for the attenuation prototype ν is enhanced to an arbitrary transformation. We present results for a non-linear spline transformation model which clearly indi- cate the superiority of the proposed reconstruction model compared to the case of no motion correction and the correction with a linear transformation model (Tab. 1). Due to the nature of the addressed SPECT reconstruction problem, the type and weight of the regularization is a crucial issue. In particular, an improved regularization model is under current investigation. As the current implementation deals with 2D reconstructions of the density of radioactivity and the attenuation, future work involves the implementation of a 3D reconstruction. Furthermore it is planned to explore different transforma- tion models, as well as different regularizations. Because of the promising and superior results of the previous section it would be interesting, to validate the proposed reconstruction model in a clinical setting. References 1. Bronnikov AV. Reconstruction of attenuation map using discrete consistency con- ditions. IEEE Trans Med Imaging. 2000;19(5):451–62. 2. Mennessier C, Noo F, Clackdoyle R, et al. Attenuation correction in SPECT using consistency conditions for the exponential ray transform. Phys Med Biol. 1999;44(10):2483–2510. 3. Censor Y, Gustafson DE, Lent A, et al. A new approach to the emission comput- erized tomography problem: simultaneous calculation of attenuation and activity coefficients. IEEE Trans Nucl Sci. 1979;26(2):2775–9. 4. Natterer F. The identification problem in emission tomography. Lect Notes com- puter Sci. 1981;8. 5. Manglos SH, Young TM. Constrained IntraSPECT reconstruction from SPECT projections. In: Proc Nucl Sci Symp Med Imaging Conf. vol. 3; 1993. p. 1605–9. 6. Natterer F. Determination of tissue attenuation in emission tomography of opti- cally dense media. Inverse Probl. 1993;9:731–6. 7. Bronnikov AV. Approximate reconstruction of attenuation map in SPECT imag- ing. IEEE Trans Nucl Sci. 1995;42(5):1483–8. 8. Welch A, Clack R, Natterer F, et al. Toward accurate attenuation correc- tion in SPECT without transmission measurements. IEEE Trans Med Imaging. 1997;16:532–41. 9. Dicken V. A new approach towards simultaneous activity and attenuation recon- struction in emission tomography. Inverse Probl. 1999;15(4):931–60. 10. Chang LT. A method for attenuation correction in radionuclide computed tomog- raphy. IEEE Trans Nucl Sci. 1978;25:638–43. 11. Natterer F, Wübbeling F. Mathematical Methods in Image Reconstruction. SIAM; 2001. 12. Nocedal J, Wright SJ. Numerical Optimization. Springer; 1999. 13. Segars WP, Sturgeon GM, Mendonca S, et al. 4D XCAT phantom for multimodal- ity imaging research. J Med Phys. 2010;37:4902–15.