Computer-aided lesion detection in positron emission tomography: A signal subspace fitting approach

PDF

Download

Open document

Flip pages

Download a page range

Download transcript

Copy asset link

Request this asset

Request accessible transcript

Transcript (if available)

Content INFORMATION TO USERS This manuscript has been reproduced from the microfilm master. UMi films the text directly from the original or copy submitted. Thus, some thesis and dissertation copies are in typewriter face, while others may be from any type of computer printer. The quality of this reproduction is dependent upon the quality of the copy submitted. Broken or indistinct print, colored or poor quality illustrations and photographs, print bleedthrough, substandard margins, and improper alignment can adversely affect reproduction. In the unlikely event that the author did not send UMI a complete manuscript and there are missing pages, these will be noted. Also, if unauthorized copyright material had to be removed, a note will indicate the deletion. Oversize materials (e.g., maps, drawings, charts) are reproduced by sectioning the original, beginning at the upper left-hand comer and continuing from left to right in equal sections with small overlaps. ProQuest Information and Learning 300 North Zeeb Road, Ann Arbor, Ml 48106-1346 USA 800-521-0600 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. C O M P U TE R -A ID E D LESION D E T E C T IO N IN POSITRON EMISSION T O M O G R A P H Y : A SIG NAL SUBSPACE F IT T IN G APPROACH by Chmig-Chioh Huang A Dissertation Presented to the FACULTY OF THE G R A D U A T E SCHOOL U N IV E R S ITY OF SO UTHERN C A LIF O R N IA In Partial Fulfillm ent of the Requirements for the Degree DO CTOR OF PH ILO SO PHY (E LE C T R IC A L EN G IN EER IN G ) December 2001 Copyright 2001 Chung-Chieh Huang Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. UMI Number: 3065797 Copyright 2001 by Huang, Chung-Chieh All rights reserved. UMI UMI Microform 3065797 Copyright 2002 by ProQuest Information and Learning Company. All rights reserved. This microform edition is protected against unauthorized copying under Title 17, United States Code. ProQuest Information and Learning Company 300 North Zeeb Road P.O. Box 1346 Ann Arbor, Ml 48106-1346 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. UNIVERSITY OF SOUTHERN CALIFORNIA The Graduate School University Park LOS ANGELES. CALIFORNIA 90089-1695 Thi s di ssertation, w ritte n b y Chung-Chieh Huang U nder the d ire c tio n o f h. . . ! S. . D issertatio n Com m i ttee, and approved b y a ll its m em bers, has been presented to an d accepted b y The G raduate School , in p a rtia l fu lfillm e n t o f requi rem ents fo r th e degree o f D O C TO R O F P H I LO S O P H Y Dean o f Graduate Studies D ate December 17, 2001 DI SSER TA T IO N C O M M IT T E E Chairperson Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Dedication To my father: Su-Jei Huang (1942-1991) my mother: Pai-Yueh Huang my wife: Li-Hsing Tsai my brother: Chung-Kwen Huang my sisters: Su-Jen Huang, Su-Fen Huang, and Su-Ju Huang my daughter: Vivian Huang my son: Ethan Huang. Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Acknowledgments I would like to express my many thanks to Professor Richard M. Leahy at SIPI and Professor Xiaoli Yu at Radiology Department who both kindly guided and supported my graduate study at I'SC’. I have benefited greatly from Dr. Leahy's rich experiences in the medical imaging and signal processing and Dr. Yu's extensive knowledge in the feature identification and target detection has made the interactions with her very rewarding. I w ill continue to benefit from these invaluable experiences at I'SC’ and extend them through my career life. I also thank to the members of my guidance committee: Dr. Peter Conti. Dr. Antonio Ortega. Dr. Manbir Singh, and Dr. Zhang Zhen for their time and suggestions. The support from the CSC' PET Imaging Center was very im portant to my research, where Dr. Peter Conti. Dr. .James Bading. Jennifer Keppler. Steve Hayles. Michelle Ernsdorff. Monica Garfield. Pamela Ewing, and Jenise Evans all played significant roles to my research. I w ill present my special thanks to Dr. .James Bading for his useful suggestion and friendly encouragement. I w ill also thank to the following people who added positively to the research and my academic experience, under Dr. Yu's research group: Piyapong Thanyasrisung. Chia-Chang Hu. .Jimbo Zheng, and Siwas Chandhrasri: under Dr. Leahy's research iii Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. group: Dr. Lance Hsu. Dr. Jinyi Qi. Dr. B ijan Tim sari. Dr. David Shattuck. John Ernier. Paola Bonetto. Evren Asina. Bing Bai. Alexei Ossadtchi. Karim .Jerbi. Quanzheng Li. and Joaquin Rapela. where special thanks go to Dr. Lance Hsu and Dr. Jinyi Qi for their valuable discussions on PE T imaging and image reconstruction. Finally, most im portantly. I would like to thank to my wife Li-Hsing Tsai for her spiritual support and her carefully taking care of my two little kids. Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Contents Dedication ii Acknowledgments iii List Of Tables viii List Of Figures ix Abstract xiv 1 Introduction 1 1.1 M o tiv a tio n ............................................................................................................ 1 1.1.1 PET Diagnosis and Computer-Aided Lesion Detection . ... 1 1.1.2 Dynamic FD G -P E T S tu d y .................................................................. 3 1.1.3 Data Processing and S im u la tio n ........................................................ 5 2 Basic Principles of PET 7 2.1 O verview ................................................................................................................ 7 2.2 Physics of P E T .................................................................................................. 7 2.2.1 Positron Emission and Photon D e te c tio n ....................................... < 2.2.2 Radiotracers Used in P E T .................................................................. 9 2.3 Statistical Model for PET Im a g in g ............................................................. 12 2.4 PET Image R econstruction.............................................................................. 14 2.4.1 Filtered Backprojection M e th o d ........................................................ 14 2.4.2 Iterative Image Reconstruction M e th o d .......................................... 14 3 PET-FDG Dynamic Data Modeling 17 3.1 O verview ................................................................................................................ 17 3.2 Tracer Kinetic T echnique.................................................................................. IS 3.2.1 Tracer S e lectio n...................................................................................... 19 3.2.2 Model S e le c tio n ...................................................................................... 22 3.2.3 Assumptions in the Tracer Kinetic M o d e lin g ................................ 23 3.3 A 4-k Com partm ental Model for the FDG Tracer K in e t ic s ................... 25 v Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 3.3.1 Homogeneous A s s u m p tio n ................................................................ 25 3.3.2 Heterogeneous A ssu m p tio n ................................................................ 29 4 Dynamic Data Formulation and Subspace Identification 34 4.1 O verview ............................................................................................................... 34 4.2 Dynamic Data Formulation .......................................................................... 34 4.3 Subspace Identification from Know Tissue T y p e ....................................... 36 4.3.1 Xon-Parametric: Singular Value D e co m p o sitio n....................... 36 4.3.2 Parametric: Least Scpiares Estim ation ........................................ 37 4.3.3 Subspace Refining: A Subspace Distance M e a s u re .................... 38 5 Matched Subspace Detector 40 5.1 O verview .............................................................................................................. 40 5.2 Hypothesized Data M o d e l............................................................................. 41 5.2.1 Replacement M o d e l............................................................................. 41 5.2.2 Superimposed M o d e l......................................................................... 42 5.3 G LRT: Spatially L’ncorrelated A s s u m p tio n .............................................. 43 5.3.1 Replacement M o d e l............................................................................. 43 5.3.2 Superimposed M o d e l......................................................................... 44 5.4 Spatially Correlated Assumption ............................................................... 46 5.4.1 Local Statistics Com putation ........................................................ 46 5.4.2 G L R T .................................................................................................... 47 5.4.3 M ulti-Pixel G LRT ............................................................................. 47 6 Experiment and Result 48 6.1 Protocol of FDG-PET Dynamic S t u d y ...................................................... 49 6.2 C linical S tu d y .................................................................................................... 50 6.2.1 Lung Cancer FD G -P E T Dynamic S t u d y .................................... 50 6.2.1.1 Subspace E s tim a tio n .......................................................... 52 6.2.1.2 G L R T ..................................................................................... 55 6.2.1.3 Subspace Refining M e th o d ................................................ 58 6.2.1.4 M AP Reconstructed Im a g e ............................................. 60 6.2.2 Breast Cancer FD G -P E T Dynamic Study ................................. 60 6.3 Phantom Dynamic S tu d y ................................................................................ 64 6.3.1 Computer Simulated Phantom Data .......................................... 64 6.3.2 Receiver Operating Characteristic (ROC) S tu d y ....................... 72 7 Conclusion 74 Appendix A Oblique P ro je c tio n ..................................................................................................... 84 vi Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Appendix B Local Statistics in PET Im a g e s ............................................................................. S C B .l Local Statistics Calculation in PET FBP im ages...................................... S O B.1.1 Variance of Pixel in an R O I ............................................................. 88 B.1.2 Covariance of Pixels in an ROI ...................................................... 89 B.2 Covariance Com putation in M AP im a g e s.................................................. 90 B.2.1 M AP R e constru ctio n.......................................................................... 91 B.2.2 Approxim ation o f Covariance ......................................................... 92 Appendix C G LR T Can Not Be Improved by Orthogonal Subspace R e fin in g .....................95 Appendix D M ulti-P ixel G LRT D e riv a tio n ................................................................................ 98 v ii Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. List Of Tables 6.1 Dynamic data acquisition protocol ............................................................. 49 v iii Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. List Of Figures 2.1 (a) A PET Siemens/CTI ECAT 953 model, (b) Two 511-keY photons generated from the annihilation of a positron w ith a nearby electron. (c) Coincidence detection of the photons by a pair of detectors, (c) Scatter coincidence (left) and random coincidence (rig h t)...................... 10 3.1 FDG structure...................................................................................................... 20 3.2 (a) The transport and metabolic reaction pathways o f FDG compared w ith glucose in tissue. FDG is transported in tissue and phosphory- lated to FDG-O-P in the same manner as glucose, however FDG-G-P is not a substrate for the further reactions in the glycolytic pathway in the metabolic compartment. (b)The sim plified compartmental model for the transport of glucose and FDG in tissue.............................. 21 3.3 A 4-k compartmental m o d e l.............................................................. 25 3.4 The compartmental model for a heterogeneous tissue.................. 30 ix Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. G.l C linical dynamic studies used in this dissertation: (a) the lung can cer image (FBP). (b) the breast cancer image (FBP). (c) one TAC example from the lung cancer dynamic study, (d) one TAC example from the breast cancer dynamic study........................................................... 51 G.2 Test image reconstructed from FBP: (a) the last frame of the dynamic FD G -P E T FBP reconstructed images, where three ROIs were chosen: L I in the big lesion. L2 in the small lesion, and BG in the background. (b) the corresponding TACs for L l. L2. and BG. respectively...................52 G.3 Comparison of accuracy and separability of extracted subspaces of TACs in lesions and normal tissue by the SVD method: (a) TAC in L l represented by the normal tissue subspace, (b) TAC in BG represented by the lesion subspace, (c) TAC in L2 represented by tlie lesion and normal tissue subspaces, respectively. Comparison of accuracy and separability of extracted subspaces of TACs in lesions and normal tissue by the LS method: (d) TAC in L l represented by the normal tissue subspace, (e) TAC in BG represented by the lesion subspace, (f) TAC in L2 represented by the lesion and background subspaces, respectively....................................................................................... 54 x Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 6.4 (a) The last frame o f the dynamic FDG-PET FBP reconstructed im ages. where a rectangle indicates a test region containing a part of the big lesion and the whole small lesion, (b) 3-D mesh o f the test region indicated in (a), (c) G LRT with replacement model for TAC subspaces obtained by the SVD method, (d) G LR T w ith replacement model for TAC subspaces obtained by the LS method, (e) G LRT with superimposed model for TAC subspaces obtained by the LS method. . 56 6.5 (a) The test region indicated by a rectangle in the lung cancer study, where the prim ary lesion, the small lesion, and the heart area were all included, (b) 3-D mesh of the original test region, (c) The result of the G L R T ......................................................................................................... 57 6.6 (a) The test region indicated by a rectangle in the lung cancer study. (b) 3-D mesh of the original test region, (c) G LR T (w hite noise), (d) G LRT (non-white), (e) M ulti-pixel G LR T........................................................58 6.7 Subspace refining using the subspace distance measure. W ithout sub space refining: (a) lesion TAC and (b) normal tissue TA C both show the fidelity of subspaces estimated from observed TAC. (c) GLRT. W ith subspace refining: (d) lesion TAC and (e) normal tissue TAC shows the fidelity of subspace refined. It is dear that the refined sub spaces m aintain the fidelity requirement, (f) G LR T performance is improved by using the refined subspaces........................................................... 59 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 6.8 Subspace refining using the orthogonal projection, (a) test region inside the lung cancer image, (b) original 3D mesh. G LRT for the superimposed model: before, (c). and after, (d). the orthogonal sub space refining. G LRT for the replacement model: before, (e). and after, (f). the orthogonal subspace refining............................................... 6.9 The M A P reconstructed lung cancer FD G -PET image, (a) The test region indicated by a rectangular for the GLRT. (b) 3-D mesh of the original test region, (c) G LRT: white Gaussian noise, (d) GLRT: non-white Gaussian noise.............................................................................. 6.10 Breast cancer study: (a) a plane containing a prim ary cancer indi cated by an arrow, (b) a plane containing a metastatic lymph node indicated by an arrow, (c) the TACs from the prim ary cancer, the metastatic lym ph node and normal breast tissues in an ROI indicated w ith a rectangle in (b)................................................................................... 6.11 Breast cancer dynamic study: (a) two selected test regions: test region 1 containing the normal breast tissue and test region 2 containing a m etastatic lym ph node in the FBP image (same as (b)). (b)-(c) 3-D mesh of test region 1 and test region 2 indicated in (d). (d) TAC of the m etastatic lym ph node represented by the lesion subspace and normal breast tissue subspace identified by the SVD method, respectively. (e)-(f) the G LRT results on test region 1 and 2. respectively. . . . Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 6.12 A G LR T result o f enhancing an "unknown" lesion in the breast cancer FD G -P E T dynam ic study. The test region containing the unknown lesion in (a) F B P and (b) MAP. The 3-D mesh o f the test region: (c) FBP and (d) M AP. GLRT (white noise): (e) FBP and (f) MAP. . . 66 6.13 Left: lesion phantom image w ith 5 lesions, middle: normal tissue phantom w ith no lesions, right: TACs for lesion and normal tissues. 67 6.14 Left column: lesion phantom image w ith 5 lesions, middle column: normal tissue phantom with no lesions, right column: TACs for lesion and normal tissues. From top row to bottom : original phantom. FBP. OSEM. and M A P reconstructed images........................................................68 6.15 The G LRT for high count FBP phantom data: (a) a test region. (b)-(c) white noise assumption, (d)-(e) non-white noise assumption. . 69 6.16 The G LRT for low count FBP data: (a) a test region, (b)-(c) white noise asumption. (d)-(e) non-white noise assum ption............................... 70 6.17 The G LR T for an M AP phantom image (low count data): (a) a test re gion indicated by a rectangular, (b): G LR T (w hite noise assumption). (c): G LRT (non-white noise assumption), (d): m ulti-pixel GLRT. . . 71 6.18 ROC curves for the low count data sim ulation: (a) FBP: G LR T (white noise, non-white noise) vs m ulti-pixel G LRT. where .-L (m ulti-pixel G LR T) > .4; (G LR T with non-white noise) = A : (G LR T w ith white noise), (b) G L R T with white noise assumption: A : (M A P) > .-C (OSEM) > A : (F B P )......................................................................................... 73 xiii Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Abstract Positron emission tomography (PET)-scanning of [1 8 ]F fiuorodeoxyglucose (FDG)- Iabeled tissues is useful for the non-invasive diagnosis of cancer. However, as w ith all nuclear medicine imaging technologies, lesion detection w ith PET-FD G is restricted by the relatively lim ited spatial resolution and low signal-to-noise ratio, both ren dering diagnosis by visual inspection difficult and potentially inaccurate, especially when the lesion diameter is small. Hence, computer-aided lesion detection algo rithm s w ith feature identification have been developed to assist visual inspection in PET. The obstacle is that spatial features such as the shape or contrast of a lesion in PET images are often varied and difficult to identify. In this dissertation, we propose that spatial and temporal metabolic features available from dynamic PET-FDG images can be distinctly identified for normal and malignant tissue and then usefully applied in statistical hypothesis tests to improve the detection of small lesions. We hypot hesize that, when characterized in terms of a physiologically based compartmental model analysis, the time activity curves (TAC) derived from dy namic PET-FDG images may have variant shapes among tumors of different size and location w ithin given organs, but the sets of exponential functions or physio logical factors which compose the tumor TACs are the same. The expansion of the xiv Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. set of exponential functions or physiological factors w ill be called the tissue subspace in this study. Using a compartmental model for homogeneous and heterogeneous tissues, the dynam ic PET-FD G images can be formed into useful temporal-spatial data matrices, from which the subspaces can be estimated through parametric and non-parametric algorithm s. To utilize the identified subspaces for unknown lesion de tection. the matched subspace statistical detection criteria, based on the generalized maximum likelihood ratio principle, are adapted and extended. Results from clinical dynamic PE T-FD G studies demonstrated that the physiological features estimated from known tumors and normal tissues are distinct and valuable for improving lesion diagnosis when combined w ith the matched subspace detection criteria. Finally, a receiver operating characteristic analysis based 011 the Monte Carlo simulation of dynamic phantom images is completed to quantitatively compare the performance of different detection methods and image reconstruction algorithms. Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Chapter 1 Introduction 1.1 M otivation 1.1.1 P E T Diagnosis and C o m p u te r-A id ed Lesion D etection Positron omission tomography (P E T ), a rapidly developing technique for producing m vivo radiotracer d istrib u tio n images, provides im portant physiological informa tion. Non-invasive diagnostic procedures using PET imaging not only enable physi cians to assess the metabolic a ctivity of lesions in vivo, but also allow patients to avoid m ultiple costly procedures and to benefit from guided surgical intervention. The diagnostic advantage of P E T is im portant because alterations in the metabolism of cells are often evidenced in disease states before structural changes can be deter mined using other anatom ically based imaging technologies, such as computerized tomography (CT) and magnetic resonance imaging (M R I). Despite a resolution infe rior to that of CT and M R I. P E T has the following capabilities: (1) to characterize mass lesions as viable tum ors. (2) to determine whether m ild ly enlarged lymph 1 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. nodes represent a tum or or a non-malignant process, and (3) to detect tumors as small as 1 cm in diameter [70]. From the recei t research [1]. static PET imaging provides a higher sensitivity (959c) for detecting lym ph node metastasis than al ternative modalities, but further improvement is needed because the remaining 5% inaccuracy results in high risk and cost. Accuracy improvement to. for example. 989? would have enormous clinical and financial im plications for treating patients w ith suspected malignancies and for avoiding unnecessary health care expenditures key issues of concern in health care. The 59c uncertainty, mostly due to the lim ited spatial resolution, low lesion-to-background contrast, and partial volume ef fects. renders conventional diagnosis of small tumors by visual inspection difficult and potentially inaccurate. Improvements in tum or detectability may be realized by combining a sophisticated computer-aided detection algorithm with conventional visual inspection. The usual approach for computer-aided lesion detection is to identify features which distinguish lesions from normal tissues. However, conventional spatial fea tures. such as size, shape, and contrast, vary from lesion to lesion and the small size and low contrast of small lesions often restrict spatial feature identification [80]. Hence, for a given type of lesion, it is im portant to determine the features with the two im portant properties: readily identifiable and realistically useful for lesion discrim ination. Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 1.1.2 D ynam ic F D G -P E T S tu dy Dynamic P E T studies have been widely applied in medical images to visualize or classify structures in the images. Principal component analysis (PCA) [49] was introduced to enhance clinically interesting inform ation in a dynamic PET imag ing sequence in the first few principal component images, but it is a data-driving technique which can not separate signals from high noise level. Factor analysis of dynamic structures (FADS) was used to estimate the parameters of a compartmen- tal model [5]. Because FADS is used to solve an under-determined problem, there are an in finite number of sets of factors for possible solutions. A procedure called information-based factor analysis in dynamic structures (IBFADS) was introduced to incorporate a priori physiological inform ation to reduce the error in the estima tion of correct model. Inform ation was added by using suitable mathematical models to describe the underlying physiological processes. However, a single physiological factor was extracted representing the particular dynamic structure of interest. By tracking how much glucose is metabolized in different areas of the body. PET scanning with 2-[fluorine-18]-fluoro-2-deoxy-D-glucose (FDG ) not only enables physicians to map glucose utilization, but also provides the a b ility of im proving the diagnosis, staging, and treatment m onitoring of a variety of human tumors [12. 62]. M ultiple reports [12. 62. 63] have shown that FDG-PET scanning is highly accu rate in detecting early localized tum or recurrence w ith about 95% sensitivity. 98% specificity, and 96% accuracy, while CT provides only 65% accuracy. Thanks to the 3 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. PET insurance reimbursement policies approved by the Health-Care Finance Ad m inistration (HCFA) in 199S. the FD G -PE T study is becoming im portant because the availability o f FD G -PET scanners progressed rapidly nationwide [63]. In this dissertation, the temporal-spatial features of time activity curves (TAC) obtained from dynamic FD G -PET studies w ill be shown to have the following prop erties [28. 82]: (a) linearly representable by a set of exponential functions or physi ological factors, (b) physiologically distinguishable as lesion and normal tissue sub spaces. (c) easily identified by using the temporal-spatial data matrices, (d) readily incorporable to a matched subspace detector for lesion detection. Our approach uses the fact that malignancies can be distinguished from normal tissues on the basis of biologically determined radiotracer accumulation or loss rate. Because malignant tumors are m etabolically active and FDG-avid on PET imaging, they metabolize glucose at a much higher rate than do most normal tissues [33]. Hence, the TAC of tum or often increases w ith time. From a physiological compartmental model analysis, a dynamic FD G -PET study can be interpreted as describing the passage of administered tracer through a finite number of independent homogeneous compartments. Each compartment is associ ated w ith a particular FDG kinetic structure and can be described by the fundamen tal exponential functions or physiological factors which w ill be called 6a.sj.s in this dissertation. Because of the lim ited resolution of a PET system, a heterogeneous com partm ental model, a m ixture of the finite number o f independent homogeneous compartments, would be more realistic to describe the kinetic data. Then a TAC 4 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. is composed of a combination of the basis which constitute the tissue subspace. Based on the physiologic compartmental model analysis, each tissue's kinetic be havior is governed by its own set of kinetic parameters, which suggests that TAC shapes may vary among different tumors w ithin a given patient, but w ill occupy the same subspace provided they have the same kinetic model structure. Therefore, it may be possible to use the subspaces identified from visible. large lesions to confirm suspected, but not unequivocally identifiable, small lesions. 1.1.3 D a ta Processing and Sim ulation By forming the dynamic FDG -PET data into a tem poral-spatial m atrix, the lesion and normal tissue subspaces can be sim ply estimated by using the signal-subspace- fittin g methods. The least square's error (LSE) algorithm is used for parametric estimation (exponential functions): while the principle component decomposition is used for non-parametric estimation (physiological factors). In order to utilize the identified subspaces. several statistical tests are derived for unknown lesion detection based on a hypothesized superimposed or replacement data model. The statistical tests are called the matched subspace detectors [55. 78] which can be derived from the generalized likelihood ratio test (GLRT) principle. The matched subspace detector is actually the extension of the rank-1 matched filter to the m ulti-rank filter. In this dissertation, both the G LR T (single-pixel) and the m ulti-pixel G LR T w ill be derived for the lesion detection. The GLRT was derived in [55] by assuming that pixels are uncorrelated spatially. Considering that image pixel values are not statistically 5 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. uncorrelated because each count contributes to all pixel values [?]. a pre-whitening procedure is required before using the G LRT. The spatial covariance m a trix required for the pre-whitening can be computed using the methods described in [7. 30. 31]. In the m ulti-pixel GLRT approach, it is assumed that the spatial inter-pixel correlation in each frame has the same structure but different energy level. In this dissertation, we w ill concentrate on the filtered backprojection (FBP) and the iterative reconstructed images. The results from clinical dynamic FD G -PET studies show that the physiological features estimated from known tumors and nor mal tissues are distinct and invariant, and are valuable for im proving lesion diagnosis by incorporating a statistical test. Finally, to quantitatively compare the sensitivity and sufficiency of different detection methods and image reconstruction algorithms, a Monte Carlo computer sim ulation of dynamic phantom images is performed and used in a receiver operating characteristic (ROC) analysis. G Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Chapter 2 Basic Principles of PET 2.1 Overview PET enables physicians to pictorially determine in vivo m etabolic activity by tracing a radioactive isotope. The isotope decays during the biological processes and can be detected by a cylindrical array of detectors in a PET scanner. PET imaging combines the early biochemical assessment of pathology achieved by nuclear medicine w ith the precise localization achieved by computerized image reconstruction algorithms. In this chapter, the physical basis for PET [8 . 32. 39. Gl] and image reconstruction methods w ill bo reviewed. 2.2 Physics of PET 2.2.1 Positron Emission and Photon D etectio n A PET Siem ens/CTI ECAT 953 model is shown in Figure 2.1 (a). In a PET study, a positron-em itting radioisotope is injected into or inhaled by a subject. The isotope, a positron em itter, then circulates through the bloodstream to reach, for example. Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. brain tissue or cardiac muscle. Positrons are em itted from the nucleus of some ra dioisotopes that are unstable because they have an excessive number of protons and a positive charge. Positron emission stabilizes the nucleus by removing a positive charge through the conversion of a proton into a neutron or a positron. A fter the decay of the mother nucleus, an emitted positron, after traveling a short distance, combines (annihilates) with an electron of a nearby atom (w ithin ~ 2 mm of the tracer). The distance that a positron travels is dependent on its energy. The re sulting annihilated energy is carried by two 511-keY(= mQ c2) photons that travel in approximately opposite directions (180” ± 0.53) to each other. The annihilation process of two photons is shown in Figure 2 .1 (b). The slight deviation from 180” is due to the residual kinetic energy of the positron. The fact that the annihilation photons are em itted simultaneously and in opposite directions is the fundamental basis for detecting and localizing the positron em itter. W ith some probability, this pair of photons, after penetrating the subject, can be recorded eoincidently by external radiation detectors. Scintillation detectors com posed of crystals, such as bismuth germinate (BG O ). and photom ultiplier tubes are posed on opposite sides of the radiation source. If both members of the detector pair receive photon signals w ithin a very narrow tim e interval (in the nanosecond range), a coincidence event w ill be detected, see Figure 2.1 (c). Coincidence de tection provides a unique detection scheme for form ing tomographic images with S Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. PET. The near collinearity of the two annihilated photons, enabling the identifica tion of the annihilation location through the coincidence detection, determines the line-of-response (LO R ) along which the annihilation occurs. To more accurately localize the annihilation event along the travel paths, the PET scanner consists of detector elements mounted on one or more rings, posi tioned so that they surround the subject. In a complete PET system w ith more than 10.000 detectors, over 20M coincidence lines are collected simultaneously, which form a large number of intersecting LORs to provide inform ation about the quan tity and spatial location of positron emitters w ithin the body. However, there are some factors, for example, the scatter and random coincidences shown in Figure 2.1 (d). which w ill affect the accuracy of PET images. Hence, some corrections (e.g.. attenuation correction, accidental coincidence correction, detector efficiency correc tion. dead tim e correction) are required before the collected data in a PET scanner can be processed by the tomographic reconstruction algorithm s. The reconstruction algorithm processes the coincidence events measured at all angular and linear posi tions to reconstruct an image that depicts the localization and concentration of the positron-em itting radioisotope w ithin a plane of the scanned organ. 2.2.2 R adio tracers Used in P E T PET, compared w ith the conventional nuclear imaging systems, offers a more useful choice of biologically significant chemical elements for labeling [39]. A ll biological activities are the result of biochemical reactions: hence, for every pathology there 9 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. (C) (cl) Figure 2.1: (a) A PET Siem ens/CTI ECAT 953 model, (b) Two 511-keV photons generated from the annihilation of a positron w ith a nearby electron, (c) Coincidence detection of the photons by a pair of detectors, (c) Scatter coincidence (left) and random coincidence (right). 10 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. is an underlying biochemical defect. The goal of medical imaging techniques is to visually identify abnormal biochemical a ctivity resulting from a given pathology and to observe the abnorm ality directly. Depending on the metabolic activity to be traced, one of four radioactive isotopes is usually chosen for PET implemen tation: ( 1 ) oxygen-15 (half-life. ~2 min: labeled metabolite: oxygen and water): for metabolism quantification, blood volume estimation: (2) nitrogen-13 (half-life. ~ 1 0 min: labeled metabolite: ammonia, various amino acids and nitrous oxide): for organ perfusion, tum or drug levels test: (3) carbon-11 (half-life. ~20 min: labeled metabolite: carbon monoxide, carbon dioxide, palm itate. glucose): to test tumor metabolism, brain permeability, blood flow, and tissue lipid content: (4) fluorine- 18 (half-life. ~ 110 min: labeled metabolite: 2-deoxy-d-glucose. 3-deo.xy-d-gluco.se. haloperidol): for brain metabolism detection, brain pharmacology, and neurorecep tor. The long-lived positron-em itting isotope of fluorine. 18F. is substituted for oxy gen in an analogue of normal metabolism or drugs, and i8 F. because of its longer half-life, is usually used for testing metabolism of the brain. The short-lived isotopes of oxygen, nitrogen, and carbon are particularly useful for tracing the metabolism of their normal counterparts but can also be incorporated into compounds of pharma cological interest. Given that their existence is so fleeting. 1:,0 . UC. and llJX must be generated close to the point of detection. To produce and handle such short-lived isotopes, it is sometimes necessary that PET systems be equipped w ith cyclotron designed for creating a variety of positron isotopes. 11 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 2.3 Statistical Model for PE T Imaging A fter a compound containing a radioactive isotope is introduced into the subject and forms an unknown em itter density A (i\ y. z) > 0. the isotope decay is gov erned by an exponential distribution [60]. Therefore, by the fact that the addition of independent exponential random variables results in a Poisson distribution, emis sions occurring in the region of interest are according to a Poisson process w ith rate A(x. y. z). The two 511-keV photons are detected coincident ly by a pair of detector elements that define a volume, usually called a detector tube cl. from which the photons orig inated. Thus, the only information acquired when a pair of detectors count a coin cidence is that somewhere along the length of the tube cl the annihilation occurred. Hence, the set o f data collected in a PET scan is the tube count [ n '( l) .........ri'(D)]. where u*(r/) is the total number of coincidences counted in the c/th detector tube and D is the total number of tubes. The problem in a PE T imaging system is to estimate the unknown emission density A(.r. ;/. c) from the observed data n“. Is - ing a counting process model, each measured data, ri'(d). can be modeled as an independent Poisson random variable w ith mean A(r/). p(n(d) = k ) = e -A (< ' > ^ . A - = 0 . 1 . 2 . . . . (2 .1 ) Also, for the purpose of display and machine com putation, it is convenient to dis cretize the density A into boxes 6 = 1 B at the outset. Thus, for each box 6 12 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. there is an associated unknown count n(b) w ith mean A(b) = E{n(b)}. b = 1 B. It is further assumed that once an emission occurs in some box b. the conditional probability that it is detected in tube d is independent of other emissions and is given by p{b . d) = ^(detected in r/|occurring inb). d — 1 D. 6 = 1 B (2.2) where p(6 . d)'s are known a priori nonnegative constants, which can be computed from the detector ring geometry. Given that n(b. d) is the number o f emissions occurring in box 6 and detected in tube d. d = 1 D. b = I B. then n(b. d)'s are Poisson random variables. independent of each other, w ith means A( 6 . d) = A(6 )p( 6 . d). d = 1 D. 6 = 1 B. (2.3) Based on this description, the measured data [n * (l) .........^ * (^ ) ] can be re-expressed as n'(d) = n {.. d) = Y , n(b . d) 6=1 = the total number of emissions detected in tube d. d = 1 D. (2.4) It is further assumed that n (6 . .) = £ n (b - d) d= i = the total number of emissions that occurred in box 6 and detected. (2.5) 13 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. then E {n (b .)} = \(b)p(b..). (2 .6 ) Given that p(b. .) is known a priori, the PET imaging reconstruction problem now is to estim ate the mean A(b). b = 1 B. or the true unobserved count n{b) for each box. from the observed data n’ [d). d = 1 D . 2.4 PET Image Reconstruction 2.4.1 F ilte re d B ackprojection M e th o d The conventional FBP method for reconstructing the PET images is simply the im plem entation of the inverse Radon transform [59]. The FBP method processes the sinogram using a linear filter, then each pixel value is derived by summing the filtered data from each angle. The FBP method is a standard method for PET image reconstruction that it is fast, easy to implement, and theoretically tractable. By the linearity of FBP. the local statistics of PET FBP images can be readily estimated. Due to the poor count statistics encountered in PET. conventional FBP tomographic image reconstruction in use today exhibits severe stripe-like artifacts. 2.4.2 Ite ra tiv e Im age Reconstruction M e th o d Maximum Likelihood (M L) Reconstruction Method: When applied to PE T image reconstruction, the M L method was shown to reduce statistical noise artifacts conventional FBP algorithms w ithout excessive smoothing 14 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. [60. 6 8 ]. This advantage is im portant because statistical noise is a m ajor lim iting factor in emission tomography. In the PET image reconstruction procedure, for each unknown density A. the observed data n* has probability or likelihood p (n '|A ). The EM algorithm , which computes an M L solution, begins w ith an in itia l estimate o f A0 and the following simple iterative procedure is used for obtaining a new estimate A,iett' from an old estimate XolJ. to obtain A*. A:=1.2....... A— (6) = , 6 = 1 ..........o. ,2.7) ,i=i £ XM {b')p{b' . d) b' - 1 This iteration w ill continue until some criterion (e.g.. maximum likelihood) is ful filled. The likelihood strictly increases at each step (unless it reaches its maximum). [){n*\\k) > p{n’ \Xk~l ). k > 1 [13]. The total number £A (fr) o f estimated count was also shown equal to the total number £ n'(d) of observed counts at each step Xk. k > 1 [13]. Finally. Xk converges as k — > oc to an estimate A70 which has maxi mum likelihood. Bayesian Reconstruction Method: The Bayesian reconstruction algorithm s [23. 24. 43] form a powerful extension of the M L reconstruction method. The Bayesian algorithms are based on a statistical model for the PET sinogram data and a probabilistic prior image model adapted from Markov random field theory. To improve the image reconstruction quality. 15 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. the Bayesian algorithm is designed to incorporate into its reconstruction formula tion prior inform ation such as smoothness constraints or partial specified topological inform ation. The Bayesian reconstruction algorithm identifies an image by maximiz ing the posterior density function, hence, it is designated the m aximum a posterior (M AP) method. Ordered Subsets Expectation Maximization (OSEM): Although the resolution is excellent compared to that of FBP images, the applica tion of M L or M A P algorithms can be computer intensive and convergence slow. Recently, a new algorithm , the ordered subsets, was proposed to accelerate iterative reconstruction [26]. The OSEM method groups the projection data into an ordered sequence of subsets (or blocks). An iteration in the OSEM algorithm is defined as a single pass through all the subsets, in each subset using the current estimate to initialize application of EM w ith that data subset. 16 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Chapter 3 PET-FDG Dynamic Data Modeling 3.1 Overview Tracer kinetic techniques [67] are generally used in physiology and biochemistry to trace dynamic processes, such as blood flow, substrate transport, and biochemical reactions. The tracer kinetic method involves a radiolabeled biologically active com pound (tracer) and a mathematical model that describes the tracer kinetics during a biological process. By the application of the tracer kinetic method. PET can ac curately estimate in vivo the tracer concentration, which provides functional and physiological variables, for example, the local cerebral m etabolic rates of glucose (LC M R G lc) [58]. Because FDG is useful in clinical PET for functional protocols (e.g.. lesion detection, cardiology, neurology, oncology), we w ill focus on FDG-PET dynamic studies in this dissertation. After being injected intravenously as a bolus. FDG concentration is measured to determine the transport rate constants of FDG in physiological compartmental models [35. 37. 58. 6 6 . 67]. The kinetic behavior of the tracer concentration in a tissue is called a time a ctivity curve (TAC). By 17 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. using a 4-k (i.e.. 4-parameter) compartmental model for a homogeneous tissue, the TAC function can be modeled as a linear combination of the fundamental exponen tial functions or physiological factors [15. 37. 44. 45. 47] which are called a basis in signal processing. For a heterogeneous tissue containing a finite number of indepen dent homogeneous compartments, the TAC can be represented as a weighted sum of the independent exponential functions or physiological factors which compose a subspace. 3.2 Tracer Kinetic Technique The tracer kinetic techniques are among the most powerful methods for measuring the rates of processes for increasing knowledge of the biochemical, transport, and metabolism of body functions. A small amount of detectable substrate (a tracer amount) can be introduced into fluid, and the speed of its passage through the system can be measured. The amount of label appearing in other chemical species w ill rise w ith some tim e delay and w ill later decrease gradually. From the amount of time delay, the rates of increase and decrease, and the concentration levels, the transport rates in a dynamic system can be estimated. However, the tracer kinetic techniques are more complicated in practical situations, e.g.. the selection of tracers, the non-steady state condition, the mixture of inhomogeneous tissues, and the design of model for the dynam ic process. Hence, it is necessary to design a tracer kinetic model to not only follow the dynamic process of interest, but also be mathematically tractable. IS Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 3.2.1 Tracer Selection In tracer kinetic techniques, an appropriate tracer is required to follow the dynamic process o f interest. Usually, a tracer needs to have the following properties: (a) structurally related to the natural substance involved in the dynamic process (i.e.. metabolic processes) or have sim ilar transport properties (i.e.. for How system), (b) measurable and distinguishable from the natural substance intended to be traced. The well-known tracer used in PET is 2 -[l8 F]fluoro-2 -deoxy-D-glucose (FDG ). which has been used to isolate the transport and phosphorylation steps from the com pli cated pathway of glucose metabolism, thus enabling the form ulation of simple tracer kinetic models for quantification of glucose utilizatio n. DG is an analog of glucose in which the hydroxyl group on the number 2 carbon has been replaced w ith hydrogen. For PET. this hydrogen is substituted by the positron em itter 18F to form FDG. see Figure 3.1 [67]. DG and FDG have been used in biochemistry to isolate the phospho rylation reaction from the rest of glycolytic metabolism. FDG is transported into tissue and phosphorylated to FDG-6 -phosphate (F D G -6 -P) in the same manner as DG or glucose. However, because of the substitution in the second carbon position. FDG -6 -P is not a substrate for the next reaction step in the metabolic pathway, see Figure 3.2 [50]. FDG-G-P does not leave the cell except through a slow hydrolysis back to free FDG. which can then be transported to plasma or be rephosphorylated. Basically, FDG behaves sim ilarly to glucose in its transport from plasma to tissue 19 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. and in its phosphorylation, the rate of which, under a steady-state condition for glu cose. is equal to the utilization rate of exogenous glucose. Therefore, the utilization rate of glucose is predictable from the kinetics o f FDG. FDG uptake for a variety of human cancers has been shown to have high tumor- background uptake ratios at 1-2 hour post i.v. injection [71]. The mechanisms for this increased FD G -6 -phosphate accumulation in many cancer cells has been shown due to [72]: • increased expression of glucose transporter molecules at the tum or cell surface • increased a ctivity of hexokinase • reduced levels o f glucose-G-phosphatase v.s. most normal tissues. The property that cancer cells avidly accumulate FD G -6 -phosphate is im portant for lesion detection and is fundamental to the method introduced in this dissertation. HO H O HO Figure 3.1: FDG structure. 20 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Metabolic Compartment G lvcogen Vascular Compartment Phosphorylase 'a' Free Space Glucose Glucose G-6-P04 F -6 -P 0 4 FDG FDG Cell Membrane Glucose FDG FDG FDG-6-P FDG (b) Figure 3.2: (a) The transport and metabolic reaction pathways of FDG compared with glucose in tissue. FDG is transported in tissue and phosphorylated to FDG-G-P in the same manner as glucose, however FDG -6 -P is not a substrate for the further reactions in the glycolytic pathway in the metabolic compartment. (b)The sim plified coinpartmental model for the transport of glucose and FDG in tissue. 21 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 3.2.2 M o d e l Selection For tracer kinetics, there exist various types of mathematical models of widely differ ent mathematical characteristics, such as deterministic versus stochastic, distributed versus non-distributed. compartmental versus non-coinpartmental. and linear versus nonlinear. In biochemical applications (e.g.. PET), linear com partm ental models are most frequently used, because of their advantageous m athem atical properties that enable straightforward solution or analysis of the model characteristics. A compartmental model, see Figure 3.3. is usually represented by a number o f compartments connected by a number of arrows indicating transport between compartments which are defined as spaces where the tracer is distributed uniformly. The amount of tracer leaving a compartment is assumed to be proportional to the tota l amount in the compartment. The arrows indicate the possible pathways the tracer can follow. The symbol k. the rate constant w ith the unit of inverse time, denotes the fraction of the total tracer that would leave the compartment per unit time. Based on the linear compartmental model assumption, the tracer kinetics of a compartmental model can be described in terms of a set of linear, first-order, constant-coefficient, ordinary differential equations. If certain conditions are fulfilled, e.g.. in itia l conditions are zero, the solutions to the differential equations w ill simply be the convolution of the input function and the system response function. The response function of a compartmental model w ill consist of a sum of exponential functions. Usually, the number of exponential components is equal to the number of 22 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. compartments in the model. The fact that the response functions o f many biological or physiological systems can be decomposed into a sum of exponential components is another indication that the com partm ental models are adequate for tracer kinetics studies. 3.2.3 Assum ptions in the Tracer K in etic M o d e lin g In tracer kinetic modeling, the following assumptions are usually made in order for the model to have tractable m athematical properties: • The tracer introduced is assumed to be in a trace amount, so th a t the measured results would reflect the effect o f the tracer introduction and not represent the original process. • the dynamic process for tracer kinetics is assumed to be in a steady state. That is. the rate of transport or reaction of the system is unchanged w ith time, and the amount of substance in any part of the system is constant during the eval uation time. For practical systems, because biological systems continuously change to adapt to the environment, there is no absolute steady state. How ever. the steady state condition is considered to be sufficient if the amount of change w ithin the tim e of evaluation is minor compared to the magnitude of the process itself. 23 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. • The contribution of different tissues to a tracer kinetic measurement is another unavoidable practical problem. W ith the advantage of PET. the concentra tion of positron-em itting in small local regions of the body can be not only quantitatively measured, but also anatom ically differentiated to a significant degree. W ith the use of radiotracers and scintillation detectors (individual detectors and scintillation cameras), the concentration of positron-einission radioactivity in small local regions of the body can be quantitatively measured with PET scanners. This capability has a large impact on the tracer kineric model, such that not only can the absolute tracer concentration in tissue be obtained, but also anatomically heterogeneous tissues can be delineated to a significant degree. Hence, a small tissue element can be reasonably regarded as a homogeneous uniform pool. • As far as the transport of tracer is concerned, the only connection between a small tissue element and other parts of the organ or body is the tracer delivery and clearance through the blood flow in the tissue element. Therefore, the tracer concentration in the supplying blood that can be measured from peripheral arterial blood samples or an ROI data drawn inside a PET heart images can be considered as the input function to the tracer kinetic system in the small tissue element. Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 3.3 A 4-k Compartmental M odel for the FDG Tracer K inetics 3.3.1 Homogeneous Assum ption The subject of a FD G -P E T dynamic study, under the homogeneous tissue assump tion. could be considered as a plasma compartment and a tissue compartment. Fig ure 3.3 shows a 4k com partm ental model for FDG. where the compartments repre sent FDG in plasma. FDG in tissue, and phosphorylated FDG in tissue (FDG-G-P). The rate constants in Figure 3.3 are defined as follows: A’* = the transport of FDG in plasma to tissue k\ = the transport of FDG in tissue to plasma A .3 = the phosphorylation of FDG to FD G -6 -P k\ = the dephosphorylation of FD G -6 -P to FDG. c*(t) C *(/) I I I ' ' Plasma Cell Tissue Cell Figure 3.3: A 4-k compartmental model 25 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. The following differential equations describe the rate of change in the concentra tions of FDG and FDG-6 -P in the homogeneous tissue: where * indicates decay-corrected FDG tracer quantities. C *(f) represents the con centration of FDG in the tissue at time t. C'ri{t) represents the concentration of FDG-G-P (and all metabolites derived from FD G -6 -P) in the same region, and C ’ (t) represents the FDG concentration in the arterial plasma. Solutions of (3.1) [58] given that the in itial conditions C *(0 ) = C *,(0 ) = 0 show that the total radioactivity for a homogeneous tissue. C '{t). is the sum of the free [lsF]FDG concentration plus the concentration of metabolites, i.e.. (3.1) c-(t) = c;(o + c,;io (3.2) where 26 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. .1 /, = — — — (a-, - A ri - k'.) (3.3) Q > — 1 and © denotes the convolution operator. The total radioactivity. Ct *tal(f). measured from a homogeneous tissue including the radioactivity in the cerebral blood can be represented as C L J t ) = vbc ;{t) + c -{ t) = Vbc ;(t) + [U te - " '1 + .\/oe'Q -'] © C ‘ p(t) (3.4) where FDG concentration contributed from blood pool is assumed to be \ lCp(t) and \ ’ b denotes the vascular space in tissue [50]. In order to characterize the complex behaviors of the realistic plasma TAC. such as the period of zero activity at the beginning due to delay from tracer delivery, the rapidly rising period, and the exponential-like decay period. Feng and Wang [16] used a compartmentalized model to represent tracer behavior in the blood circulatory system.where The FDG plasma TAC was represented by a 4th-order exponential curve w ith a pure decay and a pair of repeated eigenvalues: c;(t) = {[B ,(f - r) - B 2 - £ : 1]e - v" '- rl + £ 2e -A *“ - rl + 1 (t - r) (3.5) Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where A[. A2. and A3 are the eigenvalues of the blood circulatory system ( 1/m in ): Bi. B 2. and B-j are the coefficients (//C ,/m l): r is the delay constant (m in): and i ( i - 0 = ' °- i f , < r ‘ 1. if t > T. (3.6) A substitution of Eq. (3.5). assuming r to be 0. into Eq. (3.4) yields C'MAl(t) as C L J t ) = 4- E2e -X ‘f 4- E3e~x'-‘ 4- E ,e '^ ' + E3e-“ ‘' 4- E e e ^ 1 (3.1 where E\ — \ bB[ + A A /,B , ABB, ( A i — ci i ) — (A i — a 2) E -> = - \ l{B-, + B3) + \l\B-> - M i B i .\h(B-2 + B 3)' 4- r ~ .\E B i * M 2(B 2 + B.i)' (A1 - o 1)- A [ o [ [(A, - a , ) 2 A i — n 2 £ • :( — I t>B2 4 “ A/>B> — (X> — o i ) — (A2 — a->) r a i - r> i UxB3 , M >B, E.i — VftB;) 4----- — ---------- r 4- E, Ee — (A;j — Oi) — (A;J — a 2) M \B \ Mi{B-> + B.i) M i B, A/,B, 4 -------— --------------- r h A (Ai — 0 [ ) J — (A[ — a i ) A) — A .-j — a , M ;B i A/>(B> + Bi) M,B-> M ,B 3 - + (Ai — Ckj)“ —(A[ — Oj) A2 — C l2 A .) — Oj (3.8) For t — t[. ■ ■ ■. t \ . X is the number of total frames, the radioactivity C ’otil(t) for the homogeneous tissue can be represented as a column vector = P ap 4- T ,a i (3.9) 28 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where T t = /^O * I ^ C C C t->c~Xl'- c~Xlt- e- '-'- e-A3<- ’ p-,*If. Y p_,'2'.Y .V X 2 ap = [E .-E ^ .E , a , = [ E , . E ; ] r . (3.10) In this dissertation, the two exponential column vectors in the m atrix T i w ill be called the exponential basis for the homogeneous tissue, whereas the combination of them, resulting in the one-dimensional vector T ia ^ w ill be called the physiological factor. Both the exponential basis and the physiological factor w ill be estimated by the methods in Chapter 4. 3.3.2 Heterogeneous Assum ption Due to the lim ited spatial resolution of the PET scanner, measurement of radioac tiv ity in a homogeneous tissue is rarely, if ever achieved. For this reason, we use m ixture analysis or physiological factor analysis (PFA) to incorporate the hetero geneity effects for a heterogeneous region of interest, where the mixed tissue TAC can be represented as the mass-weighted average of physiological function in each of the finite homogeneous tissues [35. 57. 47]. For a mixed (heterogeneous) tissue 29 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. pixel containing a finite number of ./ overlapping independent homogeneous com partments. each of which is described by a physiological factor, the TAC w ill be represented as a weighted sum of the J independent physiological functions. Figure 3.4 shows a schematic diagram of a heterogeneous tissue. j c*t„ial{t) = Vv»C*p (/)+ £ “'/•£ * (0 / = 1 Figure 3.4: The compartmental model for a heterogeneous tissue. 30 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Given that a mixed tissue pixel contains J homogeneous subregions w ith the radioactivity in each subregion j described by. from Eq. (3.2). = ky [(*3j + - «1 - (Ci2j - a- - A - ) e - v ] 9 C'p(t) (3.11) then C *(f). 1^° weighted average radioactivity in a heterogeneous tissue, can be w ritten as c * ( 0 = X > j C ; ( 0 j= i j = Y. «•; [-'V " 1 '' + .1 © c;(t) (3.12) 7 = 1 where Wj is the weighting coefficient for the _/th homogeneous subregion: M \j = J '-a ,, (k jj + A 'ij ~ ° u ) aild (n '2j ~ k 'ij ~ k'\j)- tho variables I\\,. k2j. k:ij. k.ij. ciij, and a>j for the j t l i subregion in a heterogeneous tissue have the same physiological meanings for those defined in the homogeneous tissue case. Then the total radioactivity for the heterogeneous tissue can be derived by substituting the plasma function defined in Eq. (3.5): CLJt) = vb c;(t) + Y [.'A,e-^‘ + .v „e -* ‘] © c;(t) 7 = 1 = E y tc - ^ + E2c~ ^‘ + E3c~x-‘ + j j + Y VjE ojC-"11' + Y (3.13) 7=1 7=1 31 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where Ey = VbB l + Y . lL 'j J = 1 — (A i — a i j ) —(A[ — E-i = - \ \ { B 2 + Bi ) + Y . wj j =i (A[ — ayj)2 A[ — Q'i -M>jBy \[2 j{B -2 + B3) . (A [ — a->j) - A [ — ci>j E, = 1 bB) + ^ Wj j =i j E, = VbB, + Y . lL 'j j =i MyjD, M-2j B-2 - ( A j - a i j ) — (A-j — a2_,)_ M ijB , + M ijB , . - ( A:* - a i j ) — ( A3 — a 2 j ) = E<\} = A MxjBy UyjiB. + B,) M XjB2 MyjB, + + + (A i — a ^ )-’ — (A i — a i j ) \> — ctij A3 — a X j M 2j By M-ij{B-2 + By) M 2jB2 \E jB i ( A [ — c y 2 j)‘ — (Ax — Q j j ) A ’? — o .’ 2 j A;j — 0 -2 j + (3.14) Hence, the weighted average radioactivity for the heterogeneous tissue can be rep resented as a colum n vector x = Pa„ + T a (3.15) where T = [T l . T-j. • • •. T jJ.v■ . >j- T ,. the exponential basis, is defined in Eq. (3.10). a p = [ E i , E>. E-y. Ejy\T . and a = [E-^y. £>,2. • • • . E-,j. Ef)y. E r # . ■ ■ ■. E & A 1 ■ In the physi ological factor representation, x can be w ritten as x = Pa„ + T 'a ' (3.16) 32 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where T' = [T[. T'2. • • •. Tj ].% -* ./. T'. the j th physiological factor, is the combination of the column vectors in T 2. and a' are the associated coefficients. A more compact m atrix format for the radioactivity of x would be x = Ea (3.17) where E = [P. T]. if represented in the exponential basis format. E = [P. T']. if represented in the physiological factor format, and a is the corresponding coefficient vector. 33 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Chapter 4 Dynamic Data Formulation and Subspace Identification 4.1 Overview In computer-aided lesion detection w ith FD G -PET dynamic images, the TACs de fined in Chapter 3 are useful for separating lesions from norm al tissues [29]. However, how to identify the exponential subspace functions and incorporate the identified subspaces into a detector remains a problem. In this chapter, based on the dynamic temporal-spatial data structure, signal-subspace-fitting methods, such as the least squares (LS) m ethod and the singular value decomposition (SVD) method, are used to estimate param etric (exponential function) and nonparametric (physiological fac tor) subspaces, respectively. 4.2 Dynam ic Data Formulation One of the advantages in the FD G -PET dynamic study is that the dynamic images can be formed into a temporal-spatial data m atrix such that some useful signal array 34 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. processing techniques can be applied. Let vq. be the kinetic data measured by the PET scanner for the p-th pixel in a heterogeneous ROI. then yk = + n* = Eaf c + nf c (4.1) where np denotes the noise and is assumed to be Gaussian. Given that P denotes the number of pixels in the ROI. the temporal-spatial data m a trix for the ROI from a PET scanner can be expressed as Y = [xt. x-j. ■ • •. X /-] + [ni. n -_ > . • • •. n/>] = EA + N (4.2) where A = [ai.aj. • ■ • .a/>] N = [ni.n-j.---.npj. (4.3) In array signal processing, the column space of the m a trix E in (4.2) is called the signal subspace [69]. In this dissertation, the signal subspaces are called the lesion and normal tissue subspaces which are denoted as (H) and (S). respectively, where (H) is the subspace spanned by the columns of the .Y xp m a trix H. (S) is the subspace spanned by the columns of the X x t m atrix S. and p + t < X. Based on the compartmental model analysis, the subspaces (H) and (S) have the input 35 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. plasma subspaco in common, otherwise, they are linearly independent in which case the columns o f the concatenated m atrix [HS] span the (p + ^-dimensional subspace (H S). 4.3 Subspace Identification from Know Tissue T ype Two categories of subspace identification methods are used in this dissertation, namely, non-parametric and parametric methods, for solving the signal-subspace- fittin g problem. The non-parametric techniques used here simply decompose the observed data to extract the subspaces using a singular value decomposition method. The param etric method is used to estimate the parameters of interest in spectrum like functions, e.g.. least squares estimation for the exponential parameters in Eq. (3.7). 4.3.1 N o n -P ara m e tric: Singular V alu e Decom position The non-parametric methods find a m ultidim ensional signal subspace by a simple eigen-decomposition of the autocorrelation m atrix or a singular value decomposition of the data m atrix. The SYD method is used in this dissertation to estimate the physiological factors defined in Eq. (3.17). Given that the coefficient m atrix A in Eq. (4.2) is uncorrelated w ith the noise N which is assumed to be white, the autocorrelation m atrix of the data Y can be decomposed as follows: R y = Y Y 'r 3G Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. = H A A r H + a ' l = $ [ A - F < r ; I ] $ r = # sA s^ r + $ nA nC (4.4) whore we assume that A A T. is a full rank m atrix, and therefore H A A r H r can be eigendecomposed as «£sA«£if such that (H ) = ( ^ a). The r eigenvalues of the decomposition combine w ith the noise covariance to form the r x r diagonal m atrix A , = A 4- cr^I. w ith the eigenvalues in A , arranged in decreasing order. The (.V - r) x (.V — r) diagonal m atrix A „ contains the .V — r repeated eigenvalues a~. Thus the eigen-decomposition of the data autocorrelation m atrix in Eq. (4.4) results in a signal subspace < !> , and a noise-only subspace 4.3.2 P aram etric: Least Squares E stim ation Finding the subspace using the param etric algorithm is equivalent to the estimation of the exponential basis in Eq. (3.7). which is sim ilar to the well-known problem of finding the direction of arrival (DO A) in signal subspace processing [09]. In this dissertation, the least squares (LS) method is applied for the param etric method and the best fit of parameters to the observed data are chosen to constitute the subspaces. Subspace fittin g using the LS method is defined as [69]: H . A = arg min ||Y - H A | | f (4.5) H. A 37 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where || • ||/r denotes the Frobenius norm. Because the subspace fittin g problem is separable in H and A [69]. by substituting the pseudoinverse solution. A = (H r H ) _ 1H r Y . back into Eq. (4.5). the following equivalent problem is derived: H = a rg m in ||Y - H ( H r H ) - 1H r Y ||^ H = arg nun 11P 7? Y ! | H = arg max ||P //Y jfp H = arg max tr { P „ Y Y r } (4.C) H where P // = H (H r H ) - 1H r is the orthogonal projection m a trix that projects onto the column space of H and "tr" means the trace operator. In feature identification using parametric algorithm for a FDG-PET dynamic study, the only unknowns are the exponential parameters, the LS method becomes to search the parameters over a reasonable range for the tracer kinetic parameters. 4.3.3 Subspace R efining: A Subspace D istance M easure The lesion and normal tissue subspaces estimated individually by the parametric or non-parametric methods from the known types of tissue data generally capture the characteristics of lesion and normal tissues, respectively. But from the physiological compartmental model analysis, the resulting two subspaces. containing the common input plasma subspace, may be so close to each other that they can hardly be separated. Hence, to improve the lesion detection performance (see Chapter 5). a subspace correlation method is used to refine the identified subspaces. 3S Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. The subspace refining procedure selects two subsets of basis vectors (column im ization between the two candidate subsets, subject to the condition that the LSE are less than given values. By the process of refining, the common input plasma function w ill be excluded from both lesion and normal tissue subspaces. By a subspace distance measure, the subspace refining can be described as: given the identified subspaces (H ) and (S). principal correlation coefficient between K and S [21]. and f[ and e> are the thresh olds set for achieving the fidelity criterion. vectors) from the identified subspaces (H ) and (S) based on subspace distance max- H '. S' = arg max distance {H . S} H.S (4.7) where distance {H .S } = \J 1 — r'f is the subspace distance, where /•[ is the largest 39 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Chapter 5 Matched Subspace Detector 5.1 Overview A fter identifying the lesion and normal tissue subspaces using the methods described in Chapter 4. how to incorporate tlie identified subspaces into statistical decision criteria becomes critical and will be addressed in this Chapter. The generalized like lihood ratio test (GLRT) is a standard procedure for solving the detection problem [55. 75. 76. 77. 78]. The GLRT is used to derive the matched subspace detector which is the general building block o f m ulti-rank matched filter in signal processing. In sonar signal processing, the matched subspace detector is also called a matched field detector [54]. Generally, the matched subspace detectors turn out to be a ratio of generalized energy detectors. In this Chapter, three types of GLRTs w ill be derived under different assumption of voxel spatial correlation. The first G LRT directly applies the algorithm devel oped in [55] by assuming that pixels are spatially uncorrelated. Both additive and replacement noise models are used to characterize the observations of large lesions 40 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. and small lesions in this approach. Considering that image pixel values are not statistically uncorrelated because each count contributes to all pixel values [7]. the second GLRT adds a pre-whitening procedure to the first. The spatial covariance m a trix required for the pre-whitening process is to be computed using the methods described in Appendix B [7. 30. 31]. Finally, the third detector extends the G LRT criterion to a m ulti-pixel detection. In this approach, it is assumed that the spa tia l inter-pixel correlation in each frame has the same structure but different energy level. In this dissertation, depending on the size of lesions, two data models are applied for the detection hypothesis, namely, the replacement model and the superimposed model. 5.2.1 Replacem ent M o d e l Based on a replacement model, which is used to model a large lesion in FD G -PET images, two possible hypotheses regarding how the PET data was generated are needed. The null hypothesis Ho says that the data consists of a sum of normal tissue signal x 0 and noise n0. The alternate hypothesis H i says that the data consist of a sum of lesion signal X[ and noise n[. That is 5.2 Hypothesized Data Model Ho : y = X 0 + n 0 H i : y = x t + n! (5.1) 41 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where the signal x, is assumed to obey the linear subspace model x 0 = s 0 . s e n Xxt. ( p e n 1. X[ = H 9. H € f t - Vxp. O e n p. t < X - p . 5.2.2 Superim posed M o d el For the small lesion case, partial volume effects w ill cause the lesion signal to be mixed w ith the background normal tissue signal, hence, it is best represented by a superimposed model. To model the small lesion case, the superimposed model is used for the hypothesis H, which are given by + n i where the linear model for the X[ (lesion signal) and n t are the same as those defined in the replacement model and (o.3) x 00 = S 0 O . 0 O 6 n l . under hypothesis / / 0 Xqi = S 0 [. 0 [ € n l . under hypothesis H\. (5.4) 42 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 5.3 GLRT: Spatially Uncorrelated Assum ption 5.3.1 Replacem ent M o d e l Based on the spatially uncorrelated assumption, the noise n, can be modeled as normal w ith zero mean and covariance m atrix a f I. Then, the detection problem for the replacement model becomes a test of the distributions: H 0 : y = x 0 -r n 0 : A '[S 0. crjl] Hi : y = X! + n! : .Y[H0.cr2I] The likelihood ratio test can be w ritten as l{9.cr2: y ) i(y) = fry) 'cr2\ ~X/2 ■> I I_ I o ’’ M ' - i l l i o _ 2 a$J { 2a-{ 2(Tq ' 2 \ - V / '2 I 1 1 G\ \ I 1 no 1 = l i « P - w " n ‘ N ? + ^ l i „ o ! l l • Hence, the generalized likelihood ratio test (GLRT) can be derived by substituting the maximum likelihood estimate (M LE) of n,. 9 and 0: /(y) - -o\ -.V/2 * • •) 1 1 I I I J O * '> ' 0 / K 2al l a Q The M LE. n,. 0. and 0. can be w ritten as 0 = (Sr S ) - 'S r y 9 = (H r H ) ' l H r y n0 = y — S 0 = y - S(STS)~lSTy — P ^ y 43 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. h t = y - H 0 = y - H ( H r H ) _ 1H r y = P ^ y (5.8) where P j-f is the orthogonal projection onto the complement subspacc of (H ). Be cause erf is unknown and the MLE of erf is &f = [|n ,! |"|. it is more convenient to replace the G LR T by the (.Y/2)-root G LRT. Hence, the G LRT for the replacement hypothesis model becomes _ J n oi I l y r P ^ y r / \ ri( M'-’ /.V “ 0 2 M y) = / y - = ii“ i ■ ] (5.9) y r P r/y ' 5.3.2 Superim posed M odel Under the spatially uncorrelated assumption, the detection problem for superim posed model is a test of the distributions: H q ■ y = x0 o + n0 : .Y[S0o. ffjjl] Hi : y = X i + x Ui + n i .Y [H 0 -t-S0,.crflj (5.10) Hence, the G LR T can be expressed as Ky) = * ( | ) + (5.11) The M LE n0 is the same as that derived in the replacement model, i.e.. “ o = P sT - (5-12) 4-1 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. For the h [. because the MLE of 0 and < f > l are given by = [H S ]# y = ([H S]r [H S])_l [H S]r y - 1 1 X i . s ' r . H r H H r S Sr H Sr S (H r P j H ) _ 1H r P ^ (Sr P r/S ) - l S7'P r/ (5.13) and from Eq. (A .5) in Appendix A. then n i = y - H 0 = y - H (H r P ? H ) - ‘ H r P ?y - S(STP f,S ) - lSTP r{y = y - E nsy ~ E.sv/y = p r/s-y (5.14) where Ens and S^n are called the oblique projections defined in Appendix A. There fore. the (.Y /2 )-root G LRT for the superimposed model becomes £s(y) = [/(y )]2/- v = y r P g y y r P /?sy' “ 0 !n i !l (5.15) Xote that the matrices H and S used to derived the above GLRTs for both re placement and superimposed models are assumed independent (after subspace angle refining). In Appendix C. we showed that the GLRTs can not be improved by fur ther projecting the lesion subspace (H ) onto the orthogonal complement subspace of 45 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. the normal tissue subspace (S). The GLRTs for the superimposed model are exactly the same before and after using the orthogonal refined subspaces. For the replace ment model. The G L R T is degraded after using the refined orthogonal subspaces. In summary, the G LR T can not be improved if two orthogonal subspaces are defined. 5.4 Spatially Correlated Assumption 5.4.1 Local S tatistics C o m p u tatio n It is known that image pixel values are not statistically uncorrelated because each count contributes to all pixel values [7]. Hence, a decorrelation for the pre-whitening process is im portant before applying the GLRT. The decorrelation procedure re quires the com putation of a covariance m atrix in the ROI. The covariance m atrix can be estimated by sim ply averaging the samples around the ROI. However, this sample-averaging m ethod to compute the local statistics can be poor since the local statistics assumption required may be invalid [27. 79]. Calculation accuracy can be guaranteed by using the formulas developed in [7. 19. 30. 31] which exploit the Poisson data model for a PET system. The linear and non-linear covariance com putation methods for the FBP and M AP reconstructed images, respectively, are described in Appendix B. and are used here to estimate the inter-pixel covariance, frame by frame. 4G Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 5.4.2 G L R T Because of the spatially correlated assumption, the average of TACs over an ROI may enhance rather than suppress the observation noise. Thus, to m itigate this, the TAC-averaging is preceded by a spatial decorrelation operation before applying the GLRT. This spatial decorrelation is usually called a pre-whitening procedure. In this subsection, the G LRT criterion is extended to m ulti-pixel G LR T detection. In this approach, it is assumed that the spatial inter-pixel correlation in each frame has the same structure but different energy level. Let the P-dimensional vector f, = [/,( 1). / t(2 ). • • •. f , { P ) } ' denote the spatial pixel data for the i-th frame and denote F = [ft. f>. ■ • •. fv]. We assume that all the columns of the m atrix F are independent. A Hence Y = [yi-Vo. • • •. y>] = F . The derivation of the m ulti-pixel GLRT for the constrained correlation structure is based on a replacement model, i.e.. 5.4.3 M u lti-P ix e l G L R T The GLRT can be derived as (see Appendix B): 1 + (G-l/2Rs)TP£(G-1/2Rs) 1 + (G-'/-> Rs)t P ^ ( G - 1 ^Rs) (5.17) where G = RRr — (Rs)(Rs) '. R is the covariance m atrix, and s = (s s)1/Js. sT = 1 TM 1/J. The m a trix M is the correlation structure for each frame. 47 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Chapter 6 Experiment and Result In this Chapter, experiments and results for dynamic feature identification and lesion detection in FD G -P E T dynamic study are presented. Both clinical and phantom dynamic data w ill be used in this dissertation. For the clinical dynamic data, two cases w ill be demonstrated, including one lung cancer dynam ic study and one breast cancer dynamic study. Each study contains one confirm ed prim ary lesion which w ill be used for lesion subspace estim ation, and one sm all lesion which w ill be used for the subspace fidelity and G LRT demonstration. For the breast cancer dynamic study, there is one "unknown" small lesion which can not be confirmed by eye. but is successfully enhanced after applying the G LRT. For the phantom data, two phantoms, one w ith five artificial lesions inserted and the other w ith only pure normal tissue background (no lesion), are generated. The tim e a ctivity curves (TAC) of a lesion and normal tissues observed from the clinical lung cancer dynamic data are used to simulate the dynamics for the phantom. To calculate the local statistics for the pre-whitening purpose, the methods [7. 19. 30. 31] described in Appendix B will be applied for covariance computation. A receiver operating characteristic (ROC) 48 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. study for the G LR T performance w ill be presented based on the dynamic phantom study. 6.1 Protocol of FDG -PET Dynamic Study The clinical FD G -P E T dynamic data was acquired w ith a Siemens/CTI EC AT Model 953A whole-body PET scanner. This device provides 31 contiguous transaxial image planes w ith an axial field of view of 10.8 cm. The nominal intrinsic resolution of the system is 4 mm in all 3 dimensions. Consecutive detector rings are separated by tungsten septa in order to reduce scatter noise. Dynamic data was acquired from 0 to 55 minute post injection. The final 10 minute frame can be used in place of the routine static study used in purely clinical examinations. The dynamic struc ture protocol for the clinical data acquisition is listed in Table 6.1. Based on this protocol, an FD G -P E T dynamic study contains 36 frames for a fixed plane position. Scan Type Scan Times Frame Duration Dynamic scan 0 - 2 ruin 15 sec/frame Dynamic scan 3 - 5 min 30 sec/frame Dynamic scan 6 - 25 min 1 m in/fram e Dynamic scan 26 - 45 min 5 m in/fram e Static scan 45 - 55 min 10 min Table 6.1: Dynamic data acquisition protocol 49 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 6.2 Clinical Study Two clinical FD G -PET dynamic studios, one lung cancer study and one breast cancer study, w ill be used in this dissertation. Each study contains one confirmed prim ary lesion which w ill be applied for lesion subspace estim ation and one small lesion which w ill be used for the subspace fidelity and G LR T demonstration. The last frame image reconstructed from the FBP method is shown in Figure G.l: (a) and (b) for the lung cancer study and the breast cancer image, respectively. An ROI. indicated by a rectangle, can be drawn over a segment of the tissue 011 the dynamic FDG-PET images and then duplicated to ail frames, so that the activity in that tissue can be tracked over time. The advantage of the ROI method over a single pixel is that it is less sensitive to noise and typically results in a smaller error in the parameter estimates. Data obtained through ROI analysis of dynamic images produces a tissue tim e activity curve (TAC ). This curve represents the counts/second/pixol (or counts/second/ml. if calibration data is available) in a given region as a function of time. Figure G.l: (c) and (d) shows examples of TAC for the two clinical dynamic studies. 6.2.1 Lung Cancer F D G -P E T D ynam ic S tudy In the lung cancer FDG-PET dynam ic study, two lesions, one prim ary (large and clearly visualized) and the other a metastasis (smaller and barely seen), in a lung cancer patient w ith m ultiple metastases. are indicated in Figure 6.2 (a). Three 5x5 pixel ROIs. L I. L2. and BG. indicated in Figure 6.2 (a), were selected from the 50 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 6.1: Clinical dynamic studies used in this dissertation: (a) the lung cancer image (FBP). (b) the breast cancer image (FBP). (c) one TAC example from the lung cancer dynamic study, (d) one TAC example from the breast cancer dynamic study. Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. prim ary lesion, the metastasis lesion, and the normal lung tissue, respectively. Their corresponding TACs are plotted in Figure 6.2 (b) where the TACs from the L i and L2 both have the trend to increase w ith time while the BG TAC does not. This demonstrates the property of FDG kinetics to show differential uptake in lesions and normal tissue. Figure 6.2: Test image reconstructed from FBP: (a) the last frame of the dynamic FD G -PE T FBP reconstructed images, where three ROIs were chosen: L I in the big lesion. L2 in the small lesion, and BG in the background, (b) the corresponding TACs for L l. L2. and BG. respectively. 6.2.1.1 Subspace Estimation The TACs observed from the L l and BG ROIs in Figure 6.2 (b) were used to form the temporal-spatial data matrices for identifying the lesion and normal tissue subspaces. (H ) and (S). respectively. For the parametric method, the LS search algorithm was used. For this lung cancer study, the TAC subspaces o f normal lung tissue and lung cancer estimated by (b) 52 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. the LS method were found to be spanned by the exponential basis {e °l‘.e °:t.e 0i'} and {1 — e~Jlt. te~j3t}. respectively, and the parameter sets were {$i = 0.00025. 9 -2 = 0.005. 9-i = 0.05} and { 3X = 0.005. J> = 0.0223. J:! = 0.00125} for BG and L l. respectively. For the non-parametric method, we applied the SYD method, where the data correlation m atrix was estimated using the tem poral-spatial data m atrix formed from L l and BG ROIs. The eigenvectors corresponding to the first < 7 1 and d l significant eigenvalues of the two estimated correlation matrices were chosen to span the TAC subspaces of the lung cancer and normal lung tissues, respectively. In our study, e /1 = e /2 = 2 . The accuracy (fidelity) of the lesion subspace estimated using the prim ary lung tum or was examined by projecting the observed L l and L2 TACs onto the estimated lesion subspace. Also, the accuracy of the normal tissue subspace estimated using the BG ROI was examined by projecting the observed BG TAC onto the estimated normal tissue subspace. On the other hand, the separability o f the subspace was examined by projecting the L l and L2 TACs onto the normal tissue subspace. and by projecting the BG TAC onto the L l subspace. The results from the LS method and the SYD method are shown in Figure G.3. It was demonstrated that the TAC of the m etastatic lesion can be represented by the subspace features extracted from the big lesion w ith a high accuracy in both methods, and that a large separability exists between the lesion subspace and the normal lung tissue subspace. 53 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. m ; a ) >c r «• k M * t [ XT O T - JL { ) \ f I } ! H r > c (a) (b) **. ‘vk> f j f i i i f i j IX •»£ ao (‘D (r) (f) Figure 6.3: Comparison of accuracy and separability of extracted subspaces of TACs in lesions and normal tissue by the SYD method: (a) TAC in L l represented by the normal tissue subspace, (b) TAC in BG represented by the lesion subspace, (c) TAC in L 2 represented by the lesion and normal tissue subspaces, respectively. Comparison of accuracy and separability of extracted subspaces of TACs in lesions and normal tissue by the LS method: (d) TAC in L l represented by the normal tissue subspace, (e) TAC in BG represented by the lesion subspace, (f) TAC in L2 represented by the lesion and background subspaces, respectively. 54 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 6.2.1.2 GLRT A test region which includes the small lesion and some part of the prim ary lesion was selected as indicated by a rectangle in Figure 6.4 (a). The G LR T (white noise assumption) w ith either a superimposed or a replacement data hypothesis model in Eqs. (5.9) and (5.15) was used to perform a m ultiple pixel test for each 5x5 pixel ROI contained in the test region. The corresponding G LR T result is shown in Figure 6.4 (c)-(e) which indicates that the matched subspace detection incorporated w ith the lesion and normal tissue subspaces extracted by the SVD method and the LS method significantly increases the lesion-to-normal tissue contrast for the small lesion relative to the original FBP image in Figure 6.4 (b). Another test region was selected to include the heart area in order to demonstrate suppression of activity in the heart by the GLRT. Figure 6.5 (a) shows the selected test region in the lung cancer dynamic study, where the prim ary lesion, the small lesion, and the heart area were all included. Figure 6.5 (b) shows the original 3-D mesh of the test region. The G LR T of the test region is shown in Figure 6.5 (c). where the heart activity is suppressed significantly after applying the GLRT. In order to apply the whitening process for the spatially uncorrelated FBP data, we used the method proposed by Carson [7] to calculate the spatial covariance ma trix. Carson's approximation formula assumes that the ROI is smooth and the pixel variances inside the ROI are sim ilar. Therefore, the covariance m atrix is just a func tion of the filter and the location o f the pixels. The advantage of Carson's formula is that it does not need access to the original projection data and is convenient for 55 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. } - : « > * i s - *> ■ » . * <#* (a) (b ) Figure G.4: (a) The last frame o f the dynam ic FD G -PET FBP reconstructed images, where a rectangle indicates a test region containing a part of the big lesion and the whole small lesion, (b) 3-D mesh o f the test region indicated in (a), (c) G LR T with replacement model for TAC subspaces obtained by the SVD method, (d) G LRT with replacement model for TAC subspaces obtained by the LS method, (e) G LRT with superimposed model for TAC subspaces obtained by the LS method. 56 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. C C X (c) Figure 6.5: (a) The test region indicated by a rectangle in the lung cancer study, where the prim ary lesion, the small lesion, and the heart area were all included, (b) 3-D mesh o f the original test region, (c) The result of the GLRT. • ) i Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. routine use. See Appendix B for more details. Figure 6.6 shows the results of the G LR T for the FBP image after applying the estimated covariance m atrix. Figure 6.6: (a) The test region indicated by a rectangle in the lung cancer study, (b) 3-D mesh of the original test region, (c) G LRT (white noise), (d) G LRT (non-white), (e) M ulti-pixel GLRT. 6.2.1.3 Subspace Refining Method For the proposed subspace refining algorithm described by Eq. (4.7). we estimated four LS TAC bases for both lesion and normal tissue subspaces, and used the sub space refining algorithm to retain only three bases for both subspaces. Figure 6.7 58 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. (a) and (b) shows the fidelity between the estimated TAC by the identified sub spaces and the observed TAC w ith out subspace refining (four bases) and Figure C.7 (c) shows the G LRT. The corresponding results after subspace refining (three bases retained) are shown in Figure 6.7 (d). (e). and (f). respectively, where not only the fidelity was preserved, but also the G LR T performance is significantly improved. (:U (b ) (c) (d) (e) (f) Figure 6.7: Subspacc refining using the subspace distance measure. W ithout sub space refining: (a) lesion TAC and (b) norm al tissue TAC both show the fidelity of subspaces estimated from observed TAC. (c) GLRT. W ith subspace refining: (d) lesion TA C and (e) normal tissue TA C shows the fidelity of subspace refined. It is clear that the refined subspaces m aintain the fidelity requirement, (f) G LR T performance is improved by using the refined subspaces. For the orthogonal projection, we already showed, see Appendix C. that the refined subspaces by the orthogonal projection can not improve the G LRT w ith the 59 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. original estimated subspaces. The experimental results are as follows. Figure 6.8 (a) and (b) show the test region chosen and its corresponding 3-D mesh. Figure 6.8 (c) and (d) show the GLRTs for the superimposed model before and after using the orthogonal refined subspaces, where the two GLRTs were exactly the same. Then for the G LRT of the replacement model. Figure 6.8 (f) shows that the G LR T using the refined orthogonal subspaces had poor performance compared to Figure 6.8 (e). where the original subspaces were applied for the GLRT. 6.2.1.4 M AP Reconstructed Image We also tested the G LR T for the lung cancer FG D-PET dynamic study by the MAP reconstructed algorithm s [42. 18. 43. 24]. For the non-white noise case, the fast covariance calculation method proposed by Jinyi and Leahy [31] for the M AP image was used here for the decorrelation purpose. This covariance calculation method was summarized in Appendix B. We used the same subspaces estimated from the FBP image. The results of G LRT are shown in Figure 6.9. 6.2.2 B reast C ancer F D G -P E T D ynam ic S tudy In this case, we examined a clinical. 36-frame, dynamic FD G -PE T study of a patient w ith a prim ary breast cancer and a known axillary metastasis. The prim ary cancer was clearly visualized in the left breast, while the smaller axillary metastasis lymph node was m arginally visualized, see Figure 6.10 (a)-(b). TACs obtained in the ROI placed on these tumors and contralateral normal breast tissue are shown in Figure 60 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 6.S: Subspace refining using the orthogonal projection, (a) test region inside the lung cancer image, (b) original 3D mesh. G LRT for the superimposed model: be fore. (c). and after, (d). the orthogonal subspace refining. G LR T for the replacement model: before, (e). and after, (f). the orthogonal subspacc refining. Cl Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. (c) Figure 6 9: The M A P reconstructed lung cancer FDG-PET image, (a) The test region indicated by a rectangular for the GLRT. (b) 3-D mesh of the original tost region, (c) GLRT: white Gaussian noise, (d) GLRT: non-white Gaussian noise. 62 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 6.10 (c). The TACs showed that FDG radiotracer uptake w ith in the lesions increases w ith time, while the TAC decreases for normal tissue. The subspaces for the breast lesions and the norm al tissues were identified using the SYD method. The subspace fidelity and separability were demonstrated by fit ting the metastasis lymph node ROI onto the subspaces estimated from the prim ary cancer ROI and the normal breast tissue ROI. Figure 6.11 (d) shows that the sub space of the lesion TACs is approximately invariant between different tumors and the subspaces o f the normal tissue and the lesion TACs were different. For a tost region, right rectangle in Figure 6.11 (a), containing the small lesion, applications of the GLRT w ith replacement model showed that the lesion-to-normal tissue contrast was enhanced relative to the original FBP image. Figure 6.11 (f). while for a normal tissue test region, left rectangle in Figure 6.11 (a), the G LRT shows no enhancement. An interesting GLRT result for this breast cancer dynamic study was the en hancement of an unknown lesion which could not be confirmed by eye. This un known lesion was roughly located (confirmed by a later scan) inside the rectangular region shown in Figure 6.12 (a) and (b). but was undiscernible to human eyes even using a M AP reconstruction method. This patient did not have a follow-up scan un til IS months later. From the newly acquired data, the "new" lesion was confirmed. The results of G LR T (whit-noise assumption) based on the earlier FD G -PET scan of this patient showed that this unknown lesion was detected successfully in both M AP and FBP images, see Figure 6.12 (e) and (f). when compared to the 3-D mesh of the test region, see Figure 6.12 (c) and (d). 63 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. $ f i } n r J X ■ v C (a) (b) (c) Figure 6.10: Breast cancer study: (a) a plane containing a prim ary cancer indicated by an arrow, (b) a plane containing a m etastatic lym ph node indicated by an arrow, (c) the TACs from the prim ary cancer, the m etastatic lymph node and normal breast tissues in an ROI indicated w ith a rectangle in (b). 6.3 Phantom Dynamic Study 6.3.1 C o m p u ter Sim ulated P h an to m D a ta Two computer simulated phantoms, one w ith five artificial lesions at the known locations and the other w ith no lesions, were used in the phantom dynamic study, see Figure 6.13 (a)-(b). The TACs o f the normal and malignant tissues measured from the clinical lung cancer FD G -PET dynamic study were applied to simulate the lesion and normal tissue kinetics behavior, see Figure 6.13 (c). D uring the forward projection o f the phantom dynamic image, the tota l count of projection data in each frame was scaled to approach the total count of the corresponding frame in a clinical dynam ic study. This w ill be called the low count dynam ic study in this research, while the total count of projection data for a high count study is twice of the total count in a clinical dynamic study. During the 64 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure 6.11: Breast cancer dynamic study: (a) two selected test regions: test region 1 containing the normal breast tissue and test region 2 containing a metastatic lymph node in the FBP image (same as (b)). (b)-(c) 3-D mesh of test region 1 and test region 2 indicated in (d). (d) TAC of the metastatic lym ph node represented by the lesion subspace and norm al breast tissue subspace identified by the SYD method, respectively, (e)-(f) the G LR T results on test region 1 and 2. respectively. 65 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. (a) (b) (c) (d) (f) Figure 6.12: A G LR T result of enhancing an "unknown" lesion in the breast cancer FD G -PET dynamic study. The test region containing the unknown lesion in (a) FBP and (b) MAP. The 3-D mesh of the test region: (c) FBP and (d) MAP. GLRT (white noise): (e) FBP and (f) MAP. 66 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. (a) (b) (c) Figure 6.13: Left: lesion phantom image w ith 5 lesions, m iddle: normal tissue phantom w ith no lesions, right: TACs for lesion and normal tissues. forward process, the sinogram data were corrupted w ith Poisson noise and blurred kernels to simulate the noise process and system resolution. FBP. MAP. and OSEM reconstruction methods were used for dynamic phantom image reconstruction. The images (last frame) reconstructed by the FBP. MAP. and OSEM methods and the corresponding observed TACs o f normal and lesion tissues in the phantom dynamic images are presented in Figure 6.14. where the reconstructed lesion images looked like normal tissue images and the artificial lesions were hardly confirmed. The GLRT results for the FBP dynamic images are shown in Figure 6.15 (high count data) and 6.16 (low count data), and the GLRTs for the M A P low count dynamic phantom images were shown in Figure 6.17. The G LR T results for both high count data and low count data improved the visualization o f the artificial lesion locations. 67 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. (j) (k) _ (I) Figure 6.14: Left column: lesion phantom image w ith 5 lesions, middle column: normal tissue phantom w ith no lesions, right column: TACs for lesion and nor mal tissues. From top row to bottom: original phantom. FBP. OSEM. and MAP reconstructed images. GS Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. (a) Figure 6.15: The G LRT for high count FBP phantom data: (a) a test region, (b)-(c) white noise assumption, (d)-(e) non-white noise assumption. 69 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. .•>*1 • V -\V<. ► (Act: ' ' *,'■<•.'/'•. v 1 .1 v > a j o (a) * < 1 .J -X (b) (d ) m i. (c) “ ! j4 * -c ", & a £ e e Figure 6.16: The GLRT for low count FBP data: (a) a test region, (b)-(c) white noise asumption. (d)-(e) non-white noise assumption. 70 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. (b) (c) (d) Figure 6.17: The G LRT for an M AP phantom image (low count data): (a) a test region indicated by a rectangular, (b): G LR T (white noise assumption), (c): G LR T (non-white noise assumption), (d): m ulti-pixel GLRT. Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. 6.3.2 R eceiver O perating C h aracteristic (R O C ) S tu d y An ROC study was applied to the phantom dynamic study for comparing the perfor mance of detectors. F ifty lesion phantom studies (w ith 250 known lesions) and fifty normal tissue phantom dynamic studies were conducted in a Monte Carlo simula tion. The rates o f true and false positive decisions made from the three G LR T tests were counted, respectively. The ROC curves were estimated by using the software package R O C K IT 0.9B developed by Professor C.E. Metz. University of Chicago. The resulting ROC curves from the FBP reconstructed dynamic phantom images shows that the G LR T with spatial decorrelation (w ith area under curve Az=0.92) are close to that for the GLRT w ith a white-noise (Az=0.93). see Figure G.18 (a). The ROC curve from the m ulti-pixel G LR T further improved the Az (=0.98). We also compared the GLRT performance for different reconstruction algorithms. The measured Az for our ROC curves can be viewed as a computer observers for dif ferent reconstruction algorithms. A reconstruction algorithm can be claimed "bet ter" if it achieves a higher Az. Figure 6.18 (b) shows that the Az (M A P ) > Az (OSEM) > Az (FBP) which can be interpreted as the iterative reconstruction meth ods (M A P and OSEM) outperform the FBP method not only in image visualization, but also in the G LR T detection performance. Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Figure C.18: ROC curves for the low count data sim ulation: (a) FBP: G LRT (white noise, non-white noise) vs m ulti-pixel GLRT. where .4. (m ulti-pixel G LRT) > .4- (G LRT w ith non-white noise) = = ,4; (GLRT w ith white noise), (b) G LRT w ith white noise assumption: .4; (M A P ) > .4; (OSEM) > .4; (FBP). Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Chapter 7 Conclusion In order to improve early detection of small lesions using FD G -PE T. computer-aided lesion detection algorithm s have been developed in this dissertation to assist visual inspection in tum or detection using dynamic FDG-PET. In our approach, we hy pothesized that the spatial and temporal metabolic features available from dynamic FD G -PET images are identifiable in a clinically practical way and useful for lesion detection. We showed that malignancies can be distinguished from normal tissues on the basis of their rates of FDG uptake in terms of time a ctivity curves, given that lesions and normal tissues often differ in the rate of radiotracer accumulation or disappearance. Since the rates of FDG uptake in tissues are often characterized by a set of physiological factors and the time activity curves can be represented as a linear combination of the physiological factors, we used the subspace spanned by the physiological factors as a key feature to distinguish normal tissues from malig nancies. By forming the dynamic FDG-PET images into useful spatial-temporal data matrices, we were able to estimate the tum or and norm al tissue subspaces by signal subspace fittin g methods that are widely used in array signal processing. 74 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Then, applying the generalized maximum likelihood ratio principle, we adapted and generalized the matched subspace detection techniques for small lesion detection in dynamic FDG -PET. We demonstrated that w ith the matched subspace detectors, the subspace feature(s) extracted from visible, large lesions can be employed to con firm suspected, but not unequivocally identifiable, small lesions. Results from both clinical dynamic FD G -P E T studies of patients w ith breast or lung cancer and dy namic phantom data showed that the physiological subspaces are straightforward to identify and distinct for known tumors and normal tissues, and are valuable for improving lesion diagnosis when combined w ith the matched subspace detection cri teria. Therefore, the detection methods proposed in this dissertation are promising for improving early detection of small lesions. Based on these results, it is our be lief that such detection methods can be used to reduce costly m ultiple diagnostic procedures and to guide surgical intervention. I o Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Bibliography [1] L. P. Adler, .J. P. Crowe. X. K. Al-kaisi. et al.. "Evaluation of breast masses and axillary lym ph nodes with [F-18] 2-deoxy-2-fluoro-D-glucose PET." Ra diology. 187. pp. 743-750. 1993. [2] J. Aliresaie. C. Xahmias. and M. E. Jernigan. "Xeural network based seg m entation of magnetic resonance images o f the brain." IE E E Nuclear Science Symposium, 3. pp. 1397-1401. 1995. [3] X. M. A lpert. et al., "Estim ation of the local statistical noise in emission computed tomography." IE E E Transactions on Medical Imaging, vol. MI-1, no. 2. pp. 142-140. Oct. 1982. [4] S. C. A m artur. D. Piraino. and Y. Takefuji. "O ptim ization neural networks for the segmentation of magnetic resonance images." IE E E Transactions on Medical Imaging, vol. 11. no. 2. pp. 215-220. 1992. [5] D. Barber. "The use of principal components in the quantitative analysis of gamma camera dynamic studies." Physics in Medicine and Biology. 25. pp. 283-292. 1980. [6] R. T. Behrens and L. L. Scharf. 'Signal processing applications of oblique projection operators." IE E E Transactions on Signal Processing, vol. 42. no. 6. pp. 1413-1424. June 1994. [7] R. E. Carson, et al., "An approximation form ula for the variance of PET region-of-interest values." IEEE Transactions on Medical Imaging, vol. 12. no. 2, pp. 240-250. June 1993. [8] S. R. Cherry and M. E. Phelp. "Imaging brain function w ith positron emission tomography." Chapter 8 in Brain Mapping: The Methods, edited by A. \Y. Toga and J. C. M azziotta. Academic Press. 1996. [9] H. S. Choi. D. R. Haynor. and Y. Kim . "P a rtial volume tissue classification of m ultichannel magnetic resonance images - a mixed model." IE E E Trans actions on Medical Imaging, vol. 10. no. 3. Sep. 1991. 70 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. [10] P. S. Conti, et al.. "Carbon-11 labeled alpha-aminoisobutyric acid. 2-deoxy- D-glucose. and thym idine as potential imaging agents for prostatic and renal malignancies." Surgical Form. vol. 36. pp. 635-637. 1985. [11] P. S. Conti, et al.. "M ultiple radiotracers for evaluation of intracranial mass lesion using PET." Journal of Nuclear Medicine, vol. 32. pp. 954-970. 1991. [12] P. S. Conti, et al.. "PET and [18F]-FDG in oncology: a clinical update." Nuclear Medicine and Biology, vol. 23. pp. 717-735. 1996. [13] T. M. Cover. "A n algorithm for m axim izing expected log investment return." [E E E Transactions on Information Theory, vol. IT-30, no. 2. pp. 369-373. Mar. 1984. [14] L. P. Crowe. L. P. Adler. R. R. Shenk. and J. Sunshine. "Positron emission tomography and breast masses: comparison w ith clinical, mammographic. and pathological findings." Annals Surgical Oncology, vol. 1. pp. 132-140. 1994. [15] Y. J. Cunningham and T. Jones. "Spectral analysis of dynamic PET studies." Journal of Cerebral Blood Flow and Metabolism. 13. 15-23. 1993. [16] D. Feng and X. Wang. "A computer sim ulation study on the effects of input function measurement noise in tracer kinetic modeling w ith positron emission tomography (P E T )." Computer Biology Medicine, vol. 23. no. 1. pp. 57-68. 1993. [17] R. H. Huesman. "A new fast algorithm for the evaluation of regions of interest and statistical uncertainty in computed tomography." Physics in Medicine and Biology, vol. 29. no. 5. pp. 543-552. 1984. [18] E. U. Mumcuoglu. R. M. Leahy. S. R. Cherry, and Z. Zhou. "Fast gradient- based methods for Bayesian reconstruction of transmission and emission PET images." IE E E Transactions on Medical Imaging, vol. 13. no. 4. pp. 687-701. Dec. 1994. [19] J. A. Fessler. "Mean and variance of im p licitly defined biased estimators (such as penalized maximum likelihood): application to tomography." IE E E Transactions on Image Processing, vol. 5. no. 3, pp. 493-506. Mar. 1996. [20] J. A. Fessler. "Approxim ate variance images for penalized-likelihood image re construction." IE E E Nuclear Science Symposium and Medical Imaging Con ference. Albuquerque. XM. 1997. [21] G. H. Golub and C. F. Van Loan. M a trix Computation, 3rd edition, the Johns Hopkins University Press, 1996. Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. [22] G. L. Brownell. T. F. Budinger. P. C Lauterbur. and P. L. McGeer. "Positron tomography and nuclear magnetic resonance imaging." Science, pp. G19-G2G. Feb. 19S2. [23] P. J. Green. "Bayesian reconstruction from emission tom ography data using a modified EM algorithm ." IE E E Transactions on Medical Imaging, vol. MI-9, no. 1. pp. 84-93. Mar. 1990. [24] T. J. Hebert and R. M. Leahy. "A generalized EM algorithm for 3-D Bayesian reconstruction from Poisson data using Gibbs priors." IE E E Transactions on Medical Imaging, vol. M I-8, no. 2. pp. 194-202. June 1989. [25] H. M. Hudson and R. S. Larkin. "Accelerated image reconstruction using ordered subsets of projection data.” IE E E Transactions on Medical Imaging. 1994. [26] C. C. Huang and X. Yu. "Computer-aided lesion detection w ith statistical model-based features in PET images." USC-Radiology Inte rna l Report. July 1996. [27] C. C. Huang. X. Yu. J. Zheng. .J . R. Bading. and P. S. C onti. "A lgorithm to reduce the com plexity of local statistics com putation for P E T images." SPIE Proceedings on Medical Imaging, vol. 3034. pp. 1094-1111. Newport Beach. CA. Feb. 1997. [28] C. C. Huang. X. Yu. J. R. Bading. and P. S. Conti. "Computer-aided le sion detection w ith statistical model-based features in P E T images." IE EE Transactions on Nuclear Science, vol. 44. no. 6. pp. 2509-2521. Dec. 1997. [29] C. C. Huang. X. Yu. J. Bading. and P. S. Conti. "Feature extraction by subspace fitting of tim e activity curve in PET dynamic studies." IE EE Nu clear Science Symposium and Medical Imaging Conference. Albuquerque. XM . Nov. 1997. [30] R. H. Huesman. "A new fast algorithm for the evaluation o f regions of interest and statistical uncertainty in computed tomography." Physics in Medicine and Biology, vol. 29. No 5. pp. 543-552. 1984. [31] J. Qi and R. M. Leahy. "Fast com putation of the covariance o f M AP recon structions of PET images." Proceedings of SPIE. San Diego. CA. 1999. [32] A. C. Kak and M. Slaney. Principles of Computerized Tomography Imaging. IEEE Press. 1987. [33] K. Kubota, T. Matsuzawa. et al.. "Differential diagnosis o f lung tum or w ith positron emission tomography: a prospective study," Journal of Nuclear Medicine, vol. 31, pp. 1927-1933. 1990. 78 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. [34] J. Llacer. et at.. “ Results of a clinical receiver operating characteristic study comparing filtered backprojection and maximum likelihood estim ator images in FDG PET studies." Journal of Nuclear Medicine, vol. 34. pp. 1198-1203. 1993. [35] G. Lucignani. K. C. Schmidt. R. M. Moreco. et al.. "Measurement of regional cerebral glucose utilizatio n w ith fluorine-18-FDG and PE T in heterogeneous tissues: theoretical considerations and practical procedures." Journal of Nu clear Medicine, vol. 34. No. 3. 360-369. Mar. 1993. [36] S. R. Meikle and M. Dahlbom. "Positron emission tomography: application in oncology." in Nuclear Medicine in Clinical Diagnosis and Treatment, edited by I. P. C. M urray and P. J. Ell. vol. 2. pp. 1327-1337. [37] S. R. Meikle. J. C. Matthews. V. J. Cunningham, et al.. "Spectral analysis of PET projection data." IE E E Nuclear Science Symposium and Medical Imaging Conference, p. 1888-1892. 1997. [38] C. E. Metz. "Basic principles of ROC analysis." Seminars in Nuclear Medicine, vol. 4. pp. 283-298. 1978. [39] M. M. Ter-Pogossian. M. E. Raichle. and B. E. Sobel. "Positron emission tomography." Scientific American. 243 (4). pp. 171-181. 1980. [40] J. C. Mosher. P. S. Lewis, and R. M. Leahy. "M ultiple dipole modeling and lo calization from spatio-tem poral MEG data." IE E E Transactions on Biomed ical Engineering, vol. 39. no. 6. pp. 541-557. June 1992. [41] J. C. Mosher and R. M. Leahy. "Recursive MUSIC: a framework for EEG and MEG source localization." IE E E Transactions on Biomedical Engineering. vol. 45. no. 11. pp. 1342-1354. Nov. 1998. [42] E. U. Mumcuoglu. R. M. Leahy, and S. R. Cherry. "Bayesian reconstruction of PET images: quantitative methodology and performance analysis." Journal of Computer Assisted Tomography. Jan. 1994. [43] E. C . Mumcuoglu. R. M. Leahy, and S. R. Cherry. "Bayesian reconstruction of PET images: methodology and performance analysis." Physics in Medicine and Biology, vol. 41. no. 9. pp. 1777-1807. Sep. 1996. [44] K. S. X ijran and D. Barber. "Towards automatic analysis of dynamic radionu clide studies using principal-component factor analysis." Physics in Medicine and Biology. 30. pp. 1315-1325. Dec. 1985. [45] K. S. Xijran and D. Barber. "Factor analysis of dynam ic function studies using a priori physiological inform ation." Physics in Medicine and Biology. 31. pp. 1107-1117. Oct. 1986. 79 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. [46] D. Ortendahl. "The application of principal component analysis to m ultivari ate M R I data." Proceedings of the 8th Annual Conference of the Engineering in Medical and Biology Society, p. 1065-1068. Fort W orth. TX. 1986. [47] F. O'Sullivan. "Imaging radiotracer model parameters in PET: a m ixture analysis approach." IE E E Transactions Medical Imaging, vol. 12. no. 3. pp. 399-412. 1993. [48] R. M aitra and F. O'Sullivan. "Estim ating the variability of reconstructed PET data: a technique based on approximating the reconstruction filter by a sum of Gaussian kernels." IE E E Nuclear Science Symposium and Medical Imaging Conference, vol. 3. pp. 1411-1414. 1995. [49] F. Pedersen. M. Bergstrom. E. Bengtsson. and B. Langstrom. "Principal com ponent analysis of dynamic positron emission tomography images." European Journal of Nuclear Medicine, vol. 21. no. 12. pp. 1285-1292. Dec. 1994. [50] M. E. Phelps. J. C. M azziotta. and H. R. Schelbert. Positron Emission To mography and Autoradiology. Principles and Applications fo r the Brain and Heart. Raven Press. 1986. [51] I. S. Reed and X. Yu. "Adaptive m ulti-band CFAR detection of an optical pattern w ith unknown spectral distribution." IE EE Transactions on Acoustic. Speech, and Signal Processing, vol. 38. no. 10. Oct. 1990. [52] L. E. Brennan and S. Reed. "A digital adaptive antenna for JTID S ." subm it ted to the Electronic Systems Division. AFSC. USAF. Hanscom AFB. MA 01731. Dec. 6. 1989. [53] P. Santago. and H. D. Gage. "Q uantification of M R brain images by m ix ture density and partial volume modeling." IE E E Transactions on Medical Imaging, vol. 12. no. 3. Sep. 1993. [54] L. L. Scliarf. Statistical Signal Processing. MA: Addison-Wesley. 1991. [55] L. L. Scharf and B. Friedlander. "Matched subspace detectors." IE E E Trans actions on Signal Processing, vol. 42. no. 8. pp. 2146-2157. Aug. 1994. [56] B. Schmall. P. S. Conti, D. J. Schaeffer, and E. L. Kleinert. "Tum or and organ biochemical profiles determined in vivo following uptake of a combination of radiolabeled substrates: potential applications for P E T.'- American Journal of Physics Imaging, vol. 7. pp. 2-11. 1992. [57] K. Schmidt. G. Mies, and L. Sokoloff. "Model of kinetic behavior ofdeoxyglu- cose in heterogeneous tissues in brain: a reinterpretation of the significant of parameters fitted to homogeneous tissue models," Journal of Cerebral Blood Flow and Metabolism, vol. 11, pp. 10-24. 1991. 80 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. [58] S. C. Huang. M. E. Phelps. E. J. Hoffman. K. Sideris. C. J. Selin, and D. E. Kuhl. "Xoninvasive determ ination of local cerebral m etabolic rate of glucose in man." American Journal o f Physiology. E69-E82. 1980. [59] L. A. Shepp and B. F. Logan. "The Fourier reconstruction o f a head section." IE E E Transactions on Nuclear Science. NS-21. pp. 21-43. June 1974. [60] L. A. Shepp and Y. Yardi. "M axim um likelihood reconstruction for emission tomography." IE E E Transactions on Medical Imaging. M I-1, no. 2. pp. 113- 122. Oct. 1982. [61] Z. H. Cho. .J. Jones, and M. Singh. Foundations of Medical Imaging. John W iley A' Sons. Inc.. New York. 1993. [62] L. G. Strauss and P. S. Conti. "The applications of PE T in clinical oncology." Journal of Nuclear Medicine, vol. 32. no. 4. pp. 623-648. A pr. 1991. [63] T. X -L Tan. M. Y. Klein. V. T. Tabib. C. Turner, and J. M. Forde. "Oncology applications of PET: an updated review." Applied Radiology, p. 18-25. .June 2000 . [64] E. Tanaka and H. Murayama. "Properties of statistical noise in positron emission tomography." IE E E Proceedings of International Workshop Physics Engineering Medical Imaging, p. 158-164. New York. 1982. [65] G. D. Tourassi and C. E. Floyd. "A rtificia l neural networks for single pho ton emission computed tomography: a study of cold lesion detection and location." Investigative Radiology, vol. 28. no. S. pp. 671-677. Aug. 1993. [66] F. Turkheimer. R. M. Moresco. G. Lucignani. et al.. "The use o f spectral anal ysis to determine regional cerebral glucose utilizatio n w ith positron emission tomography and [1 8F]-fiuorodeoxyglueose: theory, im plem entation, and opti m ization procedures." Journal of Cerebral Blood Flow and Metabolism. 14. 406-422. 1994. [67] Let's Play PET. http://ww w.cnim p.ucla.edu/lpp. [68] Y. Yardi. L. A. Shepp. and L. Kaufman. "A statistical model for positron emission tomography." Journal of American Association, vol. 80. no. 389. 8-20, M ar. 1985. [69] M. Yiberg and B. O ttersten. "Sensor array processing based on subspace fittin g ." IE E E Transactions on Signal Processing, vol. 39, no. 5. pp. 1110- 1121, May 1991. SI Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. [70] R. L. Wahl. "Positron emission tomography: application in oncology." in Nuclear Medicine in Clinical Diagnosis and Treatment, edited by I. P. C. M urray and P. J. Ell. vol. 2. pp. 801-820. [71] R. L. Wahl. G. D. Hutchins. D. J. Buchsbaum. M. Liebert. H. B. Gross man. and S. Fisher. "18F-2-deoxy-2-fluoro-D-glucose uptake into human tu mor xenografts: feasibility studies for cancer imaging w ith positron-emission tomography." Cancer, vol. 67. pp. 15-14-1550. 1991. [72] R. L. Wahl. R. Cody. G. D. Hutchins, and E. Mudgett. "Prim ary and m etastatic breast carcinoma: in itia l clinical evaluation w ith PET w ith the radiolabeled glucose analog 2-[F-18]-fluorodeoxy-2-D-glucose(FDG)." Radi ology, vol. 179. pp. 765-770. 1991. [73] D. W. Wilson and M. W. Benjamin. "Noise properties of filtered- backprojection and M L-EM reconstructed emission tomographic images." IE E E Transactions on Nuclear Medicine, vol. 40. no. 4. pp. 1198-1203. Aug. 1993. [74] X. Yu. I. S. Reed, and A. Stocker. "Performance comparison of adaptive m ulti-band optical target detection." IE E E Transactions on Signal Process ing. Aug. 1993. [75] X. Yu. et al.. "A linear feature mapping framework for adaptive detection and recognition of targets in complex SAR data." IE E E 28th Conference on Signal. System. and Computers. Monterey. CA. Oct. 30-Xov. 2. 1994. [76] X. Yu. J. R- Bading. and P. S. Conti. "Autom ated lesion detection with feature analysis in positron emission tomography." IE E E Nuclear Science Symposium and Medical Imaging Conference. Norfolk. YA. Nov. 3-6. 1994. [77] X. Yu. A. M. Chen, and I. S. Reed. "A maximum likelihood detection of signals using feature mapping framework." Proceedings of IE E E ICASSP. D etroit. M I. May 1995. [78] X. Yu and I. S. Reed. "Adaptive detection of signals with linear feature mappings and representations." IE E E Transactions on Signal Processing, vol. 43. no. 12. Dec. 1995. [79] X. Yu, C. C. Huang. R. M. Leahy. J. R. Bading. and P. S. Conti. "Analytic com putation and sample estimation of the local statistics of positron emission tomography images: a comparative study." Proceedings of the f'ir d Annual Meeting of the Society of Nuclear Medicine. Denver. CO, June 3-5. 1996. [SO] X. Yu and C. C. Huang. "Lesion detection w ith statistical model-based prior inform ation in positron emission tomography." IE E E Nuclear Science Sym posium and Medical Imaging Conference. Anaheim. CA, Nov. 2-9. 1996. 82 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. [81] X. Yu. L. E. Hoff. I. S. Reed, et a i. "Autom ated target detection and recog n ition in m ulti-band imagery: a unified M L detection and estimation ap proach." IE E E Transactions on Image Processing, vol. G . no. 1. pp. 134-156. Jan. 1997. [82] X. Yu. C. C. Huang. J. R. Bading. and P. S. Conti. "Use of a matched subspace filter for lesion detection in dynamic positron emission tomography." A nnual Meeting of Society of Nuclear Medicine. San Antonio. TX . June 2-5. 1997. 83 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Appendix A Oblique Projection Lot (H ) denote the subspace spanned by the m atrix H . If the subspaces (H ) (/>- dimension) and (S) (f-dimension) are linearly independent, then the columns of the concatenated m atrix [HS] span a (p + O-dimensional subspace (H S). The typical orthogonal projection of a TAC y (G Tts ) onto (H S ) is denoted by P usY- where P u s is the orthogonal projection in the subspace (H S) and represented as P/,.- = [HS] ([H S]'r [H S ]) " ‘ [H S ]r . We can further re-express P //5 as (A. 1 1 H S H r H H 7 S Sr H Sr S _1 1 1 33 H 1 ST = [ h o + [ 0 s H r H H r S H r Sr H s r s s r H r H H r S -I H T Sr H STS . ST . = E hs + E. ■ ‘SH (A.2) 84 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where E HS = A H 0 H r H Hr S -1 Hr Sr H s r s . sT H r H HTS -I Hr S7H s r s ST Esii = O S The two projections. E [{S and E ^//. are called the oblique projections [0] with the property that E //s and E Sf! have respective range spaces (H ) and (S) and respective null spaces (S) and (H ). i.e.. E //5 H = H . E/f.^S = 0 . E>'//S = S. E , „ H = 0. (A.3) The sim plification version of the oblique projections E us and E sn can be w ritten as E lls = H ( H P jH ) H P j E s „ = S(S t P T[S ) - 1St P Tr (A.4) Hence. P HS — E /,5 + E s „ = H ( H r P £ H r l H r P£ + S{STP f , S ) - lSTP f,. (A .5) 85 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Appendix B Local Statistics in PET Images B .l Local Statistics Calculation in PET FBP images The analytic com putation (AC) for the local statistics of a P E T FBP image was originally developed by Huesman [17]. By the independence of observation data and the linear properties of FBP. AC estimates the (co)variance of P E T image pixels from the projection data statistics via some scaling factors. Huesman's algorithm calculates local statistics directly from the projection data w ith o u t image recon struction. Recently. Carson et al. [7] developed an approxim ation formula for the variance of PET ROI values, which accounted the radioactivity d istribu tion, atten uation. random, scatter, deadtime, detector normalization, scan length, decay, and reconstruction filter. By the assumption that the pixel variances o f any two pixels are sim ilar, and that the product of attenuation, normalization, wobbling and raw data is relatively uniform for those projection lines that substantially contribute to SG Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. the sum. they derived a simple formula for the correlation coefficient of any two pixels, which is independent of the raw data. Let the raw projection data in PET be an /-dimensional column vector p = [p f.P a PnJT- "'here I is the product o f the total number of angles (tig) and the to ta l number of rays •nr . and p# (size: n r x 1 ) is the sinogram data in the angle 0. Let the corrected data of p be p' which consists of p'o, = .-W (.\Vr(U; • po.r ~ Ro.r)-So.r). 0 = 1.2.......n0. r = 1.2 nr . (B.l) The terms A 0.r and .\'o.r in Eq. (B .l) are m ultiplicative correction terms for atten uation and detector efficiency (norm alization). U’r corrects for the fraction of time spent in each wobble position and is generally independent of angle. R0.r and So.r are the estimates of random and scatter. The terms p0,r . p'0 r. Ro.r• and S0,r have units of counts, while all other terms are dirnensionless. Therefore, for any pixel gt w ith index i at position .r,. i/t. the pixel value calcu lated by the FBP is the convolution of the corrected projection data p'0r with the reconstruction filter h. that is. n« nr 9, = ^ 2 ^ h ( . v l cos0 + y.sinfl - r)p'0 r. (B.2) 0=1r=l In practice, additional linear interpolation steps are applied and the reconstructed pixel value becomes ft* Tlr 9i = 5 Z 5 Z ho!rPo.r ( B -3 ) 0=1r=l S7 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. whore h# \ is unique for the pixel i in FBP at angle 9 and ray r. B .1.1 V arian ce o f P ixel in an R O I From Eq. (B.3). the variance of the pixel i can be sim ply w ritten as Kt) Tlr Var{<7,} = E Z K ' r ) ' • Var{ft>.r} (B.4) 0=1r=1 In a PET system, it is usually assumed that the significant source of noise in p'g r is due to the Poisson statistics in the raw projection data po,r . that is. the noise in the term Ag,r- N 0 r . U'r . Rg.r • and 5«.r is small. Also, the raw count po.r are statistically independent Poisson distribution, so are the p'0 r by the Eq. (B .l). Therefore, tlie pixel variance in Eq. (B.-l) can be approximated as 0=1r=l We can further assume that for each angle 9. the product of p0,r and [Ay r .Y,jr ll •] is relatively uniform in the ROI. i.e.. we can approximate [A#r .V |r l l '^\po.r by some 0=1 r = I n # Tlr representative of them, namely [Ajj.Yy U 'jp o • Therefore. Eq. (B.5) becomes 0=1 r= 1 "(J • 1 0 = 1 n o SS Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Po Po (B.6 ) 6=1 where ||h||-> is defined as the 2 -nonn of the filter vector h and ho is the filter vector associated w ith the angle 6. We used the observation that the 2-norm of the filter vector ho in all angle are approxim ately equal in an ROI. that is. ||h||-> % ||ho|!'_>- for all 0. Hence, the value ct = ||hj[ij can be determined a p rio ri for the convolution filte r used in the FBP. B .1 .2 Covariance o f Pixels in an R O I Following the sim ilar steps in deriving the variance of a pixel, the covariance of any two pixels can be computed as r ip n r Cov{gj.gk) = Y . Y , h d .rh o.r ' Var (Po) 0 = t r = l f i d T lr 6 = 1 r = l PO.r (B.i where g} and gk are two pixel values of interest w ith reconstruction filters li{ g j)r and hgkj.. respectively. The correlation coefficient of pixels g} and gk is given by P i 9 j - 9 k ) = Cov(gj.gk) y/VarfSj) ■ V ar(st) £ £ P O , > = t r = l L J Po., (B.S) 89 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Assume that the pixel variances of and gk in an ROI are sim ilar, and that the prod uct A j r Xg r \Vf ■ poT is relatively uniform for those projection lines that substantially contribute to the sums. Eq. (B.8 ) can be further sim plified as Z Z P (9r9k) = ' = lr= l e e ( C ) )=[ r— 1 7 J 10=1 r = I ' 7 ng n r 10= n ti Tlr , • ,i\ e e ! & » ; _ f l = l r = l ____________ rir Z E(/»fl.r)2 0=1 r= 1 = ^ < B -9 > 0=1 r=l A where c _ > = l/( c [ • n6). Thus, the correlation coefficient of any two pixels in an ROI is independent of the raw data and Eq. (B.9) can be used for a set o f pixel pairs at various positions to produce a p(d) table. B.2 Covariance Computation in M AP images The method used in this dissertation to compute the covariance o f the MAP re construction image was originally developed by Fessler [19. 20]. Recently. Jinyi and Leahy [31] proposed a fast algorithm to reduce the com putational complexity in Fessler's method. The following is a brief review of the fast covariance computation of M AP image described in [31]. 90 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. B.2.1 M A P R econstruction A PET system can be modeled as independent Poisson random variables y w ith mean y related to the emission intensity image x through an affine transform y = P x + n (B.10) where P is the detection probability m atrix and n accounts for the presence of noise in the data. The detection probability m atrix is modeled as P = D [n,]G . where D[nt] is a diagonal m atrix containing the correction factors, n,. i.e. the product of the detector norm alization, dead-time and attenuation correction factors'. G is the geometric projection m atrix representing the probability that an emission from each voxel in the image produces, in the absence of attenuation effects, a photon pair at each of the detector pairs in the system. The log-likelihood function for the Poisson data model is L{'Jk) = 5 1 & lo S f/‘ ~0i + lo S Ui- (B .ll) I The MAP reconstruction is usually described as the maximizer of the log posterior probability: x(y) — a rg m a x I(y |x ) - JU(x) (B. 12) where J is the hyper-parameter that determines the relative influence of the prior and likelihood terms. L(y\x) is the log-likelihood function for the Poisson data model. 'D[n,] denotes a diagonal matrix with the (t.t)th diagonal element equal to nt. 91 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. and i'(x) is the Gibbs energy function depending on a sparse neighborhood m atrix C B .2.2 A p p ro xim a tio n of C ovariance Based on a fixed point of the M AP objective function and using a truncated Taylor series expansion of the im plicit estimator. Fessler [19. 20] developed a closed form for the com putation of covariance in the M AP reconstruction image where F = P 'D [l/t/,]P is the Fisher inform ation m atrix. This closed form involves a com putation of the inverse of a Hessian m atrix or solving a related set of linear equations. Hence, to reduce the com putation com plexity. Jinyi and Leahy [31] pro posed a fast com putation method which uses the follow ing approxim ation for the Fisher inform ation m atrix Cov(i-) = [F -F JC'C]-‘F[F + JC'C]"1 F ^ DkG'GD, ( B . 1 4 ) where D * = D[kj] w ith £ L , a in 't/9 i K — ------- --------- ( B . 1 5 ) Then. C o v(i) £ D ^ B t x r ’ G 'G B f x r ’ D ; 1 (B. 1G) where B(x) = G'G -F JD ^C 'C D '1. (B.17) 92 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. The covariance can be further approximated using the following observations: (i) the m atrix B(x) represents a local image blurring operator since the prior energy function is defined on a local neighborhood and the blurring function G'G is also local: (ii) the correlation between voxels drops off rapidly as a function of the distance between them. The covariance w ith respect to voxel j is therefore dominated by the contribution of Kj and the following approximation can be made: Covj(j) a D~lB(x)_1G'GB(x)_lD"leJ a ^ K I J K y r ' G ' G K f J f y - r ' e j (B.18) where K (JkJ2) = G'G + Jk;-C'C. (B.19) Since G'G and C'C correspond to the shift invariant blurring operators, they have a block Toeplitz structure and can be approximately diagonalized using the 2D discrete Fourier transform , i.e.. G'G Q 'D [A,]Q C'C ^ Q 'D [//,]Q (B.20) where Q and Q '. respectively, represent the Kronecker form o f the 2D D F T m atrix and its inverse. A/s are the 2D Fourier transform of the system geometric response 93 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. at the central voxel, and //,'s are the 2D Fourier transform of the central column of C 'C . By Eq. (B.20). Eq. (B.19) can be w ritten as In cases where only a small ROI is considered, the correlation structure can be assumed to be invariant w ithin the region. Hence the covariance m atrix can be computed by one modified backprojection to compute the k/ s and one 2D FFT. One advantage o f the above derivation is that A, and /z, are independent of the data and can be pre-computed. Another advantage of Eq. (B.22) is that it is readily inverted which is an attractive property for the pre-whitening process of the non-white PET data for computer observer lesion detection. Then, the final covariance approximation formula becomes 94 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Appendix C GLRT Can Not Be Improved by Orthogonal Subspace Refining Given two identified lesion subspace (H ) and normal tissue subspace (S). the derived GLRT in chapter 5 can be w ritten as y rPjy L ,(y ) = " (superimposed model) v r P * v L r(y) = — ^—-— (replacement model). (CM) y ' P n y The refined subspace matrices by the orthogonal projection of (H ) onto the comple ment subspace of (S) can be w ritten as S' = S. H ' = P jH (C.2 ) where P £ represents the orthogonal projection onto the complement subspace of (S). Note that S' is orthogonal to H '. Then we have the following properties: H ' = P 5 H = H - P 5 H = » P // = P / / - P p,h (C.3) 95 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. P //'S ' = PH' + P^' = Ppj- H 4' Ps’ = Pp*H + P5 — PHS- S1 Superimposed model From the property in Eq. (C.4). the G LRT for superimposed model can be w ritten as My) = 0 ^ = = M y) < C 5> y ^ h's'j y ^//sy wherethe G LRTs before and after the orthogonal projection of the lesion subspace onto the orthogonal complement of the normal tissue subspace are shown to be equal. Replacement model The G LR T for the replacement model using H' and S' is i ' r ( y ) = ( C . G ) y r P fr y Then, by using the property in Eq. (C.3). we have 1 = y r P ^ y £ 'r(y) yrPjf-y _ y T [Ps + Pp.-h yrP?-y 1 , y r P p ,H y ,r . . H ------ Hence, M y) yrP5-y My) = , U,.Ky < My). lcs) M y ) + y r P j,y 90 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. Eq. (C.8 ) shows that, for the replacement model, the performance of the GLRT after applying the orthogonal projection refined subspaces is inferior to that using y T p p H y the original subspaces, because the scaler. yf p j - - . is positive. 97 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. A ppendix D M u lti-P ix e l G L R T D erivatio n Given an ROI w ith P number of pixels, let the P-dirnensional column vector f, = [ /, (1 ). / , ( 2 ). • ■■.fi{P)\1 denote the spatial pixel data for the /-tli frame and let the ROI data Fp*.v = [f\. f>. • • ■ . f.v]. We assume that all the columns of the m atrix F are independent. Let Y \ , r = [yi-y-j- • • -.y/*] = Fr. The hypothesis data model for the ROI under the replacement model assumption is Ho : Y = Xo -r No ^ ^ H x : Y = X 1 + N 1 where the noise N , is assumed to be a m ultivariate norm al (M Y X ) w ith zero mean. It is also assumed that the signal X , obeys the linear subspace model, i.e.. Xo = S$ 0 X , = (D.2) 98 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where H and S are the lesion and normal tissue subspaces, respectively, and is the corresponding coefficient m atrix. Hence, the mean of tem poral-spatial ROI data becomes _ f S $ 0. under hypothesis H 0 ^ ^ [ H ^ i. under hypothesis H We assume that the spatial covariance m atrix for each frame has the same structure but w ith different energy level [7]. That means the spatial covariance m atrix for the i-th frame can be w ritten as cr,M. i.e.. cov(f,) = cr,M. i = 1. ■ • •. A" (D.4) where M is the spatial inter-pixel correlation m atrix for all frames. Let Z = M 1 2F. then the variance m atrix for the whitened data Z is var(Z) - (T i (J y ' ‘ ' C T Y <71 (7j ‘ ’ ( 7 S O i <7> • • • ( 7 v . 1 ‘ ' . P < V :d .5) Let R = [r[. r-j. • • •. rp] = Zr. Then the covariance m atrix for each column of R is cov(r,) = cr i 0 • • • 0 0 a2 • • • 0 0 0 ■ • • crv = A. ID.G) V X .V The ROI data w ill be statistically described as a Gaussian probability density func tion (PDF) based on the temporal-spatial data m atrix R. A fte r the whitening 99 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. process, by the independence of the columns r, in R. the jo in t Gaussian PDF of R under hypothesis, say H \. becomes /i(R) = / ( r i. r - 2 . • • •. vp\Hx) = nf=1/(r,| H {) = nf=1(2-)‘ v ' - > iA1r 1/2exp [r, - E(r,)]r A"1 [r, - E(r,)]j = (2-)-vp-2!A1r P ' ie x p | ^ ^ [ r , - Efr,)]7 " A fl [r, - E(r,)j J = (-2-)- v,’,-|A,!-f’ 2ctr [R - E(R)f A f1 [R - E(R)J j = (2-)_A7’‘'!Aii“fJ 2ctr [R - E(R)][R- E(R}]rj (D.7) where "etr" denotes "exp trace", and E (R ) = E (Z 7 ') = E ( J M " 1 2F] j - E (F r M - ‘ 2) = E ( Y ) M _I'"" = H * 1M “ 1/2. (D.8 ) If we further assume that the mean o f image is the same for all frame, i.e.. (D.9) where <f}{ is the coefficient vector for all frames and l r = [1.1. • • •. 1]. Then. E (R ) = H 0 1l r M - 1 2 = b tsr (D.10) 100 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where b[ = and sr = l rM l/2. Hence, the PDF for R becomes / ( R) = (2 ;r)-vp/2 |A I | - p/2etr j ^ A f 1 (R - blSr ) (R - blSr)T} . (D .1 L) In order to apply the generalized likelihood ratio test (G LR T) for the matched subspace detector, the next step is to find the maximum likelihood estimate (MLE) o f A [ and (f>l . Based on (D. 11). the M LE of A t can be sim ply w ritten as A , = i (R - bls7')(R -b 1 sr )r . (D. 1 2 ) The M LE of (f>x can be derived from the m inim ization of jF b l|. he-. L (R ) = min |Ft,, j (D.13) t> i where F bl = ( R - b 1sr ) ( R - b 1 sT)r (D.14) Xow expand F bl as follows: F bl = RRr - R s b [ - b 1 srRr + b1 srsb[. (D.15) Let s = (sr s )l/2s and b [ = (sr s)1/' - > b 1 . so b[Sr = b[Sr . A fte r the change of variables. F6i = RRr - R s b [ - b 1 srRr + b1b[ = (b[ - Rs)(bt - Rs)r + RRr - (Rs)(Rs)r 101 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. = (b [ — Rs) (b i — Rs) ^ + G (D. 1 G ) where G = RRr - (Rs)(Rs)T and G can be shown to be invertible [52]. Hence, |F6J = |G|-|I + (b1-R s)(b l - R s ) rG - I| = |G] • |I + aarG_I| (D.17) where a = bi — Rs. Now m ultiply the m atrix D = I -t- aarG _1 by a as follows: Da = a + a(arG_1a) = (1 -+ - arG _1a)a (D.18) which shows that a is an eigenvector of D w ith the associated eigenvalue l+a7G_1a. Let v,. for ; = 1. 2. ■ • •. .V - 1. be some orthogonal set of vectors, all perpendicular to the vector G-Ia. Then. Dv, = (I-f aarG~‘)v, = v, + a (G_1a)rvIj = v, (D.19) which means v, is an eigenvector of D w ith eigenvalue equal to one. for / = 1 . 2. • • •. .Y- 1 . So Eq. (D.16) becomes |FbJ = |G| • (1 + arG_1a) = |G| • [l + (b 2 - R s ) rG -1(b1 -Rs)' 102 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. = IGI • 1 + ( H 0 i - R s ^ G '^ H ^ - R s) (D.20) where H = (srs)1/2H. The m inim um of jF ^J can be derived by using provided Q is independent of (f>x. Then we have (D.21) 4>\.mm = (HrG_lH)_1HrG_lRs. Hence. I * 0 , IG I • = IG I • = |G| * where |G| IG I 1 + [^ (H ^ -'H j-^ ^ G -'R s-R sJ G 1 H(HrG_lH)_1HrG_1Rs - Rs]} 1 + ^G~1 2H(HrG~1H)_ 1 H r G ~1 2G _ 1 2 Rs - G - 1/2Rs]' • [G-1/2H(HrG -1 H )-1HrG -1/2G -1/2Rs - G _l/-Rs]} 1 + [H(HrH )-1HrG"1/2Rs - G_1/2Rs]r • [H(HrH )-‘HrG -1/2Rs - G ~ 1 /2Rs] } H = G_1/2H 1 + [(Pft - I)G-1/2Rs]r [(Pft - I)G~l/2Rs] | 1 + [G_1/2Rs]r Pfl [G " 1 /2Rs] J 103 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission. where Pfl is an orthogonal projection m atrix onto the subspace of H. and is an orthogonal projection m atrix onto the complement subspace of H. Xow. the G LR T for the ROI data m atrix can be derived easily. Based on the assumptions made for the ROI data model in this Appendix, the detection problem becomes a test of the distributions: Ho Y: .Y S0ol r. A0 3 M Y: .Y H^lT. Ai ® M' which is equivalent to Hq : R = Y M '1'2 : .V S0ol rM -1/2.Ao ®I H 0ll rM _l/2. A[ © I (D.25) H [ : R = Y M 12: .V where 3 denotes the Kronecker product. The likelihood ratio test can be w ritten as I ^ . A n Y ) L(b,.Ai:R) /(Y) = L(0O . Aq:Y) L(b0. A0: R) ;d .2G) Hence, the G LR T can be derived by substituting the M LE o f A, and /(Y ) L(b[. A[: R) P/2 L(b0. A0: R) V IA 11 / Then it is more convenient to replace the G LRT by the (P /2 )-ro o t G LRT (0.27. A il |G| • {l + (G_1/2Rs)rP|-(G_l/2Rs)} |G| • {l + (G_1/2Rs)rP g (G '1 /2Rs)} 1 + (G-1/2Rs)rPg (G-l/2Rs) 1 + (G~1/2Rs)rP^(G~1/2Rs) ’ (D.28) 1 0 -1 Reproduced with permission of the copyright owner. Further reproduction prohibited without permission.

Asset Metadata

Creator Huang, Chung-Chieh (author)

Core Title Computer-aided lesion detection in positron emission tomography: A signal subspace fitting approach

Contributor Digitized by ProQuest (provenance)

School Graduate School

Degree Doctor of Philosophy

Degree Program Electrical Engineering

Publisher University of Southern California (original), University of Southern California. Libraries (digital)

Tag engineering, electronics and electrical,health sciences, oncology,health sciences, radiology,OAI-PMH Harvest

Language English

Advisor Leahy, Richard (committee chair), Conti, Peter (committee member), Ortega, Antonio (committee member)

Permanent Link (DOI) https://doi.org/10.25549/usctheses-c16-195082

Unique identifier UC11334888

Identifier 3065797.pdf (filename),usctheses-c16-195082 (legacy record id)

Legacy Identifier 3065797.pdf

Dmrecord 195082

Document Type Dissertation

Rights Huang, Chung-Chieh

Type texts

Source University of Southern California (contributing entity), University of Southern California Dissertations and Theses (collection)

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the au...

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus, Los Angeles, California 90089, USA