<%BANNER%>

Correction of the Respiratory Motion of the Heart from 4D Myocardial Images


PAGE 1

CORRECTION OF THE RESPIRATORY MOTION OF THE HEART FROM 4D MYOCARDIAL IMAGES By JING SUN A THESIS PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORI DA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF SCIENCE UNIVERSITY OF FLORIDA 2006

PAGE 2

Copyright 2006 by Jing Sun

PAGE 3

Dedicated to my belove d parents and my husband

PAGE 4

iv ACKNOWLEDGMENTS First of all, I would like to express my deep indebtedness and my sincere gratitude to my advisor, David Gilland, for the va luable guidance, support and continuous encouragement he has given me during the course of this work. His insight in this field and his thoughtful personality have made my study here wo rthwhile and unforgotten. I am also equally indebted to the members of my committee, Wesley Bolch and David Hintenlang, for their helpful advice, suggestio ns and encouragements Special thanks are also extended to Bernard Mair, the external committee member, for his suggestions and attendance. It is a great honor to have them on my committee. I am grateful to Willams Paul Segars, for his tremendous help with the NCAT phantom. I would like to thank Dr. Zixiong Cao for his patience and help to set up this research idea. I am also grateful to Uday Tipnis, Keelan Seabolts, and Jason Parker, who have inspired and assisted in my research study. Finally I would like to extend my sincere gratitude to my parents for their moral support from far away in my home country, China. I would like to thank my husband Dapeng for his love, understanding and moral support.

PAGE 5

v TABLE OF CONTENTS page ACKNOWLEDGMENTS.................................................................................................iv LIST OF TABLES............................................................................................................vii LIST OF FIGURES.........................................................................................................viii ABSTRACT....................................................................................................................... ..x CHAPTER 1 INTRODUCTION........................................................................................................1 1.1 An Overview...........................................................................................................1 1.2 Literature Review...................................................................................................2 1.3 Our Method.............................................................................................................2 2 RESPIRATORY MOTION COMPENSATION, IMAGE RECONSTRUCTION, AND CARDIAC WALL MOTION ESTIMATION METHODS...............................4 2.1 Introduction.............................................................................................................4 2.2 NCAT Phantom Simulation....................................................................................4 2.2.1 Myocardial Phantom Simulation..................................................................5 2.2.2 Cardiac Wall Motion Simulation..................................................................7 2.2.3 Projection Simulation...................................................................................8 2.3 Axial Center-of-Mass Motion (aCOM) Estimation and Respiratory Motion Compensation Method............................................................................................10 2.3.1 Axial Center-of-Mass (aCOM) Calculation...............................................10 2.3.2 Respiratory Motion Compensation.............................................................12 2.4 Image Reconstruction and Cardiac Wall Motion Estimation (RM) Method........13 2.4.1 RM Algorithm............................................................................................13 2.4.2 Objective Function.....................................................................................14 2.4.2.1 Negative log likelihood term............................................................15 2.4.2.2 Image matc hing term........................................................................15 2.4.2.3 Strain energy term............................................................................16 2.4.3 Reconstruction and Motion Estimation......................................................18 2.5 Image Quality and Cardiac Wall Motion Accuracy Evaluation Methods............18 2.5.1 Image Quality Evaluation...........................................................................18 2.5.2 Cardiac Wall Motion Evaluation................................................................19

PAGE 6

vi 3 EXPERIMENTAL RESULTS AND DISCUSSION.................................................20 3.1 Simulated Phantoms, True Mo tion Field, and Projections...................................20 3.1.1 Cardiac Gated NCAT Phantom..................................................................20 3.1.1.1 Reference phantom with beating heart motion only........................20 3.1.1.2 Phantom with both beating hear t motion and respiratory motion....20 3.1.2 True Motion Field of the NCAT Phantom.................................................21 3.1.3 Projection Datasets.....................................................................................22 3.2 aCOM Calculation and Respirat ory Motion Compensation Method...................23 3.2.1 aCOM Calculation......................................................................................23 3.2.2 Respiratory Motion Compensation.............................................................25 3.3 Image Reconstruction, Cardiac Wall Motion Estimation, and SSE and PME Evaluations..............................................................................................................26 3.3.1 Image Reconstruction and Wall Motion Estimation..................................27 3.3.2 SSE and PME Evaluation Method.............................................................30 4 CONCLUSIONS AND FUTURE WORK.................................................................32 4.1 Summary and Conclusions...................................................................................32 4.2 Future Work..........................................................................................................33 LIST OF REFERENCES...................................................................................................34 BIOGRAPHICAL SKETCH.............................................................................................37

PAGE 7

vii LIST OF TABLES Table page 3-1. The SSE results of OHM, ORMA, an d ORMC images quality evaluation.............30 3-2. The PME results of OHM, ORMA, a nd ORMC motion estimation accuracy evaluation.................................................................................................................31

PAGE 8

viii LIST OF FIGURES Figure page 2-1. Reference linear translation curve was given by NCAT software...........................11 3-1. Eight cardiac beating frames of a NCAT phantom with beating heart motion only. (A) Frame 1 to 8 on trans-axial or ientation. (B) Frame 1 to 8 on coronal orientation.................................................................................................................21 3-2. Beating heart frame 1 at different resp iratory phases. (A) Tr ans-axial orientation. (B) Coronal orientation............................................................................................21 3-3. True Beating Heart Motion generated by NCAT software. Column (A): enddiastole; Column (B): end-systole............................................................................22 3-4. Projection data simulated from OHM......................................................................23 3-5. Projection data simulated from ORMA...................................................................23 3-6. Projection dataset ORMU showing linear axial movement of the heart due to respiratory motion in diffe rent angles. (A) Projecti on angle 1, (B) Projection angle 38, (C) Projection angle 60.............................................................................23 3-7. Projection dataset ORMU showing no tran s-axial movement of the heart due to respiratory motion in diffe rent angles. (A) Projecti on angle 1, (B) Projection angle 38, (C) Projection angle 60.............................................................................24 3-8. The relative center of mass for each beating heart frame over one respiratory cycle.........................................................................................................................2 4 3-9. Measured movement due to respirat ory motion of each respiratory phase for different gate (phase) schema for a 5second entire respir atory cycle. (A) 5 Phases. (B) 10 Phases. (C) 20 Phases.......................................................................25 3-10. Projection datasets OHM, ORMA, and ORMC. (A) OHM. (B) ORMA. (C) ORMC......................................................................................................................26 3-11. Reconstructed images by RM algorithm a nd profile plot from projection datasets OHM, ORMA, and ORMC. (trans-axial or ientation) (A) OHM. (B) ORMA. (C) ORMC. (D) Profiles.................................................................................................27

PAGE 9

ix 3-12. Reconstructed images by RM algorithm a nd profile plot from projection datasets OHM, ORMA, and ORMC. (coronal orient ation) (A) OHM. (B) ORMA. (C) ORMC. (D) Profiles.................................................................................................28 3-13. Reconstructed images by RM algorithm a nd profile plot from projection datasets OHM, ORMA, and ORMC. (sagittal orie ntation) (A) OHM. (B) ORMA. (C) ORMC. (D) Profiles.................................................................................................29 3-14. Segmentation map of myocardial region for image quality and motion estimation accuracy evaluation................................................................................30 3-15. SSEs and PMEs as a function of number of respiratory phases...............................31

PAGE 10

x Abstract of Thesis Presen ted to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Master of Science CORRECTION OF THE RESPIRATORY MOTION OF THE HEART FROM 4D MYOCARDIAL IMAGES By Jing Sun August 2006 Chair: David Gilland Major Department: Nuclear and Radiological Engineering During myocardial nuclear medicine imaging, respiratory motion blurs the projection data and thus causes diagnostic artifacts. A method is developed and evaluated in this work to compensate fo r effects of respiratory motion in four dimensional (4D) gated cardiac Emission Co mputed Tomography (ECT). Non-uniformrational-B-splines-based CArdiac Torso (NCA T) phantoms were generated to simulate temporal frames regularly spaced during combin ed respiratory and beating heart cycles. Projection data were acquired at 60 po sitions over 180 degrees to model gated myocardial perfusion Single Photon Emi ssion Computed Tomography (SPECT) with Technetium-99m (Tc-99m) sestamibi. A 2cm respiratory rigid body translation was simulated for the axial heart movement during normal tidal breathing. The motion was tracked and estimated by calculation of th e axial center-of-mass (aCOM). Respiratory motion correction was applied to the simu lated projection data. The 4D image reconstruction and 3D cardiac wall motion estimation for all temporal frames were

PAGE 11

xi performed simultaneously by a Joint 4D R econstruction and Motion estimation method (RM). The Sum of Squared Errors (SSE) wa s calculated to evaluate the image quality, and the Phantom-matching Motion Errors (PME ) was computed to evaluate the wall motion estimation accuracy. After respiratory motion compensation, the image quality, such as spatial resolution, of corrected images was better than that of respiratory blurred images, and very close to the quality of the images obtained when no respiratory motion was present. The accuracy of cardiac wall motion estimates was also improved with respiratory motion compensation.

PAGE 12

1 CHAPTER 1 INTRODUCTION 1.1 An Overview Respiratory motion takes place during m yocardial perfusion imaging and image blurring due to patient respira tion is thus unavoidable. This respiratory effect is not negligible; it degrades the diagnostic accura cy when appropriate respiratory motion compensation is absent. Many researchers have offered methods to reduce the impact of respiratory motion on image quality. One of the common methods is to let the patients hold their breath during imaging [1]. The main disadvantage of this method is that the patient needs to be highly c ooperative in order to obtain accurate images. Also, in this method, imaging has to be performed repeat edly and during each acquisition, the patient may not hold his/her breathe at the same leve l as the previous acquisitions, which will still lead to image blurring. Using an exte rnal device, for e.g. a sensor, to set up a respiratory gating system is an alternate me thod to reduce the respiratory motion effect [2-4]. Respiratory motion is not stable or pe rfectly repeatable like beating heart motion. Therefore, it is difficult to have a system that both measures the respiratory cycle accurately such as electrocardiography (ECG) measurements of the beating heart cycle, and is acceptable to both patients and technol ogists. Currently, Center-of-Mass is the most popular method used to track and es timate respiratory motion. This method requires only the patient projection data and ha s been realized in bot h digital [5, 6] and analogy adjustments of the axial signal [7, 8].

PAGE 13

2 1.2 Literature Review There are numerous papers published on blurring due to patient motion and correction methods. Recent works have pa id more attention on Positron Emission Tomography (PET) and SPECT studies. In 1995, using MRI imaging, Wang et al. [9], determined how much the heart shifted due to respiratory motion. It was concluded that during normal tidal breathing the movement of the heart due to respiration is dominated by superior-inferior (SI) motion, which is li nearly related to the SI motion of the diaphragm. Some gating systems were developed for SPECT and PET imaging to decrease the respiratory mo tion artifacts during imaging acquisition, for e.g., Macintosh platform and LabView environment r eal time system developed by Klein et al. [3], and respiratory gated SPECT (RGS) system realized by Cho [10]. Segars et al. [11, 12] have performed simulation studies to study resp iratory motion effects on myocardial SPECT with the spline-based Mathematical Card iac Torso (MCAT) phantoms. Respiratory mechanics were modeled in MCAT phantoms for SPECT imaging, and this work led to the improved and more realistic 4D NCAT phantom. Feng et al. [13] and Bruyant et al. [14] have developed several methods to comp ensate for effects of respiratory motion by applying center-of-mass algorithms to reconstruc ted images or projection datasets. They implemented their method with simulated MCAT phantoms. However, their work focused on respiratory motion only and the hear t wall was assumed to be stationary in their SPECT myocardial perfusion study. 1.3 Our Method The goal of this study is to develop a motion compensation method using center-ofmass algorithm to apply to 4D SPECT projec tion data simulated with both beating heart motion and respiratory motion. We implement this method to improve the reconstructed

PAGE 14

3 image quality and cardiac wall motion estimatio n accuracy, for the images reconstructed and motions estimated simultaneously by RM algorithm.

PAGE 15

4 CHAPTER 2 RESPIRATORY MOTION COMPENSATION, IMAGE RECONSTRUCTION, AND CARDIAC WALL MOTION ESTIMATION METHODS 2.1 Introduction To compensate respiratory motion effect s, Non-uniform-rational-B-splines-based Cardiac Torso (NCAT) phantoms were generate d to model the realistic organs in the torso, and simulate the beating heart moti on and respiratory motion [15]. The 2D projection data of NCAT phantoms were obtained by the Radon transform. Data acquisition was simulated in SPECT by using a parallel-hole collimator geometry [16]. The respiratory motions were estimated by the axial Center-of-Mass algorithm (aCOM) [14] for the projection data. The motion co mpensations were also applied to the 4D projections. The image reconstruction a nd 3D cardiac wall motion estimation for all temporal frames were performed simulta neously by a joint 4D Reconstruction and Motion estimation method (RM) [16]. The image quality before and after compensation was evaluated by calculations of the Sum of Squared Errors (SSE), and the myocardial wall motion estimation accuracy was evaluated by computations of the Phantommatching Motion Errors (PME) [17]. 2.2 NCAT Phantom Simulation The NCAT phantoms and the cardiac wall mo tion fields were generated by NCAT C language software develope d in John Hopkins University. Motion vectors were calculated and assigned in MATLAB to each voxel in each frame from the NCAT

PAGE 16

5 software outputs. Projection datasets were al so computed by using a linear interpolation projector in MATLAB. 2.2.1 Myocardial Phantom Simulation The NCAT phantom provides a realistic source distribution that models gated myocardial perfusion imaging with Tc-99m sestamibi by using the Non-uniform Rational B-Splines (NURBS) method to define the anat omy of a human body at a series of phases over cardiac cycles and respir atory cycles. This phantom can be used to study the development of image acquisition strategi es, image processing and reconstruction algorithms and compensation methods. It can be used to resear ch the effects of anatomical variations and pa tient motions, such as the be ating heart and respiratory motions, on medical images. NCAT phantom simulation software was developed at the Johns Hopkins University a nd fully described by Segars et al [11, 12]. In this work, a 40 temporal frames NCAT phantom was generated over one entire respiratory cycle with both beating hear t motion and respiratory motion by the NCAT software. The first frame corresponded to the full exhalation of the respiratory cycle and the end-diastole of the beating heart cycle. The length of respiratory cycle was set to 5 seconds as a normal tidal breathing. The respirat ory cycle started from full exhalation. At 0.455, almost in the middle, of the entire re spiratory cycle, motion reached the full inhalation. The length of cardiac cycle was set to 1 second as a normal male heart beating. Therefore, the phantom included one completed respiratory motion cycle and five completed beating heart motion cycles. For each cardiac cycle, there are 8 temporal frames. Frame 1 was the end-diastole, and frame 4 was the end-systole. Due to the respiratory motion, the diaphr agm motion and chest Anterio r-Posterior (AP) expansion were also simulated. The extent of th e maximum diaphragm motion was set to 2

PAGE 17

6 centimeters. This motion was simulated as a translation motion along the body axis. The extent of the maximum chest AP expansion wa s set to 1.2 centimeters. This motion only effect ribcage, body and lungs. Since hear t is the only organ in cluded in the NCAT myocardial phantom, there wa s no motion along trans-axial direction modeled in this present study. The projection data will show the heart movement due to respiratory motion in chapter 3. For the heart, the lengt h of the left ventricle myocardium was 9.43 centimeters, and the average radius of the left ventricle was 2.97 centimeters. The volume of left ventricle chamber at the enddiastole was 108.36 ml, a nd at the end-systole was 42.21 ml. So for this simulated NCAT phantom, the stroke volume (the amount of blood pumped by the left ventricle of the hear t in one contraction) was 66.15 ml and the ejection fraction (the fraction of blood pumped out of a ve ntricle with each heart beat) was 61%. For each frame, 50 slices of a ma le torso in supine position were created. Within each slice, the dimension of the phantom was 96 96. The voxel size was 0.30 centimeters. The phantom provides a realis tic source distribution that models gated myocardial perfusion imaging with Tc-99m sestamibi. Activities within myocardial region (including left ventricle myocardium, right ventricle myocardium, left atrium myocardium, and right atrium myocardium) we re assigned to 75 sour ce intensity per cm3 for Tc-99m, energy in 140 keV. For this study, only myocardial source activity was included. Outside activity was 0 everywhere. The parameters determined for beating heart cycle were modeled based on tagged MR images obtained from Johns Hopkins University [15]. And the parameters dete rmined for respiratory cycle were modeled based on the Visible Human Proj ect CT dataset from the Na tional Library of Medicine [15].

PAGE 18

7 In practice, any continuous function is inst antaneous. In order to be more realistic, for each time bin, the time interval for each temporal frame, four consecutive frames were created and averaged over the time bin. In this work, the study was concentrated on the myocardial region. Also to smooth the source distribution, for each frame, the 3D cardiac phantom (96 96 50) was obtained by linearly collapsing the original torso phantom, which is generated by NCAT softwa re directly with dimension of 256 256 168 voxels for computational efficiency, and was truncated into dimension 96 96 50, which included myocardial re gion, for computation efficiency purpose. Slice 15 to slice 64 was considered and chosen as myocardial region in this presen t study. And in each slice, pixel 17 to pixel 112 on both two dime nsions was used as the final myocardial phantom. As a reference, a phantom with 8 gated frames over 1-second cardiac cycle with only beating heart motion was generated, with the same size 96 96 50 voxels for each frame. A segmentation file with intensity equals to inside myocardium and intensity outside was also generated as a map, which defined the heart region by setting a threshold along the myocardial boundary. The source distributed objects were create d as noise-free phantoms. Attenuation effects were not considered for all phantoms generated using the NCAT software in this present study. 2.2.2 Cardiac Wall Motion Simulation Myocardial wall motion can be computed by NCAT software [15]. A 3D motion field between each two successive temporal fr ames over one complete beating heart cycle was calculated. For each beating heart cycle, one NCAT software motion calculation generated 8 separate output f iles if the cardiac cycle was gated into 8 frames. Each

PAGE 19

8 output file included starting points and ending points for and only for those non-zero voxels in the first frame, but not for each non-zero voxel in that frame. Furthermore, more than one motion vector might be assigned to the same non-zero voxel. In this case, the last assigned vector would be selected as the motion vector for this voxel when motion field between these two frames was computed from the output file. For each motion vector of a voxel, the starting point in the starting frame was an integer, and the ending point in the next frame was a real number. The norm of a motion vector was calculated as the distance between the starti ng point and the ending point. The direction of a motion vector was expressed as three co mponents, namely u(r), v(r), and w(r) along axes x, y, and z, respectively, where r = ( x,y,z) denoted the voxel coordinates in the spatial domain. Those motion vectors were c onsidered true motions for each voxel from one time bin to another. Output from the software was processed in MATLAB to calculate the norm and the three components fo r each motion vector of the voxels. When evaluated the motion fields estimated by RM algorithm, those motion vectors were considered as a reference set. Therefore, the NCAT software-calculated motion fields were named True motion in this present study. 2.2.3 Projection Simulation Forward projections were performed at 60 angles over 180 degrees, which starts from patient left hand side, went in front of the patient, and to the patient right hand side. The 2D projection data were computed by usin g a linear interpolation projector for each slice per angle per frame. Effects of attenua tion, scatter, and random were not considered during projection data acquisition. Clinically, 2000 counts in the myocardial region per frame per projection angle are realistic for gated Tc-99m sestamibi imaging [18] All projection data were scaled to this

PAGE 20

9 level. Therefore, the maximum voxel intens ity was assigned to 0.237. Poisson noise was added to all of the noise-free pr ojection data after scaling. Four projection datasets were computed in this study. The first dataset was calculated from the 8-frame beating heart motion only phantom in dimension 96 50 60 8. The dataset was named Object with beating Heart Motion only (OHM) in this present study. The second dataset was comput ed from the 40-frame phantom directly, in dimension 96 50 60 40. This projection was calculated as orig inal raw data for respiratory motion compensation, and was named Object with Respiratory Motion but Uncorrected (ORMU) in this present study. A third dataset was simulated as imaging with both beating heart motion and respirat ory motion, in dimensions of 96 50 60 8. The 40-frame both motions phantom included one complete respiratory cycle and five beating heart cycles at the same time. Five cardiac cycles divide d the respiratory cycle into five phases. For example, five phases fo r each respiratory cycle could be chosen as a respiratory gate schema. In this case, there were 8 temporal frames for each beating heart cycle over one respiratory phase. Each tem poral frame in one beating heart cycle was summed over respiratory phases. Then an averaged 8-frame projection dataset was obtained. It was blurred due to the re spiratory motion. The dataset was named Object with Respiratory Motion and Averaged (ORMA) in this present study. The last dataset was computed from 40-frame both motion pha ntom, too. The projection dimension was also as 96 50 60 8. Averaging was done for each temporal frame in each beating heart cycle over five respiratory phases, but after applying aCOM respiratory motion compensation method. The dataset was named Object with Respiratory Motion but

PAGE 21

10 Corrected (ORMC) in this present study. All three pr ojection datasets in dimension 96 50 60 8 were in the same noise level. 2.3 Axial Center-of-Mass Motion (aCOM) Estimation and Respiratory Motion Compensation Method The Center-of-Mass (COM) of each projection data per temporal frame was estimated in MATLAB. By calculating the CO M of the integrated projection data per respiratory phase simulated from ORMU the movement due to respiratory motion alone body axis can be determined. To compensate for respiratory effects, each frame was shifted respectively depending on the movement in each respiratory phase to reduce the blurring in each beating heart gate. 2.3.1 Axial Center-of-Mass (aCOM) Calculation According to Feng et al [13], the COM of each frame had the same rigid-body motion as the object in that frame. The phantom was simulated a 2-centimeter translation motion along the body axis. From full exhalation to full inhalation, the heart moved down linearly with diaphragm up to 2 centimeter s. The reference linear translation curve was given by NCAT software, and is show n in Figure 2-1. Example motion curve was given for 20 respiratory phases. COM of projection data computed from ORMU for each frame was calculated to track and estimate respiratory motion. To reduce the noise, the projection dataset had been integrated over th e trans-axial direction per projection angle, over all of the projection angles, over each b eating heart cycle, a nd over each respiratory phase. The center of mass al ong the body axis was determined for any respiratory phase as given by formula Equation (2-1):

PAGE 22

11 K j NT t M i NT t M i K jf t j i p f t j i p j f aCOM1111 111 1) , , ( ) , , ( ) ( (2-1) where M was the total number of voxels in the trans-axial dimension, N was total number of projection angles, T was tota l number of beating heart gate s, K was the total number of voxels in axial dimension, and p (i, j, t, f) was the activity in pixel (i, j) in the projection acquired at angle for beating heart time frame t in respiratory phase f. 0 5 10 15 20 25 05101520Respiratory PhasesMovement of Diaphragm (mm) Figure 2-1. Reference linear translati on curve was given by NCAT software Assume within each respiratory phase, that the respiratory motion was constant. Therefore, five aCOM were calculated finally for five phases over one complete respiratory cycle. Each one indicated the axial shift magnitude for the beating heart temporal frames within that respiratory pha se. The estimated relative center of mass for each respiratory phase was in unit number of voxels. The voxel size of projection data was 0.3 centimeters.

PAGE 23

12 2.3.2 Respiratory Motion Compensation Use equation (2-1) to calcul ate the COM of the projecti on dataset for each beating heart frame. COMs of 40 temporal frames in total were computed for the five cardiac cycles, as well one entire respiratory cycle. For motion compensation, each frame of projection dataset ORMU was vertically shifted because only axial movement was consid ered in this study. Full exhalation in the respiratory cycle and end-diastole in the beat ing heart cycle was defined as the reference status with zero relative movement in this wo rk. Therefore, the last respiratory phase was chosen as the reference phase, for which the projection data were not shifted. All other projections were shifted up. Based on th e assumption that re spiratory motion was constant within each respiratory phase, the be ating heart frame within each respiratory phase will shift same amount. And the magnitude of the shifting was depended on the movement calculated by Equation (2-1) of th e phase which the beating heart frame was located. Because the axial shift magnitude wa s not an integer number of voxels, a linear interpolation was performed along the body axis to obtain an accurate shifting. After respiratory motion compensation, the COM of each shifted beating heart frames were closed to the COM of the refe rence phase. The magn itude of the errors depended on how many respiratory phases were chosen. Because, in practice, there are no very accurate devices such as ECG mon itors used for cardiac gated imaging, the narrow respiratory gates could not be obtained. Five respiratory phases were chosen because of the simplicity, so that each re spiratory phase was reasonably assigned as 1 second. Ten phases and twenty phases respir atory schema were also simulated for optimization and comparison. The aCOM meth od was also applied to the corrected

PAGE 24

13 projection dataset to check the constancy of the aCOM position during the respiratory cycle. The completed results are listed and shown in Chapter 3. Conclusively, the motion compensation al gorithm can be described as following steps: 1. Calculate the center-of-mass of each respiratory phase. 2. Choose the respiratory phase which has l east movement (the phase has smallest center-of-mass value) as the reference pha se, and calculate the movements of all other phases by giving relative to the reference phase. 3. Shift all beating heart frames in each re spiratory phase axially according to the relative movement calculated for that phase. 4. Sum each beating heart frame over all of the respiratory phases. 2.4 Image Reconstruction and Cardiac Wal l Motion Estimation (RM) Method Projections simulated from OHM ORMA and ORMC were all reconstructed by image reconstruction and cardiac wall motion esti mation (RM) algorithm. This algorithm was first discussed by Mair et al. [19] and developed in MATLAB code by Cao et al. [20]. The method was later improved and implemented in Fortran code by Gilland et al. [21] and Mair et al. [17]. 2.4.1 RM Algorithm RM algorithm was a method for gated car diac ECT that reconstructed the voxel intensities of the gated images and estimated the 3 dimensional motion of the cardiac wall for all frames simultaneously. The algorithm consisted of two steps, the image reconstruction step (R step) a nd the cardiac wall motion estimati on step (M step). It was based on an iterative algorithm for minimi zing an objective functi on that included the negative log likelihood of the data, the sta ndard optical brightness constraint, and a biomechanical model for elastic deformations of the myocardium.

PAGE 25

14 In the early version, Mair et al. [19] discussed the estima tion of motion and images for only two frames with maximum deformation, end-diastole and end-systole. During R step, the motions were fixed and the object ive function was minimized over the frames iteratively. During M step, the frames we re fixed and the objective function was minimized over the motion field iterativel y based on solving the EulerLagrange equations [22, 23]. The algorithm was im plemented in MATLAB because of its extensive library of mathematical functions However, it was rewritten in Fortran language by Gilland et al. [21] to enhance th e computational efficiency. And the most important point was that, in Mair et al. [17] the algorithm was extended from two frames to a complete cardiac cycle frames, which can be any finite number. Eight frames were chosen in this simulation study. Instead of being performed in pairs of frames, a mathematical framework, which included the entire cycle of images and motion fields, was developed to make only one single objective function. In this manner, all frames were considered as a loop. Each update of a frame involved the cu rrent estimate of the frame immediately before and after this frame. The objec tive function minimization was performed directly by using a type of se quential quadratic programming with the conjugate gradient algorithm [ 17] applied at each step. Afte r all of the improvements and optimizations, the algorithm guaranteed non-nega tivity of the reconstructed images, and each iteration decreased the value of the obj ective function. The simultaneous estimation allowed the data to influence both im age reconstruction and cardiac wall motion estimation at the same time [17]. 2.4.2 Objective Function In this study, the myocardium was consid ered as a deformable elastic material satisfying the standard constitutive relati ons in continuum m echanics [24, 25].

PAGE 26

15 2.4.2.1 Negative log likelihood term It was assumed that the small deformation between two successive frames would not result in very large st rain energy [23]. If f1(r), f2(r), . fT(r) and m1(r), m2(r), . mT(r) were denoted T frames image and T mo tion fields between each successive two frames respectively, the important non-nega tive condition on each intensity function ft will be denoted by f 0 [17]. Then the negative l og likelihood term of an objective function was defined as Equation (2-2) [17]: T t M i t t ti Hf i g i Hf f L11 ) () ( log ) ( ) (, (2-2) where M was the total number of the detector bins, the dimension of the phantom, H was the projection operator, and g1(t), g2(t), . gT(t) were gated datasets. 2.4.2.2 Image matching term When the image intensities between frames t and t + 1 perfectly matched the motion field mi between them, ft(r) should equal to ft+1(r+mt(r)) for all r [17]. In practice, image frames and motion fields were never exactly matched. Therefore, the image matching term of an objective function between frames t and t + 1, which defined as Equation (2-3) [17], was expected to be small: dr r m r f r f m f f Et t t t t t I 2 1 1)) ( ( ) ( ; , (2-3) where r was the voxel coordinate, and r = (x, y, z). If consid ered the entire cardiac cycle, the total image matching term was defined as Equation (2-4) [17]: T t t t t I Im f f E m f1 1; ;. (2-4)

PAGE 27

16 2.4.2.3 Strain energy term By treated myocardium as a deformable el astic material, material strain resulting from the cardiac motion can be modeled by th e strain energy [17]. The strain energy depended on the derivatives of the mo tion vector and material parameters and namely Lam Constants [15]. After several experiments, it was found that the material difference can be ignored in this study. Next, and were set equal to one both inside and outside the myocardium so that Youngs modulus was (3 + 2)/( + ) = 2.5 and Poissons ratio was 0.5 /( + ) = 0.25 respectively In th is study. Segmentation process by Mair et al. [17] was still included in this work but only applied for evaluation method for image quality and motion accuracy presen ted in Chapter 3. For each motion field between frame t and t + 1, mt had three components ut, vt, and wt. By applying the strain energy function in Loves classical book [15] the strain energy term of an objective function was defined as Equation (2-5) [17]: dr w v u r dr w v u r E where m E mz t y t x t z t y t x t S T t t S S 2 2 2 2 , 1) ( ) ( 2 1 (2-5) dr w v w u v u ry t z t x t z t x t y t 2 , 2 , 2 ,) ( 2 1 where ut,x, vt,y, wt,x, . were continuous first partials. e.g. x u ut x t ,. and were included here for general expres sion. They were assigned to ones in this study. Finally, the total objective function was defined as Equation (2-6) [17]: ) ( ) ; ( ) ( ) ; ( m m f E f L m fS I (2-6)

PAGE 28

17 where and were hyper-parameters, which reflected the influence that image matching, projection data and strain energy have on th e estimations. Upon several experiments by Mair et al. [17], the value of and were determined as 0.015 and 0.0035 respectively in this study. If f1 *(r), f2 *(r), . fT *(r) and m1 *(r), m2 *(r), . mT *(r) were denoted true images and true motion fields respectively, minimized function resulted as Equation (2-7) [17]: ) ; ( min arg ) ; (, 0 *m f m fm f (2-7) In RM algorithm, the initials were unifor med images and zero motion fields. The reconstructed images were updated first while motion fields were fixe d (R Step). Then the estimated cardiac wall motions were upda ted while image frames were fixed (M Step). Within R Step or M Step, the algor ithm went out of loop wh en it was converged. After one R Step and one M Step, each imag e and motion field was updated once equally. In R step, image matching and negative l og likelihood terms of the objective function were used to improve the reconstructed imag es. The expected estimate for each iterative R step was defined as Equation (2-8) [17]: ) ; ( min arg) ( 0 1 n f nm f f (2-8) In M step, updated frames in latest R step were used to improve estimated motions. The expected estimate for each iterative M st ep was defined as Equation (2-9) [17]: ) ; ( min arg) 1 ( 1m f mn m n (2-9) Experimental results of the reconstr ucted images and estimated cardiac wall motions are presented in Chapter 3.

PAGE 29

18 2.4.3 Reconstruction and Motion Estimation 100-iteration RM algorithm, with 100iteration R step a nd 100-iteration M alternately, was applied to all projection datasets ORMA ORMU and ORMC Fortran code was running on Penguin Tempest 2100 (Dua l-processor 64-bit AMD Opteron, 8 GB physical memory). It took about 4 hours to finish 100 iterations. The value of and were assigned as 0.015 and 0.0035 respectively according to Mair et al. [17]. 2.5 Image Quality and Cardiac Wall Moti on Accuracy Evaluation Methods The image quality of the reconstructed im age was evaluated using Sum of Squared Errors (SSE) method. The accu racy of the cardiac wall mo tion estimation was evaluated using images matching term of the objectiv e function. For the evaluation, phantom OHM was using as the reference, and the image matching equation was also named Phantommatching Motion Error (PME) method. 2.5.1 Image Quality Evaluation For image quality evaluation, SSE method was defined as Equation (2-10) [17]: 2 1) ( ) ( ) ( T tr t tr f r f f f SSE (2-10) where ) ( ),... ( ), (2 1r f r f r fTand ) ( ),... ( ), (2 1r f r f r fT were phantom source distributed object and reconstructed image, respectively. Only the voxels inside the myocardial region would be summed. A segmentation map was generated to define the moycardium by given a threshold intensity. In this study, only heart was included in the phantom. Therefore, the intensity threshold was assign ed as zero. In another word, any non-zero intensity voxel was considered within m yocardium. The segmentation map and the

PAGE 30

19 experimental results were shown in Chapter 3. SSE evaluation method was also used for value optimization study [17]. 2.5.2 Cardiac Wall Motion Evaluation Recall image matching equation: dr r m r f r f m f f Et t t t t t I 2 1 1)) ( ( ) ( ; (2-3) For cardiac wall motion evaluation, same equation was used but called PME method defined as Equation (2-11) [17]: 2 1 1)) ( ( ) ( ) ( T tr t t tr m r f r f m PME (2-11) where ) ( r ft was the phantom source distribut ed object at frame t, and ) ( r mt was the estimated motion field from frame t to frame t + 1. If frame t is the last frame in a sequence (e.g. the 8th frame in a 8-gate cardiac image cycle), frame t + 1 would be assigned to be the first frame in the cycle. In this evaluation st udy, only those non-zero intensity voxels were consider ed within myocardium and were computed for PME value. The segmentation map was used, too. The expe rimental results were shown in Chapter 3. PME evaluation method was used for value optimization study, too [17].

PAGE 31

20 CHAPTER 3 EXPERIMENTAL RESULTS AND DISCUSSION 3.1 Simulated Phantoms, True Motion Field, and Projections Experimental results and analysis are presente d in this Chapter. Calculated relative center-of-mass of projections before and afte r respiratory compensa tion, reconstructed images, estimated cardiac wall motion, and SSE and PME evaluation data are illustrated in the following sections. 3.1.1 Cardiac Gated NCAT Phantom 3.1.1.1 Reference phantom with beating heart motion only As a reference phantom, eight cardiac beat ing frames are shown in figure 3-1 in trans-axial and coronal orientations of a pha ntom with beating heart motion only. The dimension of the phantom was 96 96 50 for each gating frame. End-diastole and end-systole can be seen clearly at frame 1 and frame 4, respectively. 3.1.1.2 Phantom with both beating heart motion and respiratory motion To look at whether respiratory motion ha d an effect on myocardial imaging, a phantom with combined beating heart motion and respiratory motion is shown in Figure 3-2. If the length of the respiratory cycle wa s assigned as 5 seconds, and the length of the beating heart cycle was assigned as 1 second, there were 5 complete beating heart cycles during one respiratory cycle. The first be ating heart frame for each one of the five beating heart cycles at different respirat ory phases are displayed in Figure 3-2. Myocardium shifted down about 2 centimeters from full exhalation to full inhalation with the diaphragm.

PAGE 32

21 Figure 3-1. Eight cardiac beating frames of a NCAT phantom with beating heart motion only. (A) Frame 1 to 8 on trans-axial or ientation. (B) Frame 1 to 8 on coronal orientation. Figure 3-2. Beating heart frame 1 at different respiratory phases. (A) Trans-axial orientation. (B) Coronal orientation. 3.1.2 True Motion Field of the NCAT Phantom The true motion field for the beating h eart motion only phantom was superimposed on the reference object and is plotted in Figur e 3-3. Frames in the left column showed the end-diastole. All motion vectors were pointing to the direc tion where the heart contracted to during systole. Frames in th e right column show th e end-systole. All motion vectors point to the direction wher e the heart enlarged to during diastole. (A) (B) (A) (B)

PAGE 33

22 Figure 3-3. True Beating Heart Mo tion generated by NCAT software. Column (A): end-diastole; Column (B): end-systole. 3.1.3 Projection Datasets Projection dataset OHM is shown in Figure 3-4. Dataset ORMA is shown in Figure 3-5. Only projection angle 38 of each data set is shown in these two figures. It was observed that cardiac gated projection data were blurred by th e respiratory motion obviously. Hence, the respiratory motion eff ect could not be ignored during the SPECT imaging. Figure 3-6 and Figure 3-7 show the first frame for each beating heart cycle over respiratory phases (5 phases shown in this example). Projection data angle 1, 38 and 60 are shown in these two figures. Linear axia l movement can be seen vertically from Figure 3-6. Purple line indicates the heart position without respirat ory effect. No trans(A) (B)

PAGE 34

23 axial movement horizontally is found as expected how the phantom been modeled from Figure 3-7. Figure 3-4. Projection data simulated from OHM. Figure 3-5. Projection data simulated from ORMA. Figure 3-6. Projection dataset ORMU showing linear axial move ment of the heart due to respiratory motion in different angles. (A) Projection angle 1, (B) Projection angle 38, (C) Projection angle 60. 3.2 aCOM Calculation and Respiratory Motion Compensation Method The relative COM along axial dire ction of projection datasets ORMU and ORMC were presented in the following sections. The COM of each respiratory phase over the cycle was also shown here. 3.2.1 aCOM Calculation Figure 3-8 illustrated the rela tive COM of projection dataset ORMU along axial ( A ) ( B ) ( C )

PAGE 35

24 Figure 3-7. Projection dataset ORMU showing no trans-axial movement of the heart due to respiratory motion in different angles. (A) Projection angle 1, (B) Projection angle 38, (C) Projection angle 60. direction calculated by using Equation (2 -1) for each cardiac frame over one entire respiratory cycle. Compare to the 2-centimete r linear translation due to the respiratory motion, the magnitude of the movement due to the heart beating (average about 0.15 centimeters) was unobservable. Therefor e, only the displacement caused by the respiratory motion could be seen from the chart. 0 1 2 3 4 5 6 7 8 16111621263136Beating Heart Frames over one Respiratory CycleRelative Displacement of Center of Mass (number of voxels) Figure 3-8. The relative center of mass for each beating heart frame over one respiratory cycle ( A ) ( B ) ( C )

PAGE 36

25 3.2.2 Respiratory Motion Compensation Five phases (gates) over one 5-second en tire respire cycle was calculated first because of the simplicity. Ten phases and tw enty phases respiratory schema were also simulated for comparison and optimization. Al l three gate schema were shown in Figure 3-9. The height of each column indicates how much the beating heart frames within that respiratory phase need to be shifted. Figure 3-9. Measured movement due to resp iratory motion of each respiratory phase for different gate (phase) schema for a 5-second entire re spiratory cycle. (A) 5 Phases. (B) 10 Phases. (C) 20 Phases. (B) (A) -1 1 3 5 7 12345 Respiratory Phasesmovement of each respiratory phase (number of voxels) -1 1 3 5 7 12345678910 Respiratory Phasesmovement of each respiratory phase (number of voxels) -1 1 3 5 7 1234567891011121314151617181920 Respiratory Phasesmovement of each respiratory phase (number of voxels) (C)

PAGE 37

26 After respiratory motion compensation, th e relative COM for each beating heart frame should be very close to the reference phase. The ideal illustration should be a straight horizontal li ne. The error depends on which gate schema was used. The more gates used for each respiratory cycle, the le ss errors were. Twenty phases schema would be chosen according to the clinical realization by Cho et al. [10]. The projection datasets of projection angle 38 after ( ORMC ) 20-phase respiratory moti on correction is shown in Figure 3-10. Respiratory moti on blurred projection data ( ORMA ) and heart motion only projections ( OHM ) were shown in the same figure fo r comparison. The projection data can be seen from the images clearly that it was very close to the reference projections ( OHM ) after motion compensation. In practice, the respiratory motion might not be very stable. The length of the respiratory cycle might cha nge during projection data acquisition. Over sampling of the gates a nd averaging over each chosen respiratory phase would contribute to the mo tion compensation improvements. Figure 3-10. Projection datase ts OHM, ORMA, and ORMC. (A) OHM. (B) ORMA. (C) ORMC. 3.3 Image Reconstruction, Cardia c Wall Motion Estimation, and SSE and PME Evaluations. The RM reconstructed images from all three projection datasets ( OHM ORMA and ORMC ) are shown in this section. SSE and PME ev aluations are listed in the tables here. (A) (B) (C)

PAGE 38

27 3.3.1 Image Reconstruction a nd Wall Motion Estimation The reconstructed images and profile plots from all three projections ( OHM ORMA and ORMC ) for three orthogonal slice orientat ions trans-axial, coronal, and sagittal are shown in Figure 3-11, Figur e 3-12 and Figure 3-13, respectively. 20 respiratory phases schema were chosen to show result in this study according to the clinical level de termined by Cho et al. [10]. 0 0.1 0.2 0.3 0102030405060 VoxelsIntensities (counts) OHM images ORMA images ORMC images Figure 3-11. Reconstructed images by RM al gorithm and profile pl ot from projection datasets OHM, ORMA, and ORMC. (trans-axial orientation) (A) OHM. (B) ORMA. (C ) ORMC. (D) Profiles (A) (B) (C) (D)

PAGE 39

28 0 0.1 0.2 0.3 0510152025303540 VoxelsIntensities (counts) OHM images ORMA images ORMC images Figure 3-12. Reconstructed images by RM al gorithm and profile plot from projection datasets OHM, ORMA, and ORMC (coronal orientation) (A) OHM. (B) ORMA. (C ) ORMC. (D) Profiles (A) (B) (C) (D)

PAGE 40

29 0 0.1218 0.2436 0.3654 010203040 VoxelsIntensities (counts) OHM images ORMA images ORMC images Figure 3-13. Reconstructed images by RM al gorithm and profile plot from projection datasets OHM, ORMA, and ORMC (sagittal orientation) (A) OHM. (B) ORMA. (C ) ORMC. (D) Profiles From the reconstructed images in all orie ntations, it was obser ved that the blurring due to the respiratory motion had been redu ced greatly by motion compensation method. The image reconstructed from corrected projections was much closer to the image reconstructed from the reference dataset. ORMC images looked smoother than true ( OHM ) images. That was because ORMC images were summed over respiratory phases after motion compensation. From the prof ile plots, it was observed that both OHM images and ORMC images peak were narrow and had higher peak values. It demonstrated that the noise control was impr oved, and the less blurri ng on the edge of the (A) (B) (C) (D)

PAGE 41

30 myocardial region could be s een, too. Estimated wall motion would be evaluated by PME method explained in next section. 3.3.2 SSE and PME Evaluation Method The segmentation map used for SSE and PM E evaluation study is shown in Figure 3-14. Slice 27 in frame 1 was chosen to show here. Figure 3-14. Segmentation ma p of myocardial region fo r image quality and motion estimation accuracy evaluation The SSE evaluation results of OHM ORMA and ORMC images calculated by Equation (2-10) are listed in Table 31. The PME evaluation results of OHM ORMA and ORMC motion fields calculated by Equation (2-11) are listed in Table 3-2. The SSEs and PMEs values as a function of number of respiratory phases are shown in Figure 3-15. Table 3-1. The SSE results of OHM, ORMA and ORMC images quality evaluation. Sum of Squared Errors (SSE) True Object 0 OHM image 472.48 ORMA image 1503.93 ORMC image (5 phases) 643.57 ORMC image (10 phases) 560.54 ORMC image (20 phases) 521.79

PAGE 42

31 Table 3-2. The PME results of OHM, ORMA, and ORMC motion estimation accuracy evaluation. Phantom-matching Motion Errors (PME) True motion 88.25 OHM motion 183.75 ORMA motion 386.35 ORMC motion (5 phases) 307.37 ORMC motion (10 phases) 290.85 ORMC motion (20 phases) 272.91 From the SSEs and PMEs table, re spiratory motion compensation method demonstrated improved SSE error as well as improved image noise and spatial resolution characteristics. Both SSE and PME would decrease with the increased number of respiratory phases. Also accord ing to the clinical study by Cho et al. [10], 20 phases would be an ideal choice fo r respiratory-gated study. 0 200 400 600 800 0510152025 Number of Respirator y Phases SSEs and PME s SSEs PMEs Figure 3-15. SSEs and PMEs as a functi on of number of respiratory phases.

PAGE 43

32 CHAPTER 4 CONCLUSIONS AND FUTURE WORK 4.1 Summary and Conclusions A respiratory motion compensation method was discussed and developed in this study. Method was implemented by a NCAT pha ntom simulation with both beating heart motion and respiratory motion. Respiratory motion was tracked and estimated by aCOM algorithm. Motion compensation was realized by applying method on to the projections computed from the NACT phantoms with a li near interpolation ope rator. The true projection, averaged projection and the correct ed projection datasets were reconstructed by RM algorithm. The cardiac wall motion wa s estimated simultaneously with the image reconstruction. The respiratory moti on compensation method was evaluated by comparison between OHM ORMA and ORMC images and motion fields. Numerical results were explained by calcula tion of SSEs and PMEs errors. From the projections, it was observed th at cardiac gated pr ojection data were blurred by the respiratory motion. From aCOM calculation, the respiratory motion can be tracked and estimated for each respiratory phase. After respiratory motion compensation by aCOM computations, the corrected projection data are very close to th e data without resp iratory motions. It can be seen that blurring are reduced on the reconstructed images from the corrected projection data. From the evalua tion results, the values of SSE and PME demonstrated improved spatial resolution and more accuracy of wall motion estimation than that without respiratory motion compensation. Different respiratory phases schemas have been tested. More respiratory gating sampling gives more image quality improvement and more accuracy on wall

PAGE 44

33 motion estimation. According to the clinic al patient data study of other papers, 20 respiratory phases are chosen for this study. 4.2 Future Work The study included only the heart for the NCAT phantom simulation. Since for gated myocardial perfusion SPECT imaging with Tc-99m, the radiopharmaceuticals will generate intensities in the other organs, such as liver, al most same level relative to the heart. That may have big effect on the cen ter-of-mass calculation to track and estimated the respiratory motion. Future studies usi ng phantoms with more organs can make this compensation method more general. In this study, it was assumed that the re spiratory motion was c onstant within one respiratory phase. In practice, the respirator y gate can not be too narrow. This will bring errors to the simulation with both beating he art and respiratory mo tion, because within one respiratory phase (e.g. less than 5-phase schema), there are more than one beating heart frames included. One beating heart frame was the time bin which was finally implemented. Further consideration of more cardiac cycle and respiratory cycle combination may yield more accurate results. Respiratory motion e ffect will take place in almost every myocardial perfusion SPECT imaging examination. In practice, the respiratory cycle may not be fixed during imaging, and the beating heart motion and respiratory motion combination could be more complicated. Applying this compensation method to clinical patient data w ould be an important future study.

PAGE 45

34 LIST OF REFERENCES 1. Gottschalk A, Harper PV, Jimenez FF, Petasnick JP. Quantification of the respiratory scanning with the rectilinear focused collimator scanner and the gamma scintillation camera. Journal of Nuclear Medicine. 1966;7:243. 2. Zubal G, Bizais Y, Bennett GW, Brill AB Dual gated nuclear cardiac images. IEEE Transactions on Nuclear Science. 1984;31:566. 3. Klein GJ, Reutter BW, Ho MH, Reed JH Huesman RH. Real-time system for respiratory-cardiac gating in positron tomography. IEEE Transactions on Nuclear Science 1998;45(4):2139. 4. Klein GJ, Reutter BW, Huesman RH. Data-d riven respiratory gating in list mode cardiac PET. Journal of Nuclear Medicine. 1999;40(5):113P. 5. Oppenheim BE. A method using a digital compute for reducing respiratory artifacts on liver scans made with a camera. Journal of Nuclear Medicine. 1971;12:625 628. 6. Schmidlin P. Development and comparison of computer methods for organ motion correction in scintigraphy. Journal of Nuclear Medicine. 1999;40(5):113P. 7. Hoffer PB, Oppenheim BE, St terling ML, Yasillo N. A simple device for reducing motion artifacts in gamma camera images. Radiology 1972;103:199. 8. Baimel NH, Bronskill MJ. Optimization of analog-circuit motion correction for liver scintigraphy. Journal of Nuclear Medicine. 1978;19:1059. 9. Wang Y, Riederer SJ, Ehman RL. Respirat ory motion of the heart: kinematics and the implications for the spatial resolution in coronary imaging. Magnetic Resonance in Medicine. 1995;33:713. 10. Cho K, Kumiata S, Okada S, Kumazaki T. Development of respiratory gated myocardial SPECT system. Journal of Nuclear Cardiology. 1999;6:20. 11. Segars WP, Lalush DS, Tsui BMW. Mode ling respiratory mechanics in the MCAT and spline-based MCAT phantoms. IEEE Transactions on Nuclear Science. 2001;48(1):89.

PAGE 46

35 12. Segars WP, Tsui BMW. Study of the Efficacy of Respiratory Gating in Myocardial SPECT Using the New 4-D NCAT Phantom. IEEE Transactions on Nuclear Science. 2002;49(3):675. 13. Feng B, Bruyant PP, Pretorius PH, Beach RD, Gifford HC, Dey J, Gennert M, King MA. Estimation of the rigid-body motion fr om three-dimensional images using a generalized center-ofmass points approach IEEE Nuclear Science Symposium Conference Record. 2005;M07:248. 14. Bruyant PP, King MA, Pretorius PH. Correct ion of the respiratory motion of the heart by tracking of the center of mass of thresholded projections: a simulation study using the dynamic MCAT phantom IEEE Transactions on Nuclear Science 2002;49(5):2159. 15. Segars WP. Development of a new dyna mic NURBS-based cardiac torso (NCAT) phantom. Ph.D. dissertation The University of North Carolina, Chapel Hill, 2001. 16. Cao Z, Simultaneous reconstruction and 4D motion estimation for gated myocardial emission tomography. Ph.D. dissertation The University of Florida, Gainesville, 2003. 17. Mair BA, Gilland DR, Sun J, Estimation of images and non-rigid deformations in gated emission CT Transactions on Medical Imaging. Accepted in April, 2006. 18. Narayanan MV, King MA, Soares EJ, By rne CL, Pretorius PH, Wernick MN. Application of the Karhunen-Loeve transfor m to 4D reconstruc tion of cardiac gated SPECT images. IEEE Transactions on Nuclear Science 1999;46(4):1001. 19. Mair BA, Gilland DR, Cao Z. Simultaneous motion estimation and image reconstruction from gated data IEEE International Symposium on Biomedical Imaging 2002:661. 20. Cao Z, Gilland DR, Mair BA, Jaszczak RJ. Three-dimensional motion estimation with image reconstruction for gated cardiac ECT IEEE Transactions on Nuclear Science. 2003;50(3):384. 21. Gilland DR, Mair BA, Sun J. Joint 4D reconstruction and motion estimation in gated cardiac ECT. Intern Conference on Fully 3D Image Reconstruction Radiation Nuclear Medicine. 2005:303. 22. Klein J, Reutter BW, Huesman RH. Nonrig id summing of gated PET via optical flow. IEEE Transactions on Nuclear Science 1997;44(4):1509. 23. Klein GJ, Huesman RH. Elastic material model mismatch effects in deformablemotion estimation. IEEE Transactions on Nuclear Science 2000;47(3):1000.

PAGE 47

36 24. Love AEH, A Treatise on the Math ematical Theory of Elasticity. Dover Publications New York, 1944. 25. Humphrey JD, Delange SL. An introduction to biomechanics: solids and fluids, analysis and design. Springer-Verlag New York, 2004.

PAGE 48

37 BIOGRAPHICAL SKETCH Jing Sun was born in Jilin, China, on February 28th, 1976. In September 1994, she began her college education in Tsinghua Univ ersity in Beijing and obtained her Bachelor of Science degree in engineering physics in July 1999. In 2002, she obtained her Master of Science degree in nuclear engineeri ng from North Carolina State University. In August 2003, Jing was admitted to the University of Florida as a graduate student to the medical physics program in th e Department of Nucl ear and Radiological Engineering. During the last three years, she has worked with Dr. David Gilland in SPECT myocardial imaging project. Jing Sun is the only daughter of Zhongzhi Sun and Guirong Liu and is married to Dapeng Ding.


Permanent Link: http://ufdc.ufl.edu/UFE0015877/00001

Material Information

Title: Correction of the Respiratory Motion of the Heart from 4D Myocardial Images
Physical Description: Mixed Material
Copyright Date: 2008

Record Information

Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
System ID: UFE0015877:00001

Permanent Link: http://ufdc.ufl.edu/UFE0015877/00001

Material Information

Title: Correction of the Respiratory Motion of the Heart from 4D Myocardial Images
Physical Description: Mixed Material
Copyright Date: 2008

Record Information

Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
System ID: UFE0015877:00001


This item has the following downloads:


Full Text












CORRECTION OF THE RESPIRATORY MOTION OF THE HEART
FROM 4D MYOCARDIAL IMAGES















By

JING SUN


A THESIS PRESENTED TO THE GRADUATE SCHOOL
OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT
OF THE REQUIREMENTS FOR THE DEGREE OF
MASTER OF SCIENCE

UNIVERSITY OF FLORIDA


2006

































Copyright 2006

by

Jing Sun

























Dedicated to my beloved parents and my husband



A Ufta i!r RafH















ACKNOWLEDGMENTS

First of all, I would like to express my deep indebtedness and my sincere gratitude

to my advisor, David Gilland, for the valuable guidance, support and continuous

encouragement he has given me during the course of this work. His insight in this field

and his thoughtful personality have made my study here worthwhile and unforgotten. I

am also equally indebted to the members of my committee, Wesley Bolch and David

Hintenlang, for their helpful advice, suggestions and encouragements. Special thanks are

also extended to Bernard Mair, the external committee member, for his suggestions and

attendance. It is a great honor to have them on my committee.

I am grateful to Willams Paul Segars, for his tremendous help with the NCAT

phantom. I would like to thank Dr. Zixiong Cao for his patience and help to set up this

research idea. I am also grateful to Uday Tipnis, Keelan Seabolts, and Jason Parker, who

have inspired and assisted in my research study.

Finally I would like to extend my sincere gratitude to my parents for their moral

support from far away in my home country, China. I would like to thank my husband

Dapeng for his love, understanding and moral support.
















TABLE OF CONTENTS



A C K N O W L E D G M E N T S ................................................................................................. iv

L IS T O F T A B L E S ................. ........................................................ ............ ................. ... v ii

L IST O F FIG U R E S ................................................................ ........ ............... viii

A B ST R A C T .......... ..... ...................................................................................... x

CHAPTER

1 IN TR OD U CTION ............................................... .. ......................... ..

1.1 An Overview................... ........ .. ... .... ............... ....
1.2 Literature R review ..............................................................2
1.3 O ur M ethod ............................................................... ...... 2

2 RESPIRATORY MOTION COMPENSATION, IMAGE RECONSTRUCTION,
AND CARDIAC WALL MOTION ESTIMATION METHODS ............................4

2.1 Introduction..................................................... ................... .. ....... ...... 4
2.2 N CA T Phantom Sim ulation..................................................... ...................4
2.2.1 M yocardial Phantom Simulation................................. ............ ............. 5
2.2.2 Cardiac W all M otion Simulation................................. ............. ........... 7
2.2.3 Projection Sim ulation ............................................................................. 8
2.3 Axial Center-of-Mass Motion (aCOM) Estimation and Respiratory Motion
Com pensation M ethod .........................................................................................10
2.3.1 Axial Center-of-Mass (aCOM) Calculation......................... ............10
2.3.2 Respiratory M otion Compensation.................................................... ...12
2.4 Image Reconstruction and Cardiac Wall Motion Estimation (RM) Method........13
2.4.1 R M A lgorithm ........................ .. .................................. .. ........... ..13
2 .4 .2 O bjectiv e F u n action .......................................................... ..................... 14
2.4.2.1 Negative log likelihood term ............... ..........................................15
2.4.2.2 Im age m watching term ............................................... .................. 15
2.4.2.3 Strain energy term .................................... ............................ ...... 16
2.4.3 Reconstruction and M otion Estimation....................................................18
2.5 Image Quality and Cardiac Wall Motion Accuracy Evaluation Methods............ 18
2.5.1 Im age Q quality Evaluation....................................... ......... ............... 18
2.5.2 Cardiac W all M otion Evaluation..... .......... ...................................... 19









3 EXPERIMENTAL RESULTS AND DISCUSSION..............................................20

3.1 Simulated Phantoms, True Motion Field, and Projections............................. 20
3.1.1 Cardiac Gated NCAT Phantom ............................................................ 20
3.1.1.1 Reference phantom with beating heart motion only ......................20
3.1.1.2 Phantom with both beating heart motion and respiratory motion....20
3.1.2 True Motion Field of the NCAT Phantom ..............................................21
3.1.3 P rejection D atasets ............................................................................. .... 22
3.2 aCOM Calculation and Respiratory Motion Compensation Method ...................23
3.2.1 aC O M C calculation ............................................. ............................. 23
3.2.2 Respiratory M otion Compensation...........................................................25
3.3 Image Reconstruction, Cardiac Wall Motion Estimation, and SSE and PME
Evaluations..................................................... ........... ............ .. 26
3.3.1 Image Reconstruction and Wall Motion Estimation ................................27
3.3.2 SSE and PM E Evaluation M ethod .................................. ............... 30

4 CONCLUSIONS AND FUTURE WORK.....................................................32

4.1 Sum m ary and Conclusions ............................................................................32
4 .2 F u tu re W o rk .................................................................................................... 3 3

L IST O F R EFE R E N C E S ............................................................................. ............. 34

BIOGRAPH ICAL SK ETCH ...................................................... 37
















LIST OF TABLES


Table


3-1. The SSE results of OHM, ORMA, and ORMC images quality evaluation ............30

3-2. The PME results of OHM, ORMA, and ORMC motion estimation accuracy
ev alu action ......................................................................... 3 1


page















LIST OF FIGURES


Figure page

2-1. Reference linear translation curve was given by NCAT software...........................11

3-1. Eight cardiac beating frames of a NCAT phantom with beating heart motion
only. (A) Frame 1 to 8 on trans-axial orientation. (B) Frame 1 to 8 on coronal
orientation..................................... .......................... .... ..... ........ 21

3-2. Beating heart frame 1 at different respiratory phases. (A) Trans-axial orientation.
(B ) Coronal orientation. ................................................ ............................... 21

3-3. True Beating Heart Motion generated by NCAT software. Column (A): end-
diastole; Colum n (B): end-systole ......................................................................... 22

3-4. Projection data simulated from OHM. ........... ............................... ...............23

3-5. Projection data simulated from ORM A. ............ ................................................ 23

3-6. Projection dataset ORMU showing linear axial movement of the heart due to
respiratory motion in different angles. (A) Projection angle 1, (B) Projection
angle 38, (C) Projection angle 60 .................................................. ...... ......... 23

3-7. Projection dataset ORMU showing no trans-axial movement of the heart due to
respiratory motion in different angles. (A) Projection angle 1, (B) Projection
angle 38, (C) Projection angle 60. .............................. ... ..... ................ 24

3-8. The relative center of mass for each beating heart frame over one respiratory
c y c le ............................................................................. 2 4

3-9. Measured movement due to respiratory motion of each respiratory phase for
different gate (phase) schema for a 5-second entire respiratory cycle. (A) 5
Phases. (B) 10 Phases. (C) 20 Phases.................................. ......................... 25

3-10. Projection datasets OHM, ORMA, and ORMC. (A) OHM. (B) ORMA. (C)
O R M C ...............................................................................26

3-11. Reconstructed images by RM algorithm and profile plot from projection datasets
OHM, ORMA, and ORMC. (trans-axial orientation) (A) OHM. (B) ORMA. (C)
O R M C (D ) Profiles ............................................ ................. .. ...... 27









3-12. Reconstructed images by RM algorithm and profile plot from projection datasets
OHM, ORMA, and ORMC. (coronal orientation) (A) OHM. (B) ORMA. (C)
O R M C (D ) Profiles ............................................ ................. .. ...... 28

3-13. Reconstructed images by RM algorithm and profile plot from projection datasets
OHM, ORMA, and ORMC. sagittall orientation) (A) OHM. (B) ORMA. (C)
O R M C (D ) Profiles ............................................ ................. .. ...... 29

3-14. Segmentation map of myocardial region for image quality and motion
estim ation accuracy evaluation ........................................... ......................... 30

3-15. SSEs and PMEs as a function of number of respiratory phases..............................31















Abstract of Thesis Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Master of Science

CORRECTION OF THE RESPIRATORY MOTION OF THE HEART
FROM 4D MYOCARDIAL IMAGES

By

Jing Sun

August 2006
Chair: David Gilland
Major Department: Nuclear and Radiological Engineering

During myocardial nuclear medicine imaging, respiratory motion blurs the

projection data and thus causes diagnostic artifacts. A method is developed and

evaluated in this work to compensate for effects of respiratory motion in four

dimensional (4D) gated cardiac Emission Computed Tomography (ECT). Non-uniform-

rational-B-splines-based CArdiac Torso (NCAT) phantoms were generated to simulate

temporal frames regularly spaced during combined respiratory and beating heart cycles.

Projection data were acquired at 60 positions over 180 degrees to model gated

myocardial perfusion Single Photon Emission Computed Tomography (SPECT) with

Technetium-99m (Tc-99m) sestamibi. A 2-cm respiratory rigid body translation was

simulated for the axial heart movement during normal tidal breathing. The motion was

tracked and estimated by calculation of the axial center-of-mass (aCOM). Respiratory

motion correction was applied to the simulated projection data. The 4D image

reconstruction and 3D cardiac wall motion estimation for all temporal frames were









performed simultaneously by a Joint 4D Reconstruction and Motion estimation method

(RM). The Sum of Squared Errors (SSE) was calculated to evaluate the image quality,

and the Phantom-matching Motion Errors (PME) was computed to evaluate the wall

motion estimation accuracy.

After respiratory motion compensation, the image quality, such as spatial

resolution, of corrected images was better than that of respiratory blurred images, and

very close to the quality of the images obtained when no respiratory motion was present.

The accuracy of cardiac wall motion estimates was also improved with respiratory

motion compensation.














CHAPTER 1
INTRODUCTION

1.1 An Overview

Respiratory motion takes place during myocardial perfusion imaging and image

blurring due to patient respiration is thus unavoidable. This respiratory effect is not

negligible; it degrades the diagnostic accuracy when appropriate respiratory motion

compensation is absent. Many researchers have offered methods to reduce the impact of

respiratory motion on image quality. One of the common methods is to let the patients

hold their breath during imaging [1]. The main disadvantage of this method is that the

patient needs to be highly cooperative in order to obtain accurate images. Also, in this

method, imaging has to be performed repeatedly and during each acquisition, the patient

may not hold his/her breathe at the same level as the previous acquisitions, which will

still lead to image blurring. Using an external device, for e.g. a sensor, to set up a

respiratory gating system is an alternate method to reduce the respiratory motion effect

[2-4]. Respiratory motion is not stable or perfectly repeatable like beating heart motion.

Therefore, it is difficult to have a system that both measures the respiratory cycle

accurately such as electrocardiography (ECG) measurements of the beating heart cycle,

and is acceptable to both patients and technologists. Currently, Center-of-Mass is the

most popular method used to track and estimate respiratory motion. This method

requires only the patient projection data and has been realized in both digital [5, 6] and

analogy adjustments of the axial signal [7, 8].









1.2 Literature Review

There are numerous papers published on blurring due to patient motion and

correction methods. Recent works have paid more attention on Positron Emission

Tomography (PET) and SPECT studies. In 1995, using MRI imaging, Wang et al. [9],

determined how much the heart shifted due to respiratory motion. It was concluded that

during normal tidal breathing the movement of the heart due to respiration is dominated

by superior-inferior (SI) motion, which is linearly related to the SI motion of the

diaphragm. Some gating systems were developed for SPECT and PET imaging to

decrease the respiratory motion artifacts during imaging acquisition, for e.g., Macintosh

platform and LabView environment real time system developed by Klein et al. [3], and

respiratory gated SPECT (RGS) system realized by Cho [10]. Segars et al. [11, 12] have

performed simulation studies to study respiratory motion effects on myocardial SPECT

with the spline-based Mathematical Cardiac Torso (MCAT) phantoms. Respiratory

mechanics were modeled in MCAT phantoms for SPECT imaging, and this work led to

the improved and more realistic 4D NCAT phantom. Feng et al. [13] and Bruyant et al.

[14] have developed several methods to compensate for effects of respiratory motion by

applying center-of-mass algorithms to reconstructed images or projection datasets. They

implemented their method with simulated MCAT phantoms. However, their work

focused on respiratory motion only and the heart wall was assumed to be stationary in

their SPECT myocardial perfusion study.

1.3 Our Method

The goal of this study is to develop a motion compensation method using center-of-

mass algorithm to apply to 4D SPECT projection data simulated with both beating heart

motion and respiratory motion. We implement this method to improve the reconstructed









image quality and cardiac wall motion estimation accuracy, for the images reconstructed

and motions estimated simultaneously by RM algorithm.














CHAPTER 2
RESPIRATORY MOTION COMPENSATION, IMAGE RECONSTRUCTION, AND
CARDIAC WALL MOTION ESTIMATION METHODS

2.1 Introduction

To compensate respiratory motion effects, Non-uniform-rational-B-splines-based

Cardiac Torso (NCAT) phantoms were generated to model the realistic organs in the

torso, and simulate the beating heart motion and respiratory motion [15]. The 2D

projection data of NCAT phantoms were obtained by the Radon transform. Data

acquisition was simulated in SPECT by using a parallel-hole collimator geometry [16].

The respiratory motions were estimated by the axial Center-of-Mass algorithm (aCOM)

[14] for the projection data. The motion compensations were also applied to the 4D

projections. The image reconstruction and 3D cardiac wall motion estimation for all

temporal frames were performed simultaneously by a joint 4D Reconstruction and

Motion estimation method (RM) [16]. The image quality before and after compensation

was evaluated by calculations of the Sum of Squared Errors (SSE), and the myocardial

wall motion estimation accuracy was evaluated by computations of the Phantom-

matching Motion Errors (PME) [17].

2.2 NCAT Phantom Simulation

The NCAT phantoms and the cardiac wall motion fields were generated by NCAT

C language software developed in John Hopkins University. Motion vectors were

calculated and assigned in MATLAB to each voxel in each frame from the NCAT









software outputs. Projection datasets were also computed by using a linear interpolation

projector in MATLAB.

2.2.1 Myocardial Phantom Simulation

The NCAT phantom provides a realistic source distribution that models gated

myocardial perfusion imaging with Tc-99m sestamibi by using the Non-uniform Rational

B-Splines (NURBS) method to define the anatomy of a human body at a series of phases

over cardiac cycles and respiratory cycles. This phantom can be used to study the

development of image acquisition strategies, image processing and reconstruction

algorithms and compensation methods. It can be used to research the effects of

anatomical variations and patient motions, such as the beating heart and respiratory

motions, on medical images. NCAT phantom simulation software was developed at the

Johns Hopkins University and fully described by Segars et al. [11, 12].

In this work, a 40 temporal frames NCAT phantom was generated over one entire

respiratory cycle with both beating heart motion and respiratory motion by the NCAT

software. The first frame corresponded to the full exhalation of the respiratory cycle and

the end-diastole of the beating heart cycle. The length of respiratory cycle was set to 5

seconds as a normal tidal breathing. The respiratory cycle started from full exhalation. At

0.455, almost in the middle, of the entire respiratory cycle, motion reached the full

inhalation. The length of cardiac cycle was set to 1 second as a normal male heart

beating. Therefore, the phantom included one completed respiratory motion cycle and

five completed beating heart motion cycles. For each cardiac cycle, there are 8 temporal

frames. Frame 1 was the end-diastole, and frame 4 was the end-systole. Due to the

respiratory motion, the diaphragm motion and chest Anterior-Posterior (AP) expansion

were also simulated. The extent of the maximum diaphragm motion was set to 2









centimeters. This motion was simulated as a translation motion along the body axis. The

extent of the maximum chest AP expansion was set to 1.2 centimeters. This motion only

effect ribcage, body and lungs. Since heart is the only organ included in the NCAT

myocardial phantom, there was no motion along trans-axial direction modeled in this

present study. The projection data will show the heart movement due to respiratory

motion in chapter 3. For the heart, the length of the left ventricle myocardium was 9.43

centimeters, and the average radius of the left ventricle was 2.97 centimeters. The

volume of left ventricle chamber at the end-diastole was 108.36 ml, and at the end-systole

was 42.21 ml. So for this simulated NCAT phantom, the stroke volume (the amount of

blood pumped by the left ventricle of the heart in one contraction) was 66.15 ml and the

ejection fraction (the fraction of blood pumped out of a ventricle with each heart beat)

was 61%. For each frame, 50 slices of a male torso in supine position were created.

Within each slice, the dimension of the phantom was 96 x 96. The voxel size was 0.30

centimeters. The phantom provides a realistic source distribution that models gated

myocardial perfusion imaging with Tc-99m sestamibi. Activities within myocardial

region (including left ventricle myocardium, right ventricle myocardium, left atrium

myocardium, and right atrium myocardium) were assigned to 75 source intensity per cm3

for Tc-99m, energy in 140 keV. For this study, only myocardial source activity was

included. Outside activity was 0 everywhere. The parameters determined for beating

heart cycle were modeled based on tagged MR images obtained from Johns Hopkins

University [15]. And the parameters determined for respiratory cycle were modeled

based on the Visible Human Project CT dataset from the National Library of Medicine

[15].









In practice, any continuous function is instantaneous. In order to be more realistic,

for each time bin, the time interval for each temporal frame, four consecutive frames

were created and averaged over the time bin. In this work, the study was concentrated on

the myocardial region. Also to smooth the source distribution, for each frame, the 3D

cardiac phantom (96 x 96 x 50) was obtained by linearly collapsing the original torso

phantom, which is generated by NCAT software directly with dimension of 256 x 256 x

168 voxels for computational efficiency, and was truncated into dimension 96 x 96 x 50,

which included myocardial region, for computation efficiency purpose. Slice 15 to slice

64 was considered and chosen as myocardial region in this present study. And in each

slice, pixel 17 to pixel 112 on both two dimensions was used as the final myocardial

phantom.

As a reference, a phantom with 8 gated frames over 1-second cardiac cycle with

only beating heart motion was generated, with the same size 96 x 96 x 50 voxels for each

frame. A segmentation file with intensity equals to "1" inside myocardium and intensity

"0" outside was also generated as a map, which defined the heart region by setting a

threshold along the myocardial boundary.

The source distributed objects were created as noise-free phantoms. Attenuation

effects were not considered for all phantoms generated using the NCAT software in this

present study.

2.2.2 Cardiac Wall Motion Simulation

Myocardial wall motion can be computed by NCAT software [15]. A 3D motion

field between each two successive temporal frames over one complete beating heart cycle

was calculated. For each beating heart cycle, one NCAT software motion calculation

generated 8 separate output files if the cardiac cycle was gated into 8 frames. Each









output file included starting points and ending points for and only for those non-zero

voxels in the first frame, but not for each non-zero voxel in that frame. Furthermore,

more than one motion vector might be assigned to the same non-zero voxel. In this case,

the last assigned vector would be selected as the motion vector for this voxel when

motion field between these two frames was computed from the output file. For each

motion vector of a voxel, the starting point in the starting frame was an integer, and the

ending point in the next frame was a real number. The norm of a motion vector was

calculated as the distance between the starting point and the ending point. The direction

of a motion vector was expressed as three components, namely u(r), v(r), and w(r) along

axes x, y, and z, respectively, where r = (x,y,z) denoted the voxel coordinates in the

spatial domain. Those motion vectors were considered true motions for each voxel from

one time bin to another. Output from the software was processed in MATLAB to

calculate the norm and the three components for each motion vector of the voxels. When

evaluated the motion fields estimated by RM algorithm, those motion vectors were

considered as a reference set. Therefore, the NCAT software-calculated motion fields

were named True motion in this present study.

2.2.3 Projection Simulation

Forward projections were performed at 60 angles over 180 degrees, which starts

from patient left hand side, went in front of the patient, and to the patient right hand side.

The 2D projection data were computed by using a linear interpolation projector for each

slice per angle per frame. Effects of attenuation, scatter, and random were not considered

during projection data acquisition.

Clinically, 2000 counts in the myocardial region per frame per projection angle are

realistic for gated Tc-99m sestamibi imaging [18]. All projection data were scaled to this









level. Therefore, the maximum voxel intensity was assigned to 0.237. Poisson noise was

added to all of the noise-free projection data after scaling.

Four projection datasets were computed in this study. The first dataset was

calculated from the 8-frame beating heart motion only phantom in dimension 96 x 50 x

60 x 8. The dataset was named Olject i/th beating Heart Motion only (OHM) in this

present study. The second dataset was computed from the 40-frame phantom directly, in

dimension 96 x 50 x 60 x 40. This projection was calculated as original raw data for

respiratory motion compensation, and was named Olject i/th Respiratory Motion but

Uncorrected (ORMU) in this present study. A third dataset was simulated as imaging

with both beating heart motion and respiratory motion, in dimensions of 96 x 50 x 60 x

8. The 40-frame both motions phantom included one complete respiratory cycle and five

beating heart cycles at the same time. Five cardiac cycles divided the respiratory cycle

into five phases. For example, five phases for each respiratory cycle could be chosen as a

respiratory gate schema. In this case, there were 8 temporal frames for each beating heart

cycle over one respiratory phase. Each temporal frame in one beating heart cycle was

summed over respiratory phases. Then an averaged 8-frame projection dataset was

obtained. It was blurred due to the respiratory motion. The dataset was named Object

' i/th Respiratory Motion and Averaged (ORMA) in this present study. The last dataset

was computed from 40-frame both motion phantom, too. The projection dimension was

also as 96 x 50 x 60 x 8. Averaging was done for each temporal frame in each beating

heart cycle over five respiratory phases, but after applying aCOM respiratory motion

compensation method. The dataset was named Olject i th Respiratory Motion but









Corrected (ORMC) in this present study. All three projection datasets in dimension 96 x

50 x 60 x 8 were in the same noise level.

2.3 Axial Center-of-Mass Motion (aCOM) Estimation and Respiratory Motion
Compensation Method

The Center-of-Mass (COM) of each projection data per temporal frame was

estimated in MATLAB. By calculating the COM of the integrated projection data per

respiratory phase simulated from ORMU, the movement due to respiratory motion alone

body axis can be determined. To compensate for respiratory effects, each frame was

shifted respectively depending on the movement in each respiratory phase to reduce the

blurring in each beating heart gate.

2.3.1 Axial Center-of-Mass (aCOM) Calculation

According to Feng et al. [13], the COM of each frame had the same rigid-body

motion as the object in that frame. The phantom was simulated a 2-centimeter translation

motion along the body axis. From full exhalation to full inhalation, the heart moved

down linearly with diaphragm up to 2 centimeters. The reference linear translation curve

was given by NCAT software, and is shown in Figure 2-1. Example motion curve was

given for 20 respiratory phases. COM of projection data computed from ORMUfor each

frame was calculated to track and estimate respiratory motion. To reduce the noise, the

projection dataset had been integrated over the trans-axial direction per projection angle,

over all of the projection angles, over each beating heart cycle, and over each respiratory

phase. The center of mass along the body axis was determined for any respiratory phase

as given by formula Equation (2-1):









K N T M
Si y p(i, j, o, t, f)
aCOM(f) j= -=l t=l- (2-1)
Z Z p(i, j, t f)
]=1 8=1 t=l 1=1

where M was the total number of voxels in the trans-axial dimension, N was total number

of projection angles, T was total number of beating heart gates, K was the total number of

voxels in axial dimension, and p (i, j, 0, t, f) was the activity in pixel (i, j) in the

projection acquired at angle 0 for beating heart time frame t in respiratory phase f


j25
E
E 20 -

15 -

S10 -

05
10


0 0

0 5 10 15 20
Respiratory Phases

Figure 2-1. Reference linear translation curve was given by NCAT software

Assume within each respiratory phase, that the respiratory motion was constant.

Therefore, five aCOM were calculated finally for five phases over one complete

respiratory cycle. Each one indicated the axial shift magnitude for the beating heart

temporal frames within that respiratory phase. The estimated relative center of mass for

each respiratory phase was in unit number of voxels. The voxel size of projection data

was 0.3 centimeters.









2.3.2 Respiratory Motion Compensation

Use equation (2-1) to calculate the COM of the projection dataset for each beating

heart frame. COMs of 40 temporal frames in total were computed for the five cardiac

cycles, as well one entire respiratory cycle.

For motion compensation, each frame of projection dataset ORMUwas vertically

shifted because only axial movement was considered in this study. Full exhalation in the

respiratory cycle and end-diastole in the beating heart cycle was defined as the reference

status with zero relative movement in this work. Therefore, the last respiratory phase was

chosen as the reference phase, for which the projection data were not shifted. All other

projections were shifted up. Based on the assumption that respiratory motion was

constant within each respiratory phase, the beating heart frame within each respiratory

phase will shift same amount. And the magnitude of the shifting was depended on the

movement calculated by Equation (2-1) of the phase which the beating heart frame was

located. Because the axial shift magnitude was not an integer number of voxels, a linear

interpolation was performed along the body axis to obtain an accurate shifting.

After respiratory motion compensation, the COM of each shifted beating heart

frames were closed to the COM of the reference phase. The magnitude of the errors

depended on how many respiratory phases were chosen. Because, in practice, there are

no very accurate devices such as ECG monitors used for cardiac gated imaging, the

narrow respiratory gates could not be obtained. Five respiratory phases were chosen

because of the simplicity, so that each respiratory phase was reasonably assigned as 1

second. Ten phases and twenty phases respiratory schema were also simulated for

optimization and comparison. The aCOM method was also applied to the corrected









projection dataset to check the constancy of the aCOM position during the respiratory

cycle. The completed results are listed and shown in Chapter 3.

Conclusively, the motion compensation algorithm can be described as following

steps:

1. Calculate the center-of-mass of each respiratory phase.

2. Choose the respiratory phase which has least movement (the phase has smallest
center-of-mass value) as the reference phase, and calculate the movements of all
other phases by giving relative to the reference phase.

3. Shift all beating heart frames in each respiratory phase axially according to the
relative movement calculated for that phase.

4. Sum each beating heart frame over all of the respiratory phases.

2.4 Image Reconstruction and Cardiac Wall Motion Estimation (RM) Method

Projections simulated from OHM, ORMA and ORMC were all reconstructed by

image reconstruction and cardiac wall motion estimation (RM) algorithm. This algorithm

was first discussed by Mair et al. [19] and developed in MATLAB code by Cao et al.

[20]. The method was later improved and implemented in Fortran code by Gilland et al.

[21] and Mair et al. [17].

2.4.1 RM Algorithm

RM algorithm was a method for gated cardiac ECT that reconstructed the voxel

intensities of the gated images and estimated the 3 dimensional motion of the cardiac wall

for all frames simultaneously. The algorithm consisted of two steps, the image

reconstruction step (R step) and the cardiac wall motion estimation step (M step). It was

based on an iterative algorithm for minimizing an objective function that included the

negative log likelihood of the data, the standard optical brightness constraint, and a

biomechanical model for elastic deformations of the myocardium.









In the early version, Mair et al. [19] discussed the estimation of motion and images

for only two frames with maximum deformation, end-diastole and end-systole. During R

step, the motions were fixed and the objective function was minimized over the frames

iteratively. During M step, the frames were fixed and the objective function was

minimized over the motion field iteratively based on solving the Euler-Lagrange

equations [22, 23]. The algorithm was implemented in MATLAB because of its

extensive library of mathematical functions. However, it was rewritten in Fortran

language by Gilland et al. [21] to enhance the computational efficiency. And the most

important point was that, in Mair et al. [17] the algorithm was extended from two frames

to a complete cardiac cycle frames, which can be any finite number. Eight frames were

chosen in this simulation study. Instead of being performed in pairs of frames, a

mathematical framework, which included the entire cycle of images and motion fields,

was developed to make only one single objective function. In this manner, all frames

were considered as a loop. Each update of a frame involved the current estimate of the

frame immediately before and after this frame. The objective function minimization was

performed directly by using a type of sequential quadratic programming with the

conjugate gradient algorithm [17] applied at each step. After all of the improvements and

optimizations, the algorithm guaranteed non-negativity of the reconstructed images, and

each iteration decreased the value of the objective function. The simultaneous estimation

allowed the data to influence both image reconstruction and cardiac wall motion

estimation at the same time [17].

2.4.2 Objective Function

In this study, the myocardium was considered as a deformable elastic material

satisfying the standard constitutive relations in continuum mechanics [24, 25].









2.4.2.1 Negative log likelihood term

It was assumed that the small deformation between two successive frames would

not result in very large strain energy [23]. If fi(r), f2(r), fT(r) and mi(r), m2(r), .

mT(r) were denoted T frames image and T motion fields between each successive two

frames respectively, the important non-negative condition on each intensity function ft

will be denoted by f> 0 [17]. Then the negative log likelihood term of an objective

function was defined as Equation (2-2) [17]:


L(f)= 1 Hf, (i g (i)logHf,(i)], (2-2)
t=1 1=1

where M was the total number of the detector bins, the dimension of the phantom, H was

the projection operator, and gi(t), g2(t), gT(t) were gated datasets.

2.4.2.2 Image matching term

When the image intensities between frames t and t + 1 perfectly matched the

motion field mi between them, ft(r) should equal to ft+l(r+mt(r)) for all r [17]. In practice,

image frames and motion fields were never exactly matched. Therefore, the image

matching term of an objective function between frames t and t + 1, which defined as

Equation (2-3) [17], was expected to be small:


E (f, f,, ; m,) = [f (r) f, (r + m, (r))dr, (2-3)

where r was the voxel coordinate, and r = (x, y, z). If considered the entire cardiac cycle,

the total image matching term was defined as Equation (2-4) [17]:

T
ef(f;m)= _E, (f,f+, mi;). (2-4)
t=1









2.4.2.3 Strain energy term

By treated myocardium as a deformable elastic material, material strain resulting

from the cardiac motion can be modeled by the strain energy [17]. The strain energy

depended on the derivatives of the motion vector and material parameters X and [,

namely Lame Constants [15]. After several experiments, it was found that the material

difference can be ignored in this study. Next, X and [t were set equal to one both inside

and outside the myocardium so that Young's modulus was t(3k + 2[)/(k + t) = 2.5 and

Poisson's ratio was 0.5 /( X + [t) = 0.25 respectively In this study. Segmentation process

by Mair et al. [17] was still included in this work, but only applied for evaluation method

for image quality and motion accuracy presented in Chapter 3. For each motion field

between frame t and t + 1, mt had three components ut, vt, and wt. By applying the strain

energy function in Love's classical book [15], the strain energy term of an objective

function was defined as Equation (2-5) [17]:

T
Es (m) Es (m),where
t (2-5)
Es = A(r)(, u t + vt, + )2dr + J /(r)fu + v + w )dr


+--I()U
2+ pr u, \ + v )2 + 1,, + w, )2 + (^, + )r wr


where utx, vt,y, wtx, were continuous first partial. e.g. utx X and [t were
ax

included here for general expression. They were assigned to ones in this study. Finally,

the total objective function was defined as Equation (2-6) [17]:

E(f; m) = aL(f) + E, (f; m)+ Pfss (m), (2-6)









where a and p were hyper-parameters, which reflected the influence that image matching,

projection data and strain energy have on the estimations. Upon several experiments by

Mair et al. [17], the value of a and p were determined as 0.015 and 0.0035 respectively in

this study.

If fl*(r), f2*(r), ... fT*(r) and ml*(r), m2*(r), ... mT*(r) were denoted true images and

true motion fields respectively, minimized function F resulted as Equation (2-7) [17]:

(f*; m*) = argmin (f;m). (2-7)
f >,m

In RM algorithm, the initials were uniformed images and zero motion fields. The

reconstructed images were updated first while motion fields were fixed (R Step). Then

the estimated cardiac wall motions were updated while image frames were fixed (M

Step). Within R Step or M Step, the algorithm went out of loop when it was converged.

After one R Step and one M Step, each image and motion field was updated once equally.

In R step, image matching and negative log likelihood terms of the objective function

were used to improve the reconstructed images. The expected estimate for each iterative

R step was defined as Equation (2-8) [17]:

f" + = argmin (f; m n) (2-8)
f>0

In M step, updated frames in latest R step were used to improve estimated motions. The

expected estimate for each iterative M step was defined as Equation (2-9) [17]:

m"+1 = arg min E(f (+ m) (2-9)
m

Experimental results of the reconstructed images and estimated cardiac wall

motions are presented in Chapter 3.









2.4.3 Reconstruction and Motion Estimation

100-iteration RM algorithm, with 100-iteration R step and 100-iteration M

alternately, was applied to all projection datasets ORMA, ORMU, and ORMC. Fortran

code was running on Penguin Tempest 2100 (Dual-processor 64-bit AMD Opteron, 8 GB

physical memory). It took about 4 hours to finish 100 iterations. The value of a and 0

were assigned as 0.015 and 0.0035 respectively according to Mair et al. [17].

2.5 Image Quality and Cardiac Wall Motion Accuracy Evaluation Methods

The image quality of the reconstructed image was evaluated using Sum of Squared

Errors (SSE) method. The accuracy of the cardiac wall motion estimation was evaluated

using images matching term of the objective function. For the evaluation, phantom OHM

was using as the reference, and the image matching equation was also named Phantom-

matching Motion Error (PME) method.

2.5.1 Image Quality Evaluation

For image quality evaluation, SSE method was defined as Equation (2-10) [17]:


SSE(f, f) = f, (r)- f (r) (2-10)
t=1 r

where f (r), f (r),... f (r) and f (r), f (r),... f ,(r) were phantom source distributed

object and reconstructed image, respectively. Only the voxels inside the myocardial

region would be summed. A segmentation map was generated to define the moycardium

by given a threshold intensity. In this study, only heart was included in the phantom.

Therefore, the intensity threshold was assigned as zero. In another word, any non-zero

intensity voxel was considered within myocardium. The segmentation map and the









experimental results were shown in Chapter 3. SSE evaluation method was also used for

a value optimization study [17].

2.5.2 Cardiac Wall Motion Evaluation

Recall image matching equation:


E, (f,, f,; mn,) = [f, (r) fl (r + m, (r))1dr. (2-3)

For cardiac wall motion evaluation, same equation was used but called PME

method defined as Equation (2-11) [17]:


PME(m)= Y ft(r)- f +(r+m(r)) (2-11)
t=1 r

A
where f,(r) was the phantom source distributed object at frame t, and mt(r) was the

estimated motion field from frame t to frame t + 1. If frame t is the last frame in a

sequence (e.g. the 8th frame in a 8-gate cardiac image cycle), frame t + 1 would be

assigned to be the first frame in the cycle. In this evaluation study, only those non-zero

intensity voxels were considered within myocardium and were computed for PME value.

The segmentation map was used, too. The experimental results were shown in Chapter 3.

PME evaluation method was used for P value optimization study, too [17].














CHAPTER 3
EXPERIMENTAL RESULTS AND DISCUSSION

3.1 Simulated Phantoms, True Motion Field, and Projections

Experimental results and analysis are presented in this Chapter. Calculated relative

center-of-mass of projections before and after respiratory compensation, reconstructed

images, estimated cardiac wall motion, and SSE and PME evaluation data are illustrated

in the following sections.

3.1.1 Cardiac Gated NCAT Phantom

3.1.1.1 Reference phantom with beating heart motion only

As a reference phantom, eight cardiac beating frames are shown in figure 3-1 in

trans-axial and coronal orientations of a phantom with beating heart motion only. The

dimension of the phantom was 96 x 96 x 50 for each gating frame. End-diastole and

end-systole can be seen clearly at frame 1 and frame 4, respectively.

3.1.1.2 Phantom with both beating heart motion and respiratory motion

To look at whether respiratory motion had an effect on myocardial imaging, a

phantom with combined beating heart motion and respiratory motion is shown in Figure

3-2. If the length of the respiratory cycle was assigned as 5 seconds, and the length of the

beating heart cycle was assigned as 1 second, there were 5 complete beating heart cycles

during one respiratory cycle. The first beating heart frame for each one of the five

beating heart cycles at different respiratory phases are displayed in Figure 3-2.

Myocardium shifted down about 2 centimeters from full exhalation to full inhalation with

the diaphragm.
























(B)

Figure 3-1. Eight cardiac beating frames of a NCAT phantom with beating heart motion
only. (A) Frame 1 to 8 on trans-axial orientation. (B) Frame 1 to 8 on coronal
orientation.







(A)




(B)
Figure 3-2. Beating heart frame 1 at different respiratory phases.
(A) Trans-axial orientation. (B) Coronal orientation.

3.1.2 True Motion Field of the NCAT Phantom

The true motion field for the beating heart motion only phantom was superimposed

on the reference object and is plotted in Figure 3-3. Frames in the left column showed

the end-diastole. All motion vectors were pointing to the direction where the heart

contracted to during systole. Frames in the right column show the end-systole. All

motion vectors point to the direction where the heart enlarged to during diastole.





































Figure 3-3. True Beating Heart Motion generated by NCAT software.
Column (A): end-diastole; Column (B): end-systole.

3.1.3 Projection Datasets

Projection dataset OHMis shown in Figure 3-4. Dataset ORMA is shown in Figure

3-5. Only projection angle 38 of each dataset is shown in these two figures. It was

observed that cardiac gated projection data were blurred by the respiratory motion

obviously. Hence, the respiratory motion effect could not be ignored during the SPECT

imaging. Figure 3-6 and Figure 3-7 show the first frame for each beating heart cycle over

respiratory phases (5 phases shown in this example). Projection data angle 1, 38 and 60

are shown in these two figures. Linear axial movement can be seen vertically from

Figure 3-6. Purple line indicates the heart position without respiratory effect. No trans-








23

axial movement horizontally is found as expected how the phantom been modeled from

Figure 3-7.







Figure 3-4. Projection data simulated from OHM.


Figure 3-5.


-d from ORMA.


I /- I


Figure 3-6. Projection dataset ORMU showing linear axial movement of the heart due to
respiratory motion in different angles.
(A) Projection angle 1, (B) Projection angle 38, (C) Projection angle 60.

3.2 aCOM Calculation and Respiratory Motion Compensation Method

The relative COM along axial direction of projection datasets ORMU and ORMC

were presented in the following sections. The COM of each respiratory phase over the

cycle was also shown here.

3.2.1 aCOM Calculation

Figure 3-8 illustrated the relative COM of projection dataset ORMU along axial



























(A) (B) (C)
Figure 3-7. Projection dataset ORMU showing no trans-axial movement of the heart due
to respiratory motion in different angles.
(A) Projection angle 1, (B) Projection angle 38, (C) Projection angle 60.

direction calculated by using Equation (2-1) for each cardiac frame over one entire

respiratory cycle. Compare to the 2-centimeter linear translation due to the respiratory

motion, the magnitude of the movement due to the heart beating (average about 0.15

centimeters) was unobservable. Therefore, only the displacement caused by the

respiratory motion could be seen from the chart.


8
7o 7*

E X






0 **
1 6 11 16 21 26 31 36
Beating Heart Frames over one Respiratory Cycle

Figure 3-8. The relative center of mass for each beating heart frame over one respiratory
cycle










3.2.2 Respiratory Motion Compensation

Five phases (gates) over one 5-second entire respire cycle was calculated first

because of the simplicity. Ten phases and twenty phases respiratory schema were also

simulated for comparison and optimization. All three gate schema were shown in Figure

3-9. The height of each column indicates how much the beating heart frames within that

respiratory phase need to be shifted.





x 5

0' 3
o .
I>

E -1 1 2 3 4 5

Respiratory Phases

(A)

7
-@


0
El a E
o 1
E 1 2 3 4 5 6 7 8 9 10
Respiratory Phases


(B)










(C)
o >



E E -1 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
Respiratory Phases

(C)

Figure 3-9. Measured movement due to respiratory motion of each respiratory phase for
different gate (phase) schema for a 5-second entire respiratory cycle.
(A) 5 Phases. (B) 10 Phases. (C) 20 Phases.







26
After respiratory motion compensation, the relative COM for each beating heart

frame should be very close to the reference phase. The ideal illustration should be a

straight horizontal line. The error depends on which gate schema was used. The more

gates used for each respiratory cycle, the less errors were. Twenty phases schema would

be chosen according to the clinical realization by Cho et al. [10]. The projection datasets

of projection angle 38 after (ORMC) 20-phase respiratory motion correction is shown in

Figure 3-10. Respiratory motion blurred projection data (ORMA) and heart motion only

projections (OHM) were shown in the same figure for comparison. The projection data

can be seen from the images clearly that it was very close to the reference projections

(OHM) after motion compensation. In practice, the respiratory motion might not be very

stable. The length of the respiratory cycle might change during projection data

acquisition. Over sampling of the gates and averaging over each chosen respiratory

phase would contribute to the motion compensation improvements.








(A) (B) (C)

Figure 3-10. Projection datasets OHM, ORMA, and ORMC.
(A) OHM. (B) ORMA. (C) ORMC.

3.3 Image Reconstruction, Cardiac Wall Motion Estimation,
and SSE and PME Evaluations.

The RM reconstructed images from all three projection datasets (OHM, ORMA, and

ORMC) are shown in this section. SSE and PME evaluations are listed in the tables here.








27

3.3.1 Image Reconstruction and Wall Motion Estimation

The reconstructed images and profile plots from all three projections (OHM,

ORMA, and ORMC) for three orthogonal slice orientations trans-axial, coronal, and

sagittal are shown in Figure 3-11, Figure 3-12 and Figure 3-13, respectively. 20

respiratory phases schema were chosen to show result in this study according to the

clinical level determined by Cho et al. [10].


(A)


(B) (C)


--OHM mages
- ORMA images
-ORMC images


0.2
o


0.1
u,
0-


0


10 20 30 40 50
Voxels


(D)

Figure 3-11. Reconstructed images by RM algorithm and profile plot from projection
datasets OHM, ORMA, and ORMC. (trans-axial orientation)
(A) OHM. (B) ORMA. (C) ORMC. (D) Profiles


















(A) (B)


0.3 OHM images
ORMA images
ORMC Images


0.2
o



a 0.1



0
0 5 10 15 20 25 30 35 40
Voxels


(D)

Figure 3-12. Reconstructed images by RM algorithm and profile plot from projection
datasets OHM, ORMA, and ORMC. (coronal orientation)
(A) OHM. (B) ORMA. (C) ORMC. (D) Profiles



















(A) (B)


0.3654 OHM images
ORMA images
-ORMC images

5 0.2436
o


a 0.1218



0
0 10 20 30 40
Voxels

(D)
Figure 3-13. Reconstructed images by RM algorithm and profile plot from projection
datasets OHM, ORMA, and ORMC. sagittall orientation)
(A) OHM. (B) ORMA. (C) ORMC. (D) Profiles

From the reconstructed images in all orientations, it was observed that the blurring

due to the respiratory motion had been reduced greatly by motion compensation method.

The image reconstructed from corrected projections was much closer to the image

reconstructed from the reference dataset. ORMC images looked smoother than true

(OHM) images. That was because ORMC images were summed over respiratory phases

after motion compensation. From the profile plots, it was observed that both OHM

images and ORMC images peak were narrow and had higher peak values. It

demonstrated that the noise control was improved, and the less blurring on the edge of the







30

myocardial region could be seen, too. Estimated wall motion would be evaluated by

PME method explained in next section.

3.3.2 SSE and PME Evaluation Method

The segmentation map used for SSE and PME evaluation study is shown in Figure

3-14. Slice 27 in frame 1 was chosen to show here.











Figure 3-14. Segmentation map of myocardial region for image quality and motion
estimation accuracy evaluation

The SSE evaluation results of OHM, ORMA, and ORMC images calculated by

Equation (2-10) are listed in Table 3-1. The PME evaluation results of OHM, ORMA,

and ORMC motion fields calculated by Equation (2-11) are listed in Table 3-2. The SSEs

and PMEs values as a function of number of respiratory phases are shown in Figure 3-15.

Table 3-1. The SSE results of OHM, ORMA, and ORMC images quality evaluation.
Sum of Squared Errors (SSE)
True Object 0
OHMimage 472.48
ORMA image 1503.93
ORMC image (5 phases) 643.57
ORMC image (10 phases) 560.54
ORMC image (20 phases) 521.79








31

Table 3-2. The PME results of OHM, ORMA, and ORMC
motion estimation accuracy evaluation.
Phantom-matching Motion Errors (PME)
True motion 88.25
OHMmotion 183.75
ORMA motion 386.35
ORMC motion (5 phases) 307.37
ORMC motion (10 phases) 290.85
ORMC motion (20 phases) 272.91

From the SSEs and PMEs table, respiratory motion compensation method

demonstrated improved SSE error as well as improved image noise and spatial resolution

characteristics. Both SSE and PME would decrease with the increased number of

respiratory phases. Also according to the clinical study by Cho et al. [10], 20 phases

would be an ideal choice for respiratory-gated study.


, 800
LU
I 600
= 400
w 200
co 0


A SSEs
SPMEs


*0


0 5 10 15 20 25
Number of Respiratory Phases

Figure 3-15. SSEs and PMEs as a function of number of respiratory phases.














CHAPTER 4
CONCLUSIONS AND FUTURE WORK

4.1 Summary and Conclusions

A respiratory motion compensation method was discussed and developed in this

study. Method was implemented by a NCAT phantom simulation with both beating heart

motion and respiratory motion. Respiratory motion was tracked and estimated by aCOM

algorithm. Motion compensation was realized by applying method on to the projections

computed from the NACT phantoms with a linear interpolation operator. The true

projection, averaged projection and the corrected projection datasets were reconstructed

by RM algorithm. The cardiac wall motion was estimated simultaneously with the image

reconstruction. The respiratory motion compensation method was evaluated by

comparison between OHM, ORMA and ORMC images and motion fields. Numerical

results were explained by calculation of SSEs and PMEs errors.

* From the projections, it was observed that cardiac gated projection data were
blurred by the respiratory motion.

* From aCOM calculation, the respiratory motion can be tracked and estimated for
each respiratory phase.

* After respiratory motion compensation by aCOM computations, the corrected
projection data are very close to the data without respiratory motions.

* It can be seen that blurring are reduced on the reconstructed images from the
corrected projection data. From the evaluation results, the values of SSE and PME
demonstrated improved spatial resolution and more accuracy of wall motion
estimation than that without respiratory motion compensation.

* Different respiratory phases' schemas have been tested. More respiratory gating
sampling gives more image quality improvement and more accuracy on wall









motion estimation. According to the clinical patient data study of other papers, 20
respiratory phases are chosen for this study.

4.2 Future Work

The study included only the heart for the NCAT phantom simulation. Since for

gated myocardial perfusion SPECT imaging with Tc-99m, the radiopharmaceuticals will

generate intensities in the other organs, such as liver, almost same level relative to the

heart. That may have big effect on the center-of-mass calculation to track and estimated

the respiratory motion. Future studies using phantoms with more organs can make this

compensation method more general.

In this study, it was assumed that the respiratory motion was constant within one

respiratory phase. In practice, the respiratory gate can not be too narrow. This will bring

errors to the simulation with both beating heart and respiratory motion, because within

one respiratory phase (e.g. less than 5-phase schema), there are more than one beating

heart frames included. One beating heart frame was the time bin which was finally

implemented. Further consideration of more cardiac cycle and respiratory cycle

combination may yield more accurate results. Respiratory motion effect will take place

in almost every myocardial perfusion SPECT imaging examination. In practice, the

respiratory cycle may not be fixed during imaging, and the beating heart motion and

respiratory motion combination could be more complicated. Applying this compensation

method to clinical patient data would be an important future study.















LIST OF REFERENCES


1. Gottschalk A, Harper PV, Jimenez FF, Petasnick JP. Quantification of the
respiratory scanning with the rectilinear focused collimator scanner and the gamma
scintillation camera. Journal ofNuclear Medicine. 1966;7:243-251.

2. Zubal G, Bizais Y, Bennett GW, Brill AB. Dual gated nuclear cardiac images.
IEEE Transactions on Nuclear Science. 1984;31:566-569.

3. Klein GJ, Reutter BW, Ho MH, Reed JH, Huesman RH. Real-time system for
respiratory-cardiac gating in positron tomography. IEEE Transactions on Nuclear
Science. 1998;45(4):2139-2143.

4. Klein GJ, Reutter BW, Huesman RH. Data-driven respiratory gating in list mode
cardiac PET. Journal ofNuclear Medicine. 1999;40(5): 113P.

5. Oppenheim BE. A method using a digital compute for reducing respiratory artifacts
on liver scans made with a camera. Journal ofNuclear Medicine. 1971;12:625-
628.

6. Schmidlin P. Development and comparison of computer methods for organ motion
correction in scintigraphy. Journal ofNuclear Medicine. 1999;40(5): 113P.

7. Hoffer PB, Oppenheim BE, Stterling ML, Yasillo N. A simple device for reducing
motion artifacts in gamma camera images. Radiology. 1972; 103:199-120.

8. Baimel NH, Bronskill MJ. Optimization of analog-circuit motion correction for
liver scintigraphy. Journal ofNuclear Medicine. 1978; 19:1059-1066.

9. Wang Y, Riederer SJ, Ehman RL. Respiratory motion of the heart: kinematics and
the implications for the spatial resolution in coronary imaging. Magnetic
Resonance in Medicine. 1995;33:713-719.

10. Cho K, Kumiata S, Okada S, Kumazaki T. Development of respiratory gated
myocardial SPECT system. Journal ofNuclear Cardiology. 1999;6:20-28.

11. Segars WP, Lalush DS, Tsui BMW. Modeling respiratory mechanics in the MCAT
and spline-based MCAT phantoms. IEEE Transactions on Nuclear Science.
2001;48(1):89-97.









12. Segars WP, Tsui BMW. Study of the Efficacy of Respiratory Gating in Myocardial
SPECT Using the New 4-D NCAT Phantom. IEEE Transactions on Nuclear
Science. 2002;49(3):675-679.

13. Feng B, Bruyant PP, Pretorius PH, Beach RD, Gifford HC, Dey J, Gennert M, King
MA. Estimation of the rigid-body motion from three-dimensional images using a
generalized center-of-mass points approach. IEEE Nuclear Science Symposium
Conference Record. 2005;M07:248.

14. Bruyant PP, King MA, Pretorius PH. Correction of the respiratory motion of the
heart by tracking of the center of mass of thresholded projections: a simulation
study using the dynamic MCAT phantom. IEEE Transactions on Nuclear Science.
2002;49(5):2159-2166.

15. Segars WP. Development of a new dynamic NURB S-based cardiac torso (NCAT)
phantom. Ph.D. dissertation. The University of North Carolina, Chapel Hill, 2001.

16. Cao Z, Simultaneous reconstruction and 4D motion estimation for gated
myocardial emission tomography. Ph.D. dissertation. The University of Florida,
Gainesville, 2003.

17. Mair BA, Gilland DR, Sun J, Estimation of images and non-rigid deformations in
gated emission CT. Transactions on Medical Imaging. Accepted in April, 2006.

18. Narayanan MV, King MA, Soares EJ, Byrne CL, Pretorius PH, Wernick MN.
Application of the Karhunen-Loeve transform to 4D reconstruction of cardiac gated
SPECT images. IEEE Transactions on Nuclear Science. 1999;46(4):1001-1008.

19. Mair BA, Gilland DR, Cao Z. Simultaneous motion estimation and image
reconstruction from gated data. IEEE International Symposium on Biomedical
Imaging, 2002:661-664.

20. Cao Z, Gilland DR, Mair BA, Jaszczak RJ. Three-dimensional motion estimation
with image reconstruction for gated cardiac ECT. IEEE Transactions on Nuclear
Science. 2003;50(3):384-388.

21. Gilland DR, Mair BA, Sun J. Joint 4D reconstruction and motion estimation in
gated cardiac ECT. Intern Conference on Fully 3D Image Reconstruction Radiation
Nuclear Medicine. 2005:303.

22. Klein J, Reutter BW, Huesman RH. Non-rigid summing of gated PET via optical
flow. IEEE Transactions on Nuclear Science. 1997;44(4): 1509-1512.

23. Klein GJ, Huesman RH. Elastic material model mismatch effects in
deformablemotion estimation. IEEE Transactions on Nuclear Science.
2000;47(3):1000-1005.






36


24. Love AEH, A Treatise on the Mathematical Theory of Elasticity. Dover
Publications. New York, 1944.

25. Humphrey JD, Delange SL. An introduction to biomechanics: solids and fluids,
analysis and design. Springer-Verlag. New York, 2004.















BIOGRAPHICAL SKETCH

Jing Sun was born in Jilin, China, on February 28th, 1976. In September 1994, she

began her college education in Tsinghua University in Beijing and obtained her Bachelor

of Science degree in engineering physics in July 1999. In 2002, she obtained her Master

of Science degree in nuclear engineering from North Carolina State University.

In August 2003, Jing was admitted to the University of Florida as a graduate

student to the medical physics program in the Department of Nuclear and Radiological

Engineering. During the last three years, she has worked with Dr. David Gilland in

SPECT myocardial imaging project.

Jing Sun is the only daughter of Zhongzhi Sun and Guirong Liu and is married to

Dapeng Ding.