# Computational aspect of tomographic microwave imaging for biomedical applications

код для вставкиСкачатьCOMPUTATIONAL ASPECT OF TOMOGRAPHIC MICROWAVE IMAGING FOR BIOMEDICAL APPLICATIONS A Thesis Submitted to the Faculty in partial fulfillment of the requirements for the degree of Doctor of Philosophy by AMIR H. GOLNABI Thayer School of Engineering Dartmouth College Hanover, New Hampshire MAY 2012 Examining Committee: Chairman _______________________ Paul M. Meaney, Ph.D. Member ________________________ Keith D. Paulsen, Ph.D. Member ________________________ Andrea Borsic, Ph.D. Member ________________________ Elise C. Fear, Ph.D. _____________________ Brian W. Pogue Dean of Graduate Studies Copyright 2012 Amir H. Golnabi UMI Number: 3544495 All rights reserved INFORMATION TO ALL USERS The quality of this reproduction is dependent upon the quality of the copy submitted. In the unlikely event that the author did not send a complete manuscript and there are missing pages, these will be noted. Also, if material had to be removed, a note will indicate the deletion. UMI 3544495 Published by ProQuest LLC (2012). Copyright in the Dissertation held by the Author. Microform Edition © ProQuest LLC. All rights reserved. This work is protected against unauthorized copying under Title 17, United States Code ProQuest LLC. 789 East Eisenhower Parkway P.O. Box 1346 Ann Arbor, MI 48106 - 1346 Abstract Microwave imaging (MI) is based on recovering dielectric properties (permittivity and conductivity) of materials. Over the last two decades, MI has attracted increasing interests in biomedical applications, in particular, for breast cancer screening and therapy monitoring, mainly due to the significant dielectric property contrast between normal and malignant breast tissue. Compared to other conventional imaging modalities, such as X–ray mammography, the non–ionizing and non–compressive nature of MI makes it very attractive from a patient’s perspective. In addition, MI has a relatively high sensitivity to detect small tumors, and potentially high specificity to determine whether a suspicious area is malignant or benign at a significantly lower cost level compared to methods such as magnetic resonance imaging (MRI) and nuclear medicine. Several years ago, a clinical microwave imaging system was developed at Dartmouth College. During microwave data acquisition, electromagnetic fields propagate through and scatter from the tissue in a three–dimensional (3D) fashion. However, in order to reduce the computational complexity and to speed up the image reconstruction process, it is often assumed that the behavior of electromagnetic waves in 3D space can be represented as a simplified 2D model. While it benefits from less intensive computational demands, such assumption can lead into increased level of artifacts in the recovered dielectric properties. Therefore, in order to improve the accuracy and quality of reconstructed images, we have developed a viable, fast, and user-friendly 3D microwave image reconstruction. We have previously demonstrated that we can recover clinically useful microwave tomographic images for breast cancer detection without using any prior spatial information about the tissue. However, in order to enhance the quality of our previous results and recover more ii accurate dielectric property distributions, we have implemented a new approach to combine the functional information of the MI with the high spatial resolution of other imaging modalities such as MRI. This approach is based on developing a new image reconstruction strategy that exploits structural information about the object being imaged through some constraints in the microwave property contrast. iii Acknowledgements This dissertation would not have been possible without the guidance and the help of several individuals who in one way or another contributed and extended their valuable assistance in the preparation and completion of this study. First and foremost, I would like to show my sincerest gratitude to my advisor, Prof. Paul Meaney, who has supported me throughout my doctoral studies with his patience and knowledge whilst allowing me the room to grow and expand my professional and personal experiences. I should also express my gratitude to my other advisor, Prof. Keith D. Paulsen, who has constantly supported me with encouragement and enthusiasm. Despite his extremely busy schedule, Keith has always made himself available to me for any questions or concerns. In addition, his valuable comments and questions have helped me become a better scientist who is more competent in the field. In addition to my advisors, I would like to thank the rest of my thesis committee, Prof. Andrea Borsic and Prof. Elise Fear for their insightful comments and hard questions. I would also like to thank Prof. Brian Pogue, Prof. Dan Lynch, Prof. George Cybenko, Prof. Alex Hartov, Dr. Ron Lasky, and Dean Joe Helble for their advice and support during my time at Thayer. I am grateful to Amin Plaisted for his invaluable help with this dissertation, and to Dr. Hamid Ghadyani for providing me with the 2D and 3D mesh generating codes. His help has been a key to completing some of the most important parts of my research. Although Margaret Fanning and Sherri Geimer no longer work in the MIS research group, they both deserve my sincere thanks. Margaret and Sherri have been excellent role models for me as very organized and detailed oriented researchers. It was a pleasure working with them. iv I am indebted to many of my current and former colleagues who have supported me during the last five years at Dartmouth: Matt Pallone for his continuous positive attitudes, Neil Epstein and Alex Bijamov for useful discussions, and Tian Zhuo for his xcritical comments and inputs. I would like to thank Dr. Tomasz Grzegorczyk, Preston Manwaring, Dr. Qianqian Fang, and Dr. Elia Attardo for their assistance and thoughtful suggestions. Xiaoyao, Parisa, Lamia, AJ, Ashley, Alex, and Lyubomir, I am grateful for having had such a wonderful friends at Thayer with whom to share the good times working, studying, playing, or just taking a break. Also my special thanks goes to my dearest friend Niusha who has been there for me whenever I needed her. Finally, I give my deepest gratitude to my parents and my wife Marzie for supporting me unconditionally and for believing in my capabilities. Marzie has been very understanding and patient with me during my very long hours in the office including nights and weekends. Without her encouragement, I would not be where I am today without her heartening. v Dedications I dedicate this work to the loves of my life, my wife Marzie, my daughter Niusha, my mother Shahnaz, my father Jalil, my sister Elham, and my brother Mahdiar. vi Table of Contents 1. Introduction ..................................................................................................................................1 1.1. Motivation ............................................................................................................................1 1.2. Overview of Different Microwave Imaging Methods .........................................................4 1.3. Tomographic Microwave Imaging .......................................................................................8 1.4. Project Aims .......................................................................................................................10 1.5. Proposed Work ...................................................................................................................12 2. 3D Microwave Image Reconstruction .......................................................................................15 2.1. 2D vs. 3D Modeling of the Electromagnetic Field ............................................................15 2.2. Challenges With Respect to 2D Reconstruction Algorithm ..............................................19 2.3. Proposed method ................................................................................................................21 2.3.1. User Interface ...........................................................................................................21 2.3.2. 3D Data Acquisition ................................................................................................22 2.3.3. Input–File Format ....................................................................................................23 2.3.4. Optimizing the Reconstruction Algorithm...............................................................29 2.4. Results – Simulations .........................................................................................................30 2.4.1. FDTD Grid Density .................................................................................................31 2.4.2. Number of Reconstruction Nodes ............................................................................35 2.4.3. 3D images from in-plane data vs. stack of 2D images ............................................40 2.4.4. Data Selection for Reconstruction ...........................................................................43 2.4.5. Number of Iterations ................................................................................................45 2.5. Results – Phantom Experiments.........................................................................................49 2.5.1. Two Cylindrical Inclusions ......................................................................................49 vii 2.5.2. Breast-shaped Phantom Experiment with an Spherical Inclusion ...........................54 2.6. Results – Initial Clinical Data ............................................................................................63 2.6.1. Breast Imaging .........................................................................................................63 2.6.2. Bone Imaging ...........................................................................................................72 3. Incorporation of a Priori Structural Information into the Microwave Image Reconstruction Algorithm .........................................................................................................77 3.1. Introduction ........................................................................................................................77 3.2. A Mathematical Framework of Microwave Image Reconstruction Algorithm .................79 3.2.1. Overview ..................................................................................................................79 3.2.2. Soft prior Regularization..........................................................................................82 3.2.3. Hard priors – Parameter Reduction..........................................................................85 3.2.4. Error Analysis ..........................................................................................................86 3.3. Integration of Microwave Tomographic Imaging System with MRI ................................88 3.3.1. Motivation ................................................................................................................88 3.3.2. Integration of MI System into MRI .........................................................................89 3.3.3. MR Image Artifact Reduction Strategies .................................................................91 3.3.4. Multipath Signal Assessment ...................................................................................93 3.3.5. MR RF Gradient Signal Disruption of Microwave Signals .....................................95 3.4. Results – 2D .......................................................................................................................99 3.4.1. Simulation Experiments ...........................................................................................99 3.4.2. Phantom Experiments ............................................................................................131 3.4.3. Bone Imaging: Monitoring Changes in the Bone Dielectric Property Distributions............................................................................................................178 viii 3.4.4. Clinical Breast Microwave Imaging in MR ...........................................................186 3.5. Results – 3D .....................................................................................................................191 3.5.1. Simulation Experiments .........................................................................................191 3.5.2. Phantom Experiments ............................................................................................220 4. 3D Microwave Image Reconstruction GUI .............................................................................235 5. Summary and Conclusion ........................................................................................................244 ix List of Figures Figure 1.1. Tree diagram of the current MI methods for breast cancer detection ...........................4 Figure 1.2. MIST system: (a) antenna configuration on two independently-moving plates A (pink) and B (blue), (b) imaging tank, (c) exam platform, and (d) cabling and fluid reservoir underneath the table .....................................................................................9 Figure 2.1. Simulation experiment setup with a cylindrical inclusion or 1.0 cm radius ................17 Figure 2.2. Calibrated amplitude/power (top) and phase (bottom) values of the synthesized data using the 2D and 3D models as a function of the relative receivers for a single transmitter (tx = 1) ...........................................................................18 Figure 2.3. RMS difference values of the amplitude and phase ....................................................19 Figure 2.4. Schematic antenna configuration: Two interleaved arrays of 8 antennas with each sub-array (A and B) being able to move independently from the other ....................23 Figure 2.5. Schematic summary of steps to create 3D image reconstruction input files from the MIST system output data ....................................................................................24 Figure 2.6. Input File–1 .................................................................................................................25 Figure 2.7. Input File–2 .................................................................................................................25 Figure 2.8. Input File–3 .................................................................................................................26 Figure 2.9. Input File–3’ ................................................................................................................27 Figure 2.10. Input data selection GUI in MATLAB ......................................................................28 Figure 2.11. Schematic summary of reconstruction algorithm ......................................................29 Figure 2.12. Schematic of the imaging domains evaluated: The background diameter was 14 cm (antennas are positioned on a 15.2 cm diameter). The spherical inclusion of 1.5 cm radius was centered at (3, 0, 0 cm). ........................................................................30 x Figure 2.13. 1300 MHz contour slice images extracted from the reconstructed permittivity (top) and conductivity (bottom) profiles for the simulation experiment using the finite difference grid spacing cases (a) to (f) in Table 2.1. ................................................34 Figure 2.14. Reconstructed permittivity (left) and conductivity (right) profiles along with the true property values, as a function of the number of cells per wavelength (D) in the forward solution .......................................................................................................35 Figure 2.15. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) profiles for the simulation experiment on (a) 600, (b) 1587, (c) 2611, (d) 3767, (e) 4608, (f) 6095, (g) 7669, and (h) 9321 node meshes. ........................................................39 Figure 2.16. Reconstructed permittivity (left) and conductivity (right) profiles along with the true property values, as a function of the number of reconstruction nodes .................39 Figure 2.17. Total reconstruction times (min) per iteration as a function of the number of reconstruction nodes ..........................................................................................................40 Figure 2.18. (a) 2D and (b) 3D reconstructed permittivity (top row) and conductivity (bottom row) slices (in z-direction) of the simulation experiment with -100 dBm added noise at 1300 MHz, using only in-plane data ..........................................................42 Figure 2.19. 3D reconstructed permittivity (top row) and conductivity (bottom row) slices (in z-direction) of the simulation experiment with -100 dBm added noise at 1300 MHz, using (a) the full-data (all 7 planes of data, including both in-planes and cross-planes), and (b) multi-plane data with two consecutive planes ...............................44 Figure 2.20. 3D Relative errors as a function of the number of iterations from a simulation experiment........................................................................................................46 xi Figure 2.21. Relative permittivity (top) and conductivity (bottom) difference values for iterations (a) 40, (b) 60, (c) 80, and (d) 100 (i.e. for D20,40, D20,60, D20,80, and D20,100, respectively), along with the mean of the difference values ± 2 standard deviations (i.e. 95% of all difference data) ........................................................................48 Figure 2.22. Setup for the phantom experiment with two cylindrical inclusions: (a) The square-based inclusion tilted toward the right side of the imaging domain, (b) Both inclusions immersed in a coupling medium and surrounded by the antenna array ...................................................................................................................................49 Figure 2.23. 1100 MHz reconstructed permittivity values of the 3D phantom experiment with two cylindrical inclusions. The iso-surface thresholds are εr,CircInc = 10.0 and εr,SqInc = 33.0. .....................................................................................................................51 Figure 2.24. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) slices of the phantom experiment with two cylindrical inclusions using 3D reconstruction algorithm with (a) full-data set, (b) 5 multi-plane data sets with two consecutive planes, (c) 5 in-plane data sets, (d) 5 individual in-plane data sets, and (e) using 2D reconstruction algorithm with 5 individual in-plane data sets......................................................................................................................................54 Figure 2.25. Breast-shaped phantom experiment setup: (a) The spherical saline gel target inclusion suspended in the plastic breast model, (b) Rapid–prototyped plastic breast model submerged in the imaging tank ....................................................................55 Figure 2.26. 3D reconstructed (a) permittivity and (b) conductivity images, using the multi-plane data set at 1300 MHz. The iso-surface thresholds are εr,BR = 15.0, σBR = 1.0 S/m, and εr,Inc = 65.0 and σInc = 2.3 S/m. ..................................................................56 xii Figure 2.27. Extracted slices from the 3D reconstructed images in Figure 2.26: Permittivity (top row) and conductivity (bottom row) ......................................................58 Figure 2.28. 2D reconstructed dielectric properties: Permittivity (top row) and conductivity (bottom row) .................................................................................................59 Figure 2.29. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) slices of the breast-shaped phantom experiment at z = 3 cm using the 3D reconstruction algorithm with (a) full-data set, (b) 9 multi-plane data sets with two consecutive planes, (c) 5 in-plane data sets, and (d) using the 2D reconstruction algorithm with 9 individual in-plane data sets ...................................................................................61 Figure 2.30. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) slices of the breast-shaped phantom experiment at z = 0 cm using the 3D reconstruction algorithm with (a) full-data set, (b) 9 multi-plane data sets with two consecutive planes, (c) 5 in-plane data sets, and (d) using the 2D reconstruction algorithm with 9 individual in-plane data sets ...................................................................................61 Figure 2.31. 3D microwave imaging of a normal subject in the clinical MIST system ................64 Figure 2.32. Extracted slices from the 1100 MHz 3D reconstructed images using the fulldata set: Permittivity (top row) and conductivity (bottom row) ........................................65 Figure 2.33. 1100 MHz 2D reconstructed dielectric properties: Permittivity (top row) and conductivity (bottom row) .................................................................................................65 Figure 2.34. 1100 MHz reconstructed permittivity (top) and conductivity (bottom) slices of the breast-shaped phantom experiment at z = 1.5 cm using the 3D reconstruction algorithm with (a) the full-data (FD) set, (b) 6 multi-plane data sets with four consecutive planes (4CS) (two planes above and two planes below), (c) xiii 6 multi-plane data sets with two consecutive planes (2CS) (one plane above and one plane below), (d) 6 in-plane (IP) data sets, and (e) using the 2D reconstruction algorithm with 6 individual in-plane data sets. ..................................................................67 Figure 2.35. Relative permittivity (top) and conductivity (bottom) difference images with respect to the recovered dielectric property values using the full-data (FD) set: (a) full-data – multi-plane data sets with four consecutive planes (i.e. FD – 4CS) (b) full-data – multi-plane data sets with two consecutive planes (i.e. FD – 2CS), and (c) full-data – in-plane data sets (i.e. FD – IP) ..................................................................69 Figure 2.36. Relative permittivity (top) and conductivity (bottom) difference values in Figure 2.35 as a function of the reconstruction nodes, along with the mean of the difference values ± 2 standard deviations (i.e. 95% of all difference data): (a) fulldata – multi-plane data sets with four consecutive planes (i.e. FD – 4CS) (b) fulldata – multi-plane data sets with two consecutive planes (i.e. FD – 2CS), and (c) full-data – in-plane data sets (i.e. FD – IP)........................................................................71 Figure 2.37. Heel imaging in the microwave imaging system: (a) The heel extended through an aperture in the mounting plate, (b) The heel surrounded by the antenna array ...................................................................................................................................72 Figure 2.38. 2D and 3D reconstruction meshes: (a) 6.9 cm radius circular mesh: 559 nodes and 1044 triangular elements, (b) 6.9 cm radius, 5 cm height cylindrical mesh: 3132 nodes and 12707 tetrahedral elements............................................................73 Figure 2.39. Several views of the 3D reconstructed permittivity image of the heel at 1300 MHz, with an iso-surface threshold of εr = 26 ...................................................................74 xiv Figure 2.40. 1300 MHz, (a) 2D and (b) slices of the 3D reconstructed permittivity (top) and conductivity (bottom) images at 4 planes of data separated by 1 cm, starting from the ankle (P1) and moving down towards the heel apex (P4). ..................................75 Figure 3.1. Soft prior regularization matrix L with two distinct regions .......................................83 Figure 3.2. MI-MR imaging: (a) the portable data acquisition system in the control room with the auxiliary cables attached to the bulkhead connectors in the wall, (b) the imaging tank outside of the MR bore showing the semi-rigid cables extending from the tank ......................................................................................................................90 Figure 3.3. MI tall tank with a large cylinder breast phantom (green) and smaller inclusions surrounded by the monopole antenna array: (a) side view and (d) topdown view ..........................................................................................................................92 Figure 3.4. T2–weighted MR images of the phantom experiment: (a) at the center plane of the active part of the antennas and (b) at a plane 3 cm below the interface between the semi-rigid coaxial line and the antenna. For (a) and (b) the phantom included a 5.4 cm radius outer cylinder and the 2.1 and 1.0 cm radius inclusions, and the coaxial line center conductors were silver-plated copper. (c) Image at the center plane of the active part of the antennas for the 5.4 cm radius outer cylinder and just the 1.0 cm radius inclusion. In this case the coaxial cable center conductor was steel. ...........................................................................................................93 Figure 3.5. Monopole antenna array: (a) in the MI shortened tank with a breast phantom (green), and straight 5.0 cm long feedlines (short), (b) with 9.0 cm long serpentine feedlines (tank wall removed) ............................................................................................95 xv Figure 3.6. Plots of the 1300 MHz amplitudes as a function of receiver number for all 16 transmitting antennas when the imaging tank is operating in the MR bore with MR data being acquired simultaneously for the cases where (a) no low pass filters, and (b) two low pass filters installed, respectively. ................................................97 Figure 3.7. Schematic of the imaging domains evaluated. The background diameter was 14 cm (antennas are positioned on a 15.2 cm diameter). The circular inclusion’s center with radius r = 1.40 cm was located at (0, –3 cm). ...............................................100 Figure 3.8. Reconstruction parameter meshes (a) uniformly distributed 473 node mesh, (b) uniformly distributed 961 node mesh, and (c) 915 node mesh with preferential node deployment in the inclusion in Figure 3.7 ...............................................................100 Figure 3.9. Simulated 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images (a) without priors on the 473 node mesh, (b) without priors on the 961 node mesh, (c) without priors on the 915 node mesh with preferential node deployment in the inclusion, and (d) with the soft prior regularization on the 915 node mesh. ................................................................................................................101 Figure 3.10. Comparison of the 1300 MHz reconstructed permittivity (top) and conductivity (bottom) values using the no priors and soft prior regularization with different levels of added noise: (a) –110, (b) –100, (c) –90, (d) –80 dBm ......................103 Figure 3.11. Comparison of 1300 MHz reconstructed permittivity (top row) and conductivity (bottom row) profiles for the no priors and the soft prior regularization: (a) no permittivity contrast (left column), (b) no conductivity contrast (right column). ....................................................................................................105 xvi Figure 3.12. (a) Schematic of the imaging domains evaluated. The arbitrarily shaped inclusion was located in the upper part of the imaging domain. The reconstructed values of the simulation experiment with arbitrarily shaped inclusion were extracted at 30 points evenly distributed along the line x = –2 cm. (b) Customized soft prior reconstruction mesh comprised of 1725 nodes and 3248 triangular elements. ..........................................................................................................................106 Figure 3.13. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images from a phantom experiment with arbitrarily-shaped target inclusion for (a) no priors on the 473 node mesh, and (b) soft priors on the customized 1725 node mesh in Figure 3.12(b) .....................................................................................................107 Figure 3.14. Comparison of the 1300 MHz reconstructed (a) permittivity and (b) conductivity profiles (along the line x = –2 cm) in Figure 3.13 phantom experiment with arbitrarily-shaped target inclusion for the no priors and soft prior regularization ...................................................................................................................108 Figure 3.15. Customized soft prior reconstruction meshes used in 1 and 3–inclusion simulation experiments: (a) 402 node mesh with one inclusion and two distinct regions (RBK and RI1), (b) 1203 node mesh with one inclusion and two distinct regions (RBK and RI1), (c) 407 node mesh with three inclusions and two distinct regions (RBK and RI1), (d) 1203 node mesh with three inclusions and two distinct regions (RBK and RI1), (e) 407 node mesh with three inclusions and four distinct regions (RBK, RI1, RI2, and RI3), and (f) 1203 node mesh with three inclusions and four distinct regions (RBK, RInc1, RInc2, and RInc3) .............................................................109 xvii Figure 3.16. Soft prior reconstructed (a) permittivity and (b) conductivity values using two meshes with different number of nodes, along with the exact properties of background and inclusion at 1300 MHz...........................................................................110 Figure 3.17. Maximum soft prior recovered (a) permittivity and (b) conductivity values along with the exact properties of the background and target inclusions at 1300 MHz using the customized reconstruction meshes in Figures 3.15(c) – (f) .....................112 Figure 3.18. Maximum soft prior recovered (a) permittivity and (b) conductivity values along with the exact properties of the background and target inclusions at 1300 MHz using the customized reconstruction meshes in Figures 3.15(e) and (f) .................115 Figure 3.19. Comparison between the soft prior reconstructed (a) permittivity and (b) conductivity values of I1 with and without additional targets (three inclusions and one inclusion cases, respectively) ....................................................................................117 Figure 3.20. Customized reconstruction mesh comprised of 1196 nodes and 2215 triangular elements, accounting for two inclusion regions: a false circular region of interest (in red) with radius of 1.4 cm centered at (0, 3 cm), and the actual circular target zone (in green) with radius of 1.4 cm centered at (0, -3 cm) ....................119 Figure 3.21. Transect plots of the 1300 MHz reconstructed soft prior permittivity (top) and conductivity (bottom) profiles along the y-axis for three simulation cases with a false region of interest. Dielectric properties of the actual inclusion were set to those of the (a) 40:60, (b) 55:45 and (c) 70:30 mixtures of glycerin and water, reported in Table 3.3 ........................................................................................................121 Figure 3.22. Relative errors between the maximum soft prior recovered properties and the actual values in Figure 3.21 .............................................................................................121 xviii Figure 3.23. (a) Patient’s breast MR image (b) The customized soft prior mesh composed of 1465 nodes and 2706 triangular elements containing only internal structure of the breast – adipose (orange) and fibroglandular (blue) tissue, in (a). Reconstructed dielectric property values were extracted along the dark blue line across the breast and fibroglandular tissues. ....................................................................123 Figure 3.24. Simulated 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images with -100 dBm added noise, using the soft priors (left) and no priors (right) for cases (a) – (d) in Table 3.4....................................................................125 Figure 3.25. Comparison of the 1300 MHz recovered permittivity (top) and conductivity (bottom) values along an arbitrary line (illustrated in Figure 3.23(b)) in Figure 2.24: (a) – (d) correspond to different target values in case (a) – (d) in Table 3.4 ..........127 Figure 3.26. Comparison of the 1300 MHz soft prior reconstructed permittivity (top) and conductivity (bottom) values with the standard (Stnd w) and new (w*0.1) scaling factor along the same arbitrary line shown in Figure 3.23(b): (a) – (d) correspond to different target values in case (a) – case (d) in Table 3.4. ...........................................130 Figure 3.27. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images from a phantom experiment with a single circular inclusion for (a) no priors on the 473 node mesh and (b) soft priors on the 915 node mesh in Figure 3.8 ....................132 Figure 3.28. Comparison of the 1300 MHz reconstructed permittivity (left) and conductivity (right) profiles extracted along the y-axis in Figure 3.26 phantom experiment using the no priors and soft prior regularizations. ........................................133 Figure 3.29. Comparison of the reconstructed permittivity (left), and conductivity (right) values from a phantom experiment using different spatial prior coefficient λ xix values at (a) 900, (b) 1100, (c) 1300, (d) 1500, (e) 1700, (f) 1900, (g) 2100, (h) 2300, and (i) 2500 MHz ...................................................................................................138 Figure 3.30. Weighted (a) εr and (b) σ errors for a phantom experiment over a range of frequencies from 900 to 2500 MHz using six different soft prior coefficients: λ = 0.01, 0.1, 1, 10, 100, and 1000 .........................................................................................139 Figure 3.31. Customized soft prior reconstruction meshes used for the phantom experiments with a circular inclusion of radius: (a) 0.65 cm – 1329 nodes and 2530 triangular elements, (b) 1.4 cm – 1037 nodes and 1943 triangular elements, and (c) 2.1 cm – 813 nodes and 1498 triangular elements ...............................................141 Figure 3.32. Comparison of the reconstructed permittivity (left), and conductivity (right) values from three phantom experiments with different sized inclusions (of radius 0.65, 1.4, and 2.1 cm) using the no prior (LM) and soft prior (SP) regularizations at (a) 900, (b) 1100, (c) 1300, (d) 1500, and (e) 1700 MHz ............................................143 Figure 3.33. Customized soft prior reconstruction meshes used for the phantom experiments with: (a) circular inclusion – 1329 nodes and 2530 triangular elements, (b) square-shaped inclusion – 1096 nodes and 2096 triangular elements, and (c) diamond-shaped inclusion – 1040 nodes and 1984 triangular elements .............144 Figure 3.34. Comparison of the reconstructed permittivity (left), and conductivity (right) values from a phantom experiment with three different shaped but similar sized inclusions (1.3 cm diameter circular, 1.3 cm side length square, and 1.3 cm side length diamond-shaped cylinders) using the no prior (LM) and soft prior (SP) regularizations at (a) 900, (b) 1100, (c) 1300, (d) 1500, and (e) 1700 MHz ....................147 xx Figure 3.35. Comparison of the soft prior reconstructed permittivity (left) and conductivity (right) values of five different contrast-level target-inclusions (50:50, 60:40, 70:30, 90:10 and 100:0 mixture of glycerin:water) at (a) 900, (b) 1100, (c) 1300, (d) 1500, and (d) 1700 MHz, along with linear least-square fits ...........................150 Figure 3.36. Comparison of the soft prior reconstructed permittivity (left), and conductivity (right) values of a target inclusion using several soft prior customized meshes with different node densities: case (a) number of nodes in the inclusion region (NInc) changes, case (b) number of nodes in the background region (NBK) changes, and case (c) number of nodes in both regions (NTotal) changes proportionally .....................................................................................................154 Figure 3.37. 1300 MHz soft prior reconstructed permittivity (top) and conductivity (bottom) images of a phantom experiment using the reconstruction mesh with a false region of interest in Figure 3.20 ..............................................................................156 Figure 3.38. Comparison of the 1300 MHz soft prior reconstructed permittivity (left), and conductivity (right) profiles from a phantom experiment using the reconstruction mesh with a false region of interest in Figure 3.20 ..........................................................157 Figure 3.39. 1100 MHz reconstructed permittivity (top) and conductivity (bottom) images of a phantom experiment using different prior inclusion-sized meshes: (a) 0.8 cm, (b) 1.0 cm, (c) 1.2 cm, (d) 1.4 cm (the true size), (e) 1.6 cm, (f) 1.8 cm, and (g) 2.0 cm were considered as the radius of the inclusion in the reconstruction mesh. The actual size (radius of 1.4 cm) and location (0, – 3 cm) of the inclusion is outlined by a circle on the reconstructed images. ..........................................................................159 xxi Figure 3.40. Comparison of the extracted 1100 MHz reconstructed permittivity (left), and conductivity (right) values in Figure 3.39 along with the exact property distributions for a phantom experiment using different inclusion–sized soft prior meshes ..............................................................................................................................160 Figure 3.41. Comparison of the extracted (a) 1300, (b) 1500, and (c) 1700 MHz reconstructed permittivity (left), and conductivity (right) profiles along with the exact property distributions for a phantom experiment using different inclusion– sized meshes.....................................................................................................................162 Figure 3.42. 1100 MHz reconstructed permittivity (top) and conductivity (bottom) images of a phantom experiment using different prior inclusion–located meshes: x = 0 and y–coordinates of (a) – 3.9, (b) – 3.0 (the true location), (c) – 2.1, (d) – 1.2, (e) – 0.3, (f) 0.6, and (g) 1.5 cm were considered as the center of the inclusion in the reconstruction mesh. The actual size (radius of 1.4 cm) and location (0, – 3 cm) of the inclusion is outlined by a circle on the reconstructed images. ...................................164 Figure 3.43. Comparison of the extracted (a) 1100, (b) 1300, (c) 1500, and (d) 1700 MHz reconstructed permittivity (left), and conductivity (right) values along with the exact property distributions for a phantom experiment using different inclusion– located meshes .................................................................................................................166 Figure 3.44. Real breast–shaped model: (a) breast mold, (b) rapid–prototyped breast model................................................................................................................................168 Figure 3.45. 3D optical scanning of the model ............................................................................169 Figure 3.46. Exported 3D model (a) 3D surface mesh, (b) Formation of a solid object and extraction of 2D slices in the SolidWorks software package ..........................................170 xxii Figure 3.47. 2D customized soft prior meshes used in a breast–shaped phantom experiment with circular target inclusions of (a) 0.65 cm and (b) 1.0 cm radius, as well as a polyline along which the actual and reconstructed dielectric property distributions were extracted. The mesh in (a) is comprised of 1778 nodes and 3468 triangular elements, whereas the one in (b) contains 1788 nodes and 3488 triangular elements. ..........................................................................................................171 Figure 3.48. Comparison of the 1300 MHz reconstructed permittivity (left) and conductivity (right) profiles along an arbitrary polyline for the breast–shaped phantom experiments with (a) small inclusion of 0.65 cm radius and (b) medium– sized inclusion of 1 cm radius, using the no priors and soft prior regularization ............172 Figure 3.49. T2–weighted MR images of the phantom experiment at the center plane of the active part of the antennas ..........................................................................................174 Figure 3.50. 2D reconstruction meshes used for the phantom experiment in MR: (a) No priors, with 559 uniformly distributed nodes and 1044 triangular elements, and (b) soft priors, with 1034 nodes and 1821 triangular elements .............................................175 Figure 3.51. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images of the multi–region phantom case using the (a) no priors (circles indicate the exact object locations) and (b) soft prior regularization............................................................176 Figure 3.52. Experiment flow diagram for evaluating the correlation between bone mineral and dielectric properties of porcine trabecular bone specimens .........................179 Figure 3.53. Microwave data acquisition of the bone sample .....................................................180 Figure 3.54. Soft prior microwave image reconstruction mesh comprised of 1107 nodes and 2092 triangular elements ...........................................................................................181 xxiii Figure 3.55. The 1100 MHz reconstructed permittivity (top) and conductivity (bottom) images of the (a) first, (b) second, and (c) final (5th) microwave scans...........................181 Figure 3.56. Relationship between the dielectric properties of saline–saturated bones at 1100 MHz and BVF of dry bone specimens ....................................................................182 Figure 3.57. X–ray CT images of the patient’s heels for the (a) normal weight–bearing and (b) injured leg (White zones within the calcaneous bones are the radiologists’ labels.) ..............................................................................................................................184 Figure 3.58. 2D customized reconstruction meshes of the (a) normal weight–bearing and (b) affected leg comprised of 1196 nodes and 2264 triangular elements, and 815 nodes and 1519 triangular elements, respectively ...........................................................185 Figure 3.59. 1300 MHz soft prior permittivity (top) and conductivity (bottom) images for the heel of a (a) normal and an (b) injured leg utilizing the CT images as prior information during the reconstruction process ................................................................186 Figure 3.60. Photograph of a patient in the combined MI-MRI system (a) before entering the MR bore and (b) inside the bore ................................................................................187 Figure 3.61 (a) MR T1-weighted image of the patient’s left breast at the plane of the active part of the monopole antennas, and (b) the corresponding segmented soft prior mesh ........................................................................................................................188 Figure 3.62. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images of the patient’s left breast utilizing the (a) no priors and (b) soft prior regularizations, respectively ............................................................................................189 xxiv Figure 3.63. 1300 MHz 3D reconstructed permittivity (top) and conductivity (bottom) images of the simulation experiment with – 100 dBm added noise using (a) the soft prior regularization and (b) no prior spatial information ..........................................192 Figure 3.64. Extracted 3D recovered permittivity (left) and conductivity (right) values along with the exact property values, as a function of the number of nodes in the soft prior reconstruction mesh .........................................................................................195 Figure 3.65. 1300 MHz 3D reconstructed permittivity (top) and conductivity (bottom) profiles for the simulation experiment using different soft prior coefficients λ = (a) 0.01, (b) 0.1, (c) 1.0, (d) 10, and (e) 100..........................................................................198 Figure 3.66. Extracted 3D recovered permittivity (left) and conductivity (right) values along with the exact property distributions, as a function of the base – 10 logarithm of the soft prior coefficient λ ...........................................................................198 Figure 3.67. Relative RMS errors of the recovered dielectric property distributions as a function of the base – 10 logarithm of the soft prior regularization coefficient λ ...........199 Figure 3.68. Extracted soft prior 3D recovered permittivity (left) and conductivity (right) values along with the exact property distributions, using different inclusion–sized meshes as prior information .............................................................................................201 Figure 3.69. Extracted soft prior 3D recovered permittivity (left) and conductivity (right) values along with the exact property distributions, using different inclusion– located meshes as prior information ................................................................................204 Figure 3.70. 1300 MHz 3D reconstructed images using the corresponding soft prior meshes with a false region of interest in Table 3.18 ........................................................208 xxv Figure 3.71. 3D soft prior recovered (a) permittivity and (b) conductivity values extracted at the center of each region, along with the exact properties, using customized reconstruction meshes in Table 3.20 ................................................................................211 Figure 3.72. Soft prior 3D recovered (a) permittivity and (b) conductivity values extracted at the center of each region, along with the exact properties, using 4region customized reconstruction meshes in Table 3.20 .................................................214 Figure 3.73. Comparison between the soft prior reconstructed (a) permittivity and (b) conductivity values of I1 with and without additional targets (i.e. multi-inclusion and single inclusion cases, respectively) .........................................................................216 Figure 3.74. Recovered permittivity (top) and conductivity values at the center of the (a) cylindrical and (b) spherical inclusions as a function of the inclusion size, using different reconstruction approaches (2D and 3D, with no priors, soft priors, and hard priors) .......................................................................................................................219 Figure 3.75. Complex-shaped phantom experiment with a cylindrical and a cup-shaped gelatin inclusion: (a) cylindrical inclusion (CylInc) inserted into the cup-shaped gelatin inclusion (GelInc), (b) both inclusions submerged in the imaging tank, and (c) top view of the phantom setup in the microwave imaging tank .................................221 Figure 3.76. MRI Imaging of the phantom experiment with a cylindrical and a cupshaped gelatin inclusions: (a – b) pictures of the phantom placed in an empty microwave imaging tank inside the MR bore, and (c) top view of the schematics of the phantom in the MRI scanner..................................................................................223 xxvi Figure 3.77. (a) One of the corresponding MR images of the inclusions, (b) the corresponding 3D soft prior mesh composed of 3393 nodes and 17027 tetrahedral elements ...........................................................................................................................224 Figure 3.78. A vertical slice through (a) soft prior and (b) no prior reconstructed permittivity (top) and conductivity (bottom) images of the phantom experiment with a cylindrical and a cup-shaped gelatin inclusions at 1300 MHz ..............................225 Figure 3.79. Breast–shaped phantom experiment with two arbitrarily shaped target inclusions: (a) Two arbitrary-shaped gelatin target inclusions suspended in the plastic breast model, (b) rapid–prototyped plastic breast model submerged in the imaging tank, and (c) schematic configuration of the breast phantom in the microwave imaging tank. .................................................................................................227 Figure 3.80. MRI Imaging of the breast phantom experiment with two inclusions: (a – b) pictures of the phantom placed in an empty microwave imaging tank inside the MR bore, and (c) top view of the phantom configuration in the MRI scanner................229 Figure 3.81. (a) Stack of binary-segmented MRI images of the breast phantom experiment with two target inclusions, (b) the corresponding 3D soft prior mesh composed of 7540 nodes and 34591 tetrahedral elements...............................................230 Figure 3.82. A vertical slice through the (a) soft prior and (b) no prior reconstructed permittivity (top) and conductivity (bottom) images of the breast phantom experiment with two target inclusions at 1100 MHz .......................................................231 Figure 4.1. Input file selection: (a) file format and (b) raw measurement file selection windows ...........................................................................................................................235 xxvii Figure 4.2. Sample antenna configuration plot where 6 planes of data with a slice thickness of 1 cm have been acquired ..............................................................................236 Figure 4.3. Output file format selection window .........................................................................237 Figure 4.4. 3D data selection window .........................................................................................238 Figure 4.5. Example of a multi-plane data selection with 2 consecutive planes .........................239 xxviii List of Tables Table 2.1. Finite difference grid spacing (in mm) used to compute the forward solution for the simulation experiment, along with the corresponding reconstruction times per iteration ........................................................................................................................32 Table 2.2. Parameter reconstruction meshes used to solve the inverse problem for the simulation experiment........................................................................................................36 Table 2.3. Independently measured dielectric properties of the phantom experiment with two cylindrical inclusions ..................................................................................................50 Table 2.4. Independently measured dielectric properties of the background coupling medium (BK), breast model (BR), and target inclusion (Inc) at 1300 MHz.......................55 Table 3.1. Relative (a) permittivity and (b) conductivity errors for the multiple inclusions simulation experiment using different reconstruction meshes with two and four distinct regions .................................................................................................................113 Table 3.2. Relative (a) permittivity and (b) conductivity errors for the multiple inclusion simulation experiment with different property contrasts using the soft prior reconstruction meshes with four distinct regions in Figures 3.15(e) and (f) ...................116 Table 3.3. Dielectric properties of different mixtures of glycerin and water used for the simulation experiments with a false region of interest ....................................................118 Table 3.4. Dielectric property values of the background medium, breast, and fibroglandular region, in four breast simulation experiments, cases (a) – (d), at 1300 MHz .........................................................................................................................124 Table 3.5. Soft prior relative dielectric property errors in the fibroglandular region, breast-shaped simulation experiments with two different scaling factors .......................130 xxix Table 3.6. Independently measured dielectric properties of the background medium (BK) and target inclusion (Inc) over the range of frequencies evaluated .................................134 Table 3.7. Properties of the soft prior customized meshes used for the phantom experiment with a single target inclusion of 1.4 cm radius..............................................152 Table 3.8. Weighted soft prior εr and σ errors for a phantom experiment using a reconstruction mesh with a false region of interest for different soft prior coefficients: λ = 0.01, 0.1, and 1.0 at 1300 MHz .............................................................158 Table 3.9. Characteristics of the soft prior reconstruction meshes used for a phantom experiment with imperfect prior size of the target inclusion ...........................................159 Table 3.10. Characteristics of the soft prior customized meshes used for a phantom experiment with imperfect prior location of the target inclusion ....................................163 Table 3.11. Independently measured dielectric properties of the background medium (BK), breast (BR), and inclusions (Inc) at 1300 MHz ......................................................169 Table 3.12. (a) The actual and averaged recovered property values of the fibroglandular and tumor inclusion regions in the phantom experiment in the MR bore using the no priors, soft priors, and hard priors, (b) The corresponding relative dielectric property errors ..................................................................................................................177 Table 3.13. Pearson coefficient for the recovered dielectric properties at various frequencies .......................................................................................................................183 Table 3.14. The permittivity and conductivity RMS errors associated with the inclusion region in Figure 3.63 ........................................................................................................193 Table 3.15. Customized reconstruction meshes used for the simulation experiment with soft prior regularization....................................................................................................194 xxx Table 3.16. Characteristics of the 3D soft prior customized meshes used for a simulation experiment with imperfect prior size of the target inclusion ...........................................200 Table 3.17. Characteristics of the 3D soft prior customized meshes used for a simulation experiment with imperfect prior location of the target inclusion ....................................203 Table 3.18. Characteristics of the 3D soft prior customized meshes used for a simulation experiment with a false target inclusion (Bk = background, TI = True Inclusion, and FI = False Inclusion) .................................................................................................206 Table 3.19. Recovered dielectric properties (εr and σ) of the true (TI) and false inclusion (FI) extracted at the center of each region .......................................................................209 Table 3.20. Soft prior reconstruction meshes used for the simulation experiment with three spherical target inclusions .......................................................................................210 Table 3.21. Relative (a) permittivity and (b) conductivity errors for the multi-inclusion simulation experiment using different reconstruction meshes with two and four distinct regions. The negative sign indicates that the reconstructed property value is underestimated. (i.e. less than the exact property values) ...........................................212 Table 3.22. Relative (a) permittivity and (b) conductivity errors for the multi-inclusion simulation experiment with different property contrasts using two soft prior reconstruction meshes with 4 distinct regions .................................................................215 Table 3.23. Characteristics of the customized 2D and 3D meshes used for the simulation experiments with different sized and shaped target inclusions ........................................217 Table 3.24. Independently measured dielectric properties of the background medium (Bk), cylindrical inclusion (CylInc), and cup-shaped gelatin inclusion (GelInc) at 1300 MHz .........................................................................................................................222 xxxi Table 3.25. Soft-prior and no prior RRMS errors of the recovered properties in each region of the breast phantom experiment with two target inclusions ..............................226 Table 3.26. Independently measured dielectric properties of the background medium (Bk), breast model (Br), and the target inclusions (Inc1 and Inc2) at 1100 MHz ............228 Table 3.27. Soft prior and no prior RRMS errors of the recovered properties in each region of the breast phantom experiment with two target inclusions ..............................232 Table 3.28. The computed contrast enhancements of the target inclusions with respect to the breast region, when the prior structural information of the phantom is used during the reconstruction .................................................................................................233 Table 4.1. A list of inputs/3D reconstruction parameters along with a brief description of them..................................................................................................................................240 xxxii List of Acronyms 2D Two-dimensional 3D Three-dimensional ACS American Cancer Society Adp Adipose Bk Background Br Breast BVF Bone Volume Fraction CS Consecutive planes CT Computed Tomography DHMC Dartmouth Hitchcock Medical Center EIS Electrical Impedance Spectroscopy Err Error FD full-data FDTD Finite Difference Time Domain Fg Fibroglandular GUI Graphical User Interface Inc Inclusion IP In-plane LM Levenberg Marquardt MI Microwave Imaging MIS Microwave Imaging Spectroscopy MRE Magnetic Resonance Elastography xxxiii MRI Magnetic Resonance Imaging NIR Near Infrared NP Number of Planes RF Radiofrequency Field RX Receiver SA Spectrum Analyzer SP Soft Prior TSAR Tissue Sensing Adaptive Radar TK Tikhonov TM Transverse Magnetic Tu Tumor TX Transmitter UWB Ultrawide Band VHF Very High Frequency xxxiv 1. Introduction 1.1. Motivation Microwave imaging (MI) is based on recovering dielectric properties (permittivity and conductivity) of materials. Over the last two decades MI has attracted increasing interest in biomedical applications. In particular, interest in breast cancer screening and therapy monitoring. Excluding skin cancer, breast cancer is the most commonly diagnosed cancer in the United States, and the second leading cause of cancer death in women, exceeded only by lung cancer [1]. In fact, it was estimated that over 230,000 new cases of invasive breast cancer would be diagnosed among women in the U.S. in 2011 [1]. The American Cancer Society (ACS) estimates that a woman in the United States has a 1 in 8 chance of developing invasive breast cancer during her lifetime [1]. Despite the overall decline in breast cancer mortality over the past decade, about 40,420 deaths were expected in 2011 alone [1]. It has been shown that early detection of breast cancer is among the most effective ways to improve patient’s long-term survival [2]. The most common method for detection in clinical practice is X-ray mammography, which is generally effective for the broad population of women over 50 years of age. However, screening mammography has substantial limitations, primarily a high false-positive rate (up to 29%), which can result in unnecessary and costly surgical interventions [3]. Breast cancer detection is a particularly challenging problem in younger women and those with radiographically dense breasts. In these cases, the increased levels of fibroglandular tissue can easily obscure small tumors or masquerade as an abnormality because of the tissue overlap on plane film, and as a result, the overall diagnostic performance of mammography can be 1 significantly degraded [4]. Mammography has other drawbacks from the patient’s perspective, including uncomfortable and painful breast compression and exposure to ionizing radiation [5]. Other clinical standards, such as ultrasound and magnetic resonance imaging (MRI), have also been used to detect breast cancer. While both can achieve high spatial resolution, neither can provide information about the molecular-level changes occurring in breast tissue at the present time [6, 7]. In response to these shortcomings, alternative and/or complementary medical imaging modalities are being developed to improve both the sensitivity and specificity of current imaging tools, and to supply more functional information about the tissue properties comprising the breast. These properties can include electrical or optical characteristics, temperature, or tissue elasticity [5]. At Dartmouth, four different modalities, namely, microwave imaging spectroscopy (MIS), electrical impedance spectroscopy (EIS), near infrared imaging (NIR), and magnetic resonance elastrography (MRE), are being developed for breast cancer detection [8]. These alternative imaging modalities may prove to be reliable detection tools because they provide significant and consistent contrast between normal and malignant breast tissue [9-12]. Early microwave studies for biomedical applications showed a significant dielectric property contrast between normal and malignant breast tissues [13-15]; however, more recent data reported by Lazbnik et al. [16] has indicated that the properties of the normal breast are more variable than originally thought and that the contrast may not be as great for some types of breast tissue. This is particularly true for radiographically denser breast with higher concentrations of fibroglandular tissue [16]. Notwithstanding, early clinical microwave imaging studies on patients with suspected tumors has demonstrated significant discrimination between those with malignant cancers versus those with benign lesions and other normal tissues [17, 18]. 2 Compared to other conventional imaging modalities, such as X–ray mammography, the non–ionizing and non–compressive nature of MI makes it very attractive from the patient’s perspective [19]. In addition, MI has a relatively high sensitivity for detecting small tumors, and potentially high specificity to determine whether a suspicious area is malignant or benign at a significantly lower cost level compared to methods such as MRI and nuclear medicine [19]. 3 1.2. Overview of Different Microwave Imaging Methods At radio and microwave frequencies, several methods have been investigated for MI, particularly for breast cancer detection. These methods include: passive, hybrid, and active MI. These techniques are summarized in a tree diagram illustrated in Figure 1.1. Figure 1.1. Tree diagram of the current MI methods for breast cancer detection Passive MI is based on using microwave radiometry to measure temperature differences in the tissue. The principle of operation of this method is that there is a measurably greater temperature in tumor region compared with healthy breast tissue, which can be potentially detected using radiometers [20-25]. In hybrid methods, also called microwave–acoustic imaging, the principle of operation is to use microwave energy to illuminate the tissue and measure the pressure waves generated by the higher conductivity tumor regions using an ultrasound transducer [26-28]. In this case, the dielectric contrast selectively chooses the tumor for enhance absorption and vibration, which is subsequently detected by the ultrasound. Finally, the basic idea in the active MI is to illuminate the tissue with microwaves and measure the scattered or reflected signals around the tissue [29-33]. 4 Due to the dielectric property contrasts between normal and malignant tissue, in the presence of a tumor, microwave signals traveling through the tissue encounter a change in electrical properties, which causes the incident wave to scatter and create some variation in phase and the amount of detected energy at the receivers [19]. This information can be ultimately used to obtain the associated dielectric property maps of the tissue at microwave frequencies. Active MI approach may be classified as tomographic and confocal. Tomographic MI is a transmission–reflection imaging technique where the complete spatial dielectric property distributions of the object are obtained from the transmitted (incident) and scattered (received) electromagnetic fields [34]. Confocal MI, on the other hand, is a time-domain imaging technique that exploits the principles of radar synthetic focusing, and instead of recovering the complete dielectric property profiles of the object, only locations of significant scatterers are reconstructed [35]. The main idea behind the confocal MI is to create a map of microwave scattering arising from the contrast in dielectric properties within the imaging domain. This approach is very closely related to optical confocal systems; however, electromagnetic signals at microwave frequencies benefit from a greater penetrability in tissue [36]. After transmitting an ultra-short pulse into the tissue, the backscatter waveforms are collected, filtered, and time–shifted to create synthetic focal points. Returning waveforms from the scattering tissue add coherently at the focal points, while those surrounding the confocal points add incoherently [30, 37, 38]. Confocal microwave imaging was first presented by Hagness et al. in 1998 [37]. Traditionally, radar imaging is only a reflection technique and utilizes monostatic radar, where the signal is transmitted and received by the same antenna. In addition, due to the structural 5 complexity of the living tissue, more measurements are desired and as a result, sensor arrays have been introduced to increase the amount of measurement data [39, 40]. Radar-based imaging systems include microwave imaging via space time beamforming [41, 42], and tissue sensing adaptive radar (TSAR) [40, 43, 44]. In the former, the patient lies supine with the antennas scanned over the naturally flattened breast, and some advanced clutter reduction algorithms are used to calculate the location of the scatterer(s) within the imaging domain. On the other hand, in the TSAR systems the patient lies prone and the antennas are scanned around the breast, which is suspended through a hole in the examination table. In this case, simpler clutter reduction algorithms are used to create the corresponding TSAR images [40, 43, 44]. Transmission-based confocal imaging is a relatively new radar imaging technique where the transmitted signals are collected and used to exploit the varying attenuation induced by the electrical property differences within the tissue. Craddock et al at the University of Bristol has successfully developed a clinical radar–based UWB microwave system for breast cancer detection [45]. Their system uses a real aperture array of UWB antennas positioned on a section of a hemi-sphere (conforming to the breast’s curvature) and operates in a multi-static mode [45]. Microwave tomography commonly refers to imaging by cross-sectioning or slicing through an inhomogeneous object/tissue. The word tomography is derived from the Greek tomē ("a cutting") or tomos ("a cut" or "section") and graphein ("to write") [46]. This method aims at the complete reconstruction of the dielectric property distributions (permittivity and conductivity) of the tissue being imaged. It should be noted that despite recovering the complete 3D map of the actual dielectric properties of the tissue at once (instead of 2D slices), the three-dimensional MI is frequently referred to as 3D tomography. 6 Several research groups have been investigating the 2D and 3D tomographic MI in biomedical applications. Researchers at Duke have developed a fixed 3D array of antennas positioned throughout the surface of the imaging chamber and connected to the source and measurement devices through an RF switching network [47]. They use a fast volume integral equation method (diagonal tensor approximation, DTA) to accelerate the computation of the forward solution [48-50]. The 3D imaging systems and methods at Duke have been tested on several phantoms from experimental data, and their results indicate that the detection and imaging of inclusions with dielectric contrast is possible [51]. Similarly, other researchers have been studying microwave tomography as a diagnostic modality for non-invasive assessment of functional and pathological conditions of biological tissues. They have developed both 2D and 3D imaging systems [52-54], in addition to the image reconstruction methods [55-57] for experimental imaging of animals [53, 54, 58] and functional imaging of extremity soft tissues [59]. However, the prototypical imaging system developed at Dartmouth College [9] is the only MI system thus far which has been successfully progressed to clinical trial stage. Details of this system are provided in the following section. 7 1.3. Tomographic Microwave Imaging Several research groups have been investigating tomographic MI for breast cancer detection [60, 61], but to the best of our knowledge, we at Dartmouth College, are the only group to this date, who has developed an actual clinical system [9, 29]. Figure 1.2 shows our latest clinical microwave imaging system (MIST), located at Dartmouth Hitchcock Medical Center (DHMC). The imaging array consists of 16 monopole antennas positioned on a 15.2 cm diameter circle. They are placed on one of the two independently–moving plates, A or B, as shown in Figure 1.2(a). During data acquisition, the antenna array scans through a number of vertical positions in 0.5, 0.75, or 1.0 cm increments. This configuration enables us to collect both in-plane (i.e. when all antennas are at the same height) and cross-plane (i.e. when the two sets of antennas are at different heights) data. However, due to a time limitation and the lack of a viable, fast, and user– friendly three–dimensional (3D) image reconstruction procedure, only in-plane patient data is currently collected and the corresponding images are only reconstructed in two dimensions (2D). Sequentially, at each antenna position, that antenna acts as a transmitter while the response is measured at each of the remaining 15 antennas. The system is designed to operate over the frequency range from 500 MHz to 3 GHz. The patient lies prone on the bed with the breast to be examined suspended through an aperture in the top of the tank, as shown in Figure 1.2(b). The tank is filled with a coupling liquid, closely mimicking the average constitutive parameters of the breast [62], maximizing the amount of microwave energy coupled into the tissue. The imaging platform (c) and the space underneath (d) are shown in Figure 1.2. This space is used to house the data acquisition instrumentation, illumination tank modules, and the coupling medium reservoir. 8 (a) (b) (c) (d) Figure 1.2. MIST system: (a) antenna configuration on two independently moving plates A independently-moving (pink) and B (blue), (b) imaging tank, ((c) exam platform, and (d)) cabling and fluid reservoir underneath the table 9 1.4. Project Aims Question 1: Can a viable, fast, and user–friendly image reconstruction procedure be developed for 3D microwave imaging? To recover representative electrical property profiles of tissue from measurement data, matching the numerical model to the physical situation is one of the most important factors to ensure the correct interpretation of the data (by the model) and algorithmic convergence. Modeling of 3D electromagnetic wave propagation in a complex medium is generally computationally expensive, even with modern computational resources. Using 2D methods to model this inherently 3D phenomenon may save significant computational time; however, it likely imposes simplifications, which may introduce image artifacts. Our goal is to develop a practical 3D image reconstruction procedure, which is computationally feasible and balances the tradeoff between the accuracy and efficiency of the model. Other issues such as compatibility with previously developed reconstruction techniques should also be taken into account. In order to evaluate the performance of the reconstruction algorithm, a number of simple and complex– shaped 3D simulation and phantom experiments will be performed. Question 2: Will a priori structural information from other imaging modalities, such as MRI, improve the quality of microwave images and enrich the accuracy of the reconstructed property distributions? This question will be answered through a series of 2D and 3D simulations and phantom experiments. Moreover, different analyses including: frequency, noise level, shape, size and number of regions of interest, imperfect or false priors, and contrast levels will be performed to evaluate how prior spatial information of the object being imaged can enhance the quality and accuracy of the recovered property distributions. We will also extend the technique of using prior 10 structural information to 2D patient data. Our ultimate goal is to collect simultaneously 3D microwave data and the corresponding anatomical information from MRI. We believe this additional information will improve the accuracy of the reconstructed property distributions significantly, enabling us to characterize smaller regions of interest more precisely. 11 1.5. Proposed Work Currently our imaging system operates primarily in the 2D mode, where multiple vertical planes of data (in-plane configuration) are collected, and the corresponding 2D slice images are reconstructed. In this mode, the image reconstruction is based on the assumption that the scattering problem can be reasonably represented as a 2D problem [63]. Using 2D methods to model a 3D phenomenon significantly reduces the computation time; however, it may impose excessive simplifications, which can introduce image artifacts. This has been certainly used as a rationale for why other 2D imaging efforts by other researchers failed [64]. In the past, we have performed studies showing that the match between the 2D model and 3D actual measurements is quite good and is not a deterrent for reconstructing 2D images, especially for the phase [63, 65]. Certainly, this has been part a key reason we have been able to recover useful 2D images to this point. Previously, a 3D microwave image reconstruction procedure based on a finite–difference time–domain (FDTD) algorithm was developed at Dartmouth College [66]. The original implementation of this reconstruction algorithm was written in FORTRAN and a few simple– geometry simulations and phantom experiments were performed to evaluate the efficiency and accuracy of the results [67]. One of the aims of this work is to modify and redesign the existing algorithm and develop a not only viable and fast, but also a user–friendly 3D image reconstruction module which is computationally feasible and balances the tradeoff between the accuracy and efficiency of the model. In chapter 2 a full discussion about the 3D microwave image reconstruction, including the proposed redesigning and optimization methods, along with the results of several simulations, phantom experiments, and preliminary patient data are presented. 12 At the core of active MI is the inverse reconstruction algorithm. Since the inverse electromagnetic problem is inherently non–linear, iterative schemes are usually used in the image reconstruction process [68, 69]. Moreover, due to the ill–posed nature of the problem, regularization procedures are required to impose additional constraints on the reconstructed images [69, 70]. This is often accomplished by introducing prior information about the object being imaged, which is helpful to ensure convergence of the reconstruction algorithm to the correct electromagnetic property distribution [71]. We have previously demonstrated that we could recover clinically useful microwave tomographic images for breast cancer detection without using any prior spatial information about the tissue [9]. However, in order to enhance the quality of our previous results and recover more accurate dielectric property distributions, we propose a new approach combining the functional information of the MI (property contrast) with the high spatial resolution of other imaging modalities such as MRI and X–ray CT. This approach is based on developing new image reconstruction/regularization strategies that exploit structural information about the object being imaged through some constraints in the microwave property contrast [72]. In chapter 3 the mathematical framework of our microwave imaging reconstruction algorithm along with the new regularization schemes are presented. Then, our prototyped MI-MRI multi-modality imaging system is described in detail, and finally, the new regularization approaches are evaluated through several 2D and 3D simulations, phantom experiments, and initial clinical data. The final outcome of this project is two graphical user interfaces (GUI) for the 2D and 3D FDTD microwave image reconstructions. These GUIs are designed to enable the user to control and adjust some of the important parameters in the reconstruction procedure, and to visualize the 13 final reconstructed images more conveniently. A detailed description of the developed GUIs is presented in chapter 4. 14 2. 3D Microwave Image Reconstruction 2.1. 2D vs. 3D Modeling of the Electromagnetic Field In our image reconstruction algorithm, the numerical modeling of the propagation of electromagnetic (EM) waves is governed by Maxwell's equations. Solutions to Maxwell’s equations in 2D cylindrical waveguides can be approximated as linear combinations of the Bessel functions, known as Hankel functions of the first and second kind of nth order with an argument hr, Hn(hr), where h2 = γ2 + k2, γ is the propagation constant, k is the wave number, and r is the distance from the source to the observation point in the cylindrical coordinate system [73]. For the far field solution: π nπ 2 ± j hr− 4 − 2 lim ( H n (hr)) = e r→∞ π hr (2.1) and as a result, the field E is proportional to [73]: () E r ∝ 1 jhr e πr (2.2) On the other hand, in antenna and radiation problems in 3D space, solutions to Maxwell’s equations can be approximated using dyadic Green’s function G [74], which in turn can be expressed in terms of the scalar Green’s function g as: 1 G r, r′ = I + 2 ∇∇ g r, r′ k ( ) ( ) (2.3) where r and r′ are the Cartesian coordinates of the observation point and the source, respectively. Using the radiation field approximation, the electric field E is proportional to [74]: 15 () E r ∝ 1 jhr e πr (2.4) Based on the Poynting’s theorem, the magnitude of the electromagnetic radiation is 2 2 proportional to E [73]. From (2.2) and (2.4), E ∝ 2 E ∝ 1 in the 2D case, whereas for the 3D case, r 1 . In both cases, the phase varies as a function of r, the distance from the source to the r2 observation point, i.e. receiver location. Based on simulation and phantom experimentation, we have noticed a strong correlation between the magnitude of the electric field and the recovered conductivity values. We have observed a similar relationship between the phase and the recovered permittivity values. As discussed above, the phase varies as a function of the distance r in both 2D and 3D models, whereas the magnitude varies differently in the two models (1/r in 2D vs. 1/r2 in 3D). Therefore, we hypothesize that there exists a greater mismatch between the magnitudes of electric field modeled in 2D and 3D, as compared to the corresponding phase values. In addition, since the electromagnetic wave propagation is a 3D phenomenon, this may explain why we have been able to recover more accurate permittivity values using our 2D model. On the other hand, due to the greater mismatch between the actual and 2D modeled magnitudes of the electric field, overall the 2D reconstructed conductivity images contain more artifacts. In order to validate our hypothesis, we performed a simulation experiment as shown in Figure 2.1. A cylindrical inclusion of 1.0 cm radius with dielectric properties of εr = 40 and σ = 2 S/m was centered at (3, 0 cm). The permittivity and conductivity values of the surrounding background medium was set to εr = 24.0 and σ = 1.13 S/m, respectively. For the 2D simulated data, a cross section of the setup in Figure 2.1 was used, whereas for the 3D case, the entire 16 cylinder of 3.0 cm height was considered. In both cases, only one plane of data (16 transmitters × 15 receivers) was synthesized at 1100 MHz. In addition, the data was calibrated with respect to the homogeneous background. Figure 2.1. Simulation experiment setup with a cylindrical inclusion or 1.0 cm radius Figure 2.2 shows the calibrated amplitude/power (top) and phase (bottom) values of the synthesized data using the 2D and 3D models as a function of the relative receivers for a single transmitter (tx = 1). In both power and phase plots, the mismatch between the two models appears in the regions where the projection of the actual object (cylinder/circle) is recorded (i.e. between the receiver # 5 and 13). 17 Figure 2.2. Calibrated amplitude/power (top) and phase (bottom) values of the synthesized data using the 2D and 3D models as a function of the relative receivers for a single transmitter (tx = 1) In order to quantify the mismatch for all receivers and compare the 2D and 3D models, the RMS difference values of the amplitude and phase were calculated at each transmitter tx as: 15 D = V tx (V 2D rx − Vrx3D ) ∑ max(V ) − min(V ) 15 (2.5) rx=1 where V is the amplitude or phase, rx is the receiver number, and the superscript 2D and 3D correspond to the amplitude or phase values synthesized in the 2D and 3D model, respectively. The calculated RMS difference values are plotted in Figure 2.3. 18 Figure 2.3. RMS difference values of the amplitude and phase The results in Figure 2.3 confirm our hypothesis that there is a greater mismatch between the magnitudes of electric field modeled in 2D and 3D, as compared to the corresponding phase values. In fact, the amplitude RMS differences are between 2-3 times larger than the counterpart phase values. 2.2. Challenges With Respect to 2D Reconstruction Algorithm During microwave data acquisition, electromagnetic fields propagate through and scatter from the tissue in a three–dimensional (3D) fashion [75]. However, in order to reduce the computational complexity and to speed up the image reconstruction process, it is often assumed that the behavior of electromagnetic waves in 3D space can be represented by a simplified 2D model. While it benefits from less intensive computational demands, such assumptions, including the field confinement to only the transverse magnetic (TM) mode, can lead to an increased level of artifacts in the recovered dielectric properties. Moreover, since only in-plane data is collected in 2D imaging, if the region of interest is small enough to fall between two 19 consecutive imaging slices, the 2D reconstruction algorithm may not detect the target accurately. Therefore, in order to improve the accuracy and quality of reconstructed images, a viable 3D microwave image reconstruction scheme is desired. Several years ago, a 3D microwave reconstruction algorithm was developed by Fang et al [66]. The algorithm was based on the FDTD modeling of electromagnetic wave propagation and it was implemented in FORTRAN. In order to demonstrate its feasibility, the algorithm was tested on several simulated data and phantom experiments with simple–geometry configurations [66, 67]. Despite the numerical efficiency of the reconstruction scheme, the lack of transparent and flexible input–data selection process and the lack of a user–friendly interface prevented further development of a practicable 3D image reconstruction procedure. 20 2.3. Proposed method 2.3.1. User Interface MATLAB (The MathWorks, Inc., Natick, MA) is a high–performance language for technical computing, which integrates computation, visualization, and programming in a user– friendly environment. The MATLAB graphical user interface (GUI) is also a powerful tool exploiting advanced computational and graphical functions of MATLAB. However, for computationally complex problems such as our 3D microwave image reconstruction, MATLAB is inherently less efficient than other prominent programming languages, such as FORTRAN, C, or C++. Nonetheless, MATLAB enables us to evaluate our codes using “performance profiling” and take advantage of MEX–files to improve speed of computationally intensive portions of our program. MEX–files are MATLAB EXecutables, which provide an interface between MATLAB and subroutines written in FORTRAN, C, or C++. When compiled, MEX–files can be run from within MATLAB in the same way as any other built–in functions. To support the development of MEX–files, MATLAB offers external interface functions to facilitate data transfer between MEX–files and MATLAB, and the ability to call MATLAB functions from FORTRAN, C, or C++ code [76]. Finally, MATLAB offers the Parallel Computing Toolbox to parallelize high– level constructs, such as for–loops, using multi–core processors or computer clusters. Parallel computing is very appealing in 3D microwave image reconstruction, especially to solve the forward model, where the electromagnetic field has to be calculated for each source (in our case, for each transmitter), as well as calculating the Jacobian matrix [77]. The forward problem is certainly the most time consuming portion of the image reconstruction algorithm, but since it can be solved for each source independently, the field calculation corresponding to each antenna can 21 be parallelized. In an optimal situation, the forward model for all sources can be calculated concurrently by using as many core–processors/cluster–nodes as the number of sources. 2.3.2. 3D Data Acquisition Based on the current design of our imaging system as described in Section 1.3, the 16 monopole antennas are grouped into two interleaved arrays of 8 antennas, and they are placed on two independently–moving plates, which in turn are attached to separate pairs of opposing computer controlled motors positioned underneath the tank. Since all antennas on the same mounting plate move together and have the same height, the collected measurement data contains some redundancies. For example, as is illustrated in Figure 2.4, when antennas 2 and 4 on plate A act as transmitter–receiver pairs, redundant data is collected at all possible imaging planes of plate B. In theory, all redundant measured amplitudes and phases are supposed to be equal, simply because the same antennas at the same locations transmit and receive the signal. However, there are minor differences between these measurements in practice, mainly due to the contribution of multipath signals related to other neighboring antennas. In order to handle redundant data, several approaches can be considered. One of them is that each redundant transmitter–receiver pair can be considered a new pair by assigning a different Source ID. The main advantage of this approach is that no measured data is discarded and the contribution of external sources is also taken into account. However, this assumption could make the number of forward solutions larger, and as a result, the computational time could increase significantly. In addition, due to the fact that the difference between redundant data is minimal (usually less than 0.1 dBm in amplitude and 1 degree in phase,) this approach does not appear to be beneficial. 22 Figure 2.4. Schematic antenna configuration: Two interleaved arrays of 8 antennas with each sub-array (A and B) being able to move independently from the other In order to account for all measured data and also minimize the number of unnecessary forward solutions, we use a more efficient approach. The average of redundant amplitudes and phases are calculated and used for each repeated transmitter–receiver pair. 2.3.3. Input–File Format In order to develop an efficient, flexible, and transparent image reconstruction procedure, we designed a new input file format. Figure 2.5 shows a schematic summary of steps to create 3D image reconstruction input files from the MIST system output data. 23 Figure 2.5. Schematic summary of steps to create 3D image reconstruction input files from the MIST system output data The first step is data subtraction and phase unwrapping, where calibrated data is formed by subtracting calibration data from raw measurement data, and then, the phase is unwrapped. In the second step, MIST system output files are converted into a single and flexible format: Files 1, 2, and 3. Based on the measurement configuration, a unique Source ID is assigned to each possible antenna position. These IDs along with their corresponding antenna positions are recorded in File–1. (Figure 2.6) File–2 records the order in which measurements were acquired at each frequency. As illustrated in Figure 2.7, transmitter (TX) and receiver (RX) sequences corresponding to measurements are listed in a table, in which a unique Measurement ID is assigned to each row. It should be noted that there are redundancies in (TX, RX) pairs (circled in Figure 2.7), which need to be identified and handled as described in the previous section. 24 Figure 2.6. Input File–1 Figure 2.7. Input File–2 File–3 consists of a table containing calibrated and phase unwrapped data for each frequency. As shown in Figure 2.8, each row of this table includes a measurement ID, TX and RX IDs, along with their corresponding data in the last two columns. Similar to File–2, File–3 contains redundant measurements (circled pairs in Figure 2.8). In order to compare these redundancies, the second and third columns in File–3 are sorted in ascending order, thus a list of all data corresponding to each (TX, RX) pair appears in a block. 25 Figure 2.8. Input File–3 Files 1, 2, and 3 are a basis for all possible subsequent reconstructions. The next step is to generate input files that are free of redundancies and are ready to be used in the image reconstruction algorithm. These input files include File–1’, 3’ and 4. As illustrated in Figure 2.5 and Figure 2.9, File–3’ consists of two 2–dimensional matrices (one for amplitude and one for phase), which are formed based on File–2 and File–3. The size of each matrix is n_Sources × n_Sources, where n_Sources is the number of selected sources in File–1’. Entry ampi,j or phasei,j of each matrix indicates a unique amplitude or phase transmitted by Source ID Si and received at Source ID Sj. At this stage, data redundancy is eliminated, as described in Section 2.3.2. Moreover, if no data has been collected for a given (TX, RX) pair – an impossible pair, the corresponding entry of matrices in File–3’ are empty. For example, since each antenna cannot act as both transmitter and receiver simultaneously, the diagonals of the matrices are always empty. 26 Figure 2.9. Input File–3’ File–4 consists of a n_Sources × n_Sources matrix of flags (0’s or 1’s), which are linked to matrices in File–3’ (Figure 2.9). RX’s associated with 1 are used for the inversion in the reconstruction algorithm, whereas those associated with 0 are not. The flagged matrix entries corresponding to empty cells in File–3’ (diagonal and impossible pairs) always remain zero, simply because no data has been collected at those (TX, RX) pairs. Prior to image reconstruction, the user can define what data will be used in the reconstruction algorithm by activating (switching the flags on) or deactivating (switching the flags off) receivers. This process is done through a MATLAB GUI, where each checked box represents a flag corresponding to (TX, RX) pair (i, j), as shown in Figure 2.10. Diagonal and impossible pairs of the matrix are disabled so the user cannot switch them on. A few buttons, including SELECT ALL, CLEAR ALL, SELECT IN–PLANE, and MULTI PLANES are also designed to help the user make selections faster and more efficiently. Based on the number of 27 planes used to acquire measurements, the user can select a specific set of in-plane data (for example P1, P2, …) by checking/un–checking the corresponding checkboxes, in Figure 2.10. Once the transmitter/receiver selection process is done, the user presses CONTINUE and the flag file (File–4) will be generated. This file is used in the reconstruction algorithm to select specific transmitter/receiver pairs for the inversion process. Based on the selected transmitter/receiver pairs in File–4, some antennas may not be selected to act as either transmitter or receiver. Therefore, in order to reduce the size of input data and avoid unnecessary forward solution calculation, the unused antennas (sources) in File–4 are removed from the antenna list in File–1 and the new list of sources with their spatial coordinates is saved into File–1’. Accordingly, Files–3’ and 4 are updated with the new source IDs. Figure 2.10. Input data selection GUI in MATLAB 28 2.3.4. Optimizing the Reconstruction Algorithm As described in Section 2.3.1, the forward solution for each source can be computed independently. In order to speed up this process, they are calculated in parallel. Once forward solutions are obtained, the Jacobian matrix is formed. The components of the Jacobian matrix are partial derivatives (or the sensitivities) of the computed electric field with respect to the complex contrast in the object. Depending upon the number of sources and the resolution of reconstruction mesh, the computation of the Jacobian matrix can be intensive in MATLAB. In order to speed up this process, we use MEX–files to populate the Jacobian matrix. Figure 2.11 illustrates a schematic summary of steps in the image reconstruction algorithm. In the current study, we evaluate several important aspects of 3D microwave imaging in simulations, phantom experiments, and some preliminary patient data. Figure 2.11. Schematic summary of reconstruction algorithm 29 2.4. Results – Simulations Simulated 3D measurement data is generated by our FDTD forward solver for different shaped and different sized target inclusions. In order to simulate the system noise, –100 dBm synthetic noise is also introduced to the data. In addition, in order to avoid so called “inverse crime,” different meshes are used for the forward and inverse solvers [78]. More specifically, the finite difference grid used for generating the simulated data is 50% denser than the one used for the reconstruction. Figure 2.12. Schematic off the imaging domains evaluated: The background diameter was 14 cm (antennas are positioned on a 15.2 cm diameter). The spherical inclusion of 1.5 cm radius was centered at (3, 0, 0 cm). Figure 2.12 shows a schematic setup of the simulation experiment experiments studied in this section. 112 monopole antennas are configured in 7 evenly–spaced circles of 15.2 cm diameter and 1 cm separated from each other. A spherical inclusion of 1.5 cm radius with dielectric properties of εr, 30 Inc = 40.0 and σInc = 2.0 S/m and centered at (x, y, z) = (3, 0, 0 cm) is embedded in a background medium with dielectric properties of εr, BK = 22.4 and σBK = 1.23 S/m. 2.4.1. FDTD Grid Density We use a dual-mesh approach for the image reconstruction procedure, where the forward solution is calculated on a rectangular uniform FDTD grid, while the distribution of the constitutive parameters are reconstructed on a tetrahedral-element mesh [67]. It is generally assumed that a resolution of 10-20 grid points per wavelength is sufficient to guarantee acceptable convergence for forward solution [79]. In this section, we study the effects of the finite difference (FD) grid density on the recovered dielectric property distributions. The permittivity (ε) and permeability (µ) determine the speed of propagation of a wave (c) through a given medium by c= 1 (2.6) εµ The wavelength (λ) is determined by λ= 2π c ω = 2π (2.7) ω 2εµ where ω is the angular frequency. To determine the size of the FD grid for the forward solution, we use the already known background medium permittivity value to calculate the wavelength. The resolution of the FD grid (∆x, ∆y, and ∆z) is then calculated by 31 ∆x = ∆y = ∆z = λ (2.8) D where D is the number of cells per wavelength. Table 2.1 contains a list of FD grid spacing used to compute the forward solution for the simulation experiment described in Section 2.4. Table 2.1. Finite difference grid spacing (in mm) used to compute the forward solution for the simulation experiment, along with the corresponding reconstruction times per iteration Case D ∆x = ∆y = ∆z (mm) Reconstruction Time (min)/ iteration (a) (b) (c) (d) (e) (f) 30 25 20 15 10 5 1.52 1.83 2.29 3.05 4.58 9.16 16 9 4.5 2.5 1.7 0.8 Figure 2.13 shows the corresponding reconstructed images at 1300 MHz on a mesh composed of 3594 uniformly distributed nodes and 15095 tetrahedral elements. 32 (a) D = 30 (b) D = 25 (c) D = 20 (d) D = 15 33 (e) D = 10 (f) D = 5 Figure 2.13. 1300 MHz contour slice images extracted from the reconstructed permittivity (top) and conductivity (bottom) profiles for the simulation experiment using the finite difference grid spacing cases (a) to (f) in Table 2.1. In all 6 cases above, both permittivity and conductivity profiles of the target inclusion are successfully recovered from those of the background medium. In order to compare the recovered properties with their true values, the reconstructed dielectric profiles in Figure 2.13 were extracted at the center of the inclusion (3, 0, 0 cm), as well as the center of the imaging domain (0, 0, 0 cm) – i.e. background region, and they are plotted as a function of the number of cells per wavelength (D in equation 2.8) in Figure 2.14. 34 Figure 2.14. Reconstructed permittivity (left) and conductivity (right) profiles along with the true property values, as a function of the number of cells per wavelength (D) in the forward solution For D = 10, the conductivity values of the target in Figure 2.14 does not seem to have reached an asymptote – i.e. it is still getting better, albeit slowly, for even larger D. However, neither reconstructed permittivity nor conductivity values of the background medium and the target inclusion appear to be significantly sensitive to the size of the FDTD grid spacing when over 10 cells per wavelength (D > 10) are used to calculate the forward solution. Therefore, in order to guarantee that the FDTD grid size does not affect the recovered property distributions, without significantly increasing the total reconstruction time, we use 20 cells per wavelength in the forward solver. 2.4.2. Number of Reconstruction Nodes The total reconstruction process time is directly proportional to the number of nodes (N) in the property reconstruction mesh. More specifically, N is the size of the solution of the inverse problem, where dielectric property distributions are calculated. Therefore, the larger N is, the more computationally intensive the image reconstruction becomes. On the other hand, as N increases, dielectric property values are calculated at more nodes within the imaging domain. As 35 a result, the reconstructed image has a higher quality, and the property distributions are potentially more accurate. In order to quantitatively study the sensitivity of the recovered property profiles to the parameter reconstruction mesh density, we reconstructed 3D dielectric property profiles of the simulation experiment described in the previous section on different reconstruction meshes. A list of these meshes (including the number of nodes and elements) is shown in Table 2.2. Table 2.2. Parameter reconstruction meshes used to solve the inverse problem for the simulation experiment Case Number of Nodes Number of Elements (a) 600 2275 (b) 1587 7545 (c) 2611 11378 (d) 3767 17761 (e) 4608 22184 (f) 6095 29796 (g) 7669 39372 (h) 9321 47978 In order to exclude the effects of other factors, all parameters except the reconstruction mesh were kept the same in this experiment. Figure 2.15 shows the reconstructed property profiles for cases (a) through (h) in Table 2.2. 36 (a) 600 Nodes (b) 1587 nodes (c) 2611 Nodes (d) 3767 nodes 37 (e) 4608 Nodes (f) 6095 nodes (g) 7669 Nodes (h) 9321 nodes 38 Figure 2.15. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) profiles for the simulation experiment on (a) 600, (b) 1587, (c) 2611, (d) 3767, (e) 4608, (f) 6095, (g) 7669, and (h) 9321 node meshes. Although both permittivity and conductivity profiles of the target inclusion are effectively detected from the background medium in all cases, the recovered dielectric properties appear to be improved as the number of nodes in the reconstruction mesh increases. In order to verify this observation, the reconstructed dielectric profiles in Figure 2.15 were extracted at the center of the inclusion (3, 0, 0 cm), as well as at the center of the imaging domain (0, 0, 0 cm) – i.e. background region. Figure 2.16 shows the extracted property values along with the exact solutions as a function of the number of reconstruction nodes in Table 2.1. Furthermore, for each case, the total reconstruction time per iteration was recorded and they are plotted as a function of the number of reconstruction nodes in Figure 2.17. Figure 2.16. Reconstructed permittivity (left) and conductivity (right) profiles along with the true property values, as a function of the number of reconstruction nodes 39 Figure 2.17. Total reconstruction times (min) per iteration as a function of the number of reconstruction nodes Results in Figure 2.16 confirm that the recovered dielectric properties are more accurate as the number of reconstruction nodes increases. However, the changes in the property values are not significant when reconstruction meshes with over 3000 nodes are used. Moreover, the reconstruction times in Figure 2.17 increase quadratically with the number of reconstruction nodes. Therefore, for a configuration of this size, it appears to be sufficient to use a mesh composed of approximately 3000-4000 nodes. 2.4.3. 3D images from in-plane data vs. stack of 2D images During the data acquisition, when the transmitting antenna is at the same height as the 15 receiving antennas – i.e. both plates A and B in Figure 1.2(a) are at the same level, in-plane data is collected. On the contrary, when they are at different levels, cross-plane data is collected. A complete set of data (full-data) includes both in-plane and cross-plane data at multiple vertical positions. In 3D imaging, data can be acquired in either in-plane or a combination of in-plane and cross-plane forms. 40 On the other hand, only in-plane data is collected and reconstructed for 2D imaging. Therefore, the only difference between a 3D reconstructed image from in-plane data and a stack of 2D images reconstructed from the same set of data is the forward model and the size of the respective inversion problems. Propagation of electromagnetic waves in a complex medium is a 3D phenomenon, and using 2D methods to model such a phenomenon is under the assumption that the scattering problem can be reasonably represented as a 2D problem. Nonetheless, this assumption may impose over simplifications on the reconstruction process that consequently introduce image artifacts. In order to study the effects of using a 2D model for a 3D problem, we have used the same set of in-plane synthetic data in Section 2.4 to reconstruct images in 2D and 3D on uniformly distributed triangular and tetrahedral element meshes composed of 473 and 4608 nodes, respectively. (a) 2D 41 (b) 3D Figure 2.18. (a) 2D and (b) 3D reconstructed permittivity (top row) and conductivity (bottom row) slices (in z-direction) of the simulation experiment with -100 dBm added noise at 1300 MHz, using only in-plane data Figure 2.18(a) and (b) show 7 slices of the corresponding 2D and 3D reconstructed images, respectively. In terms of the location and recovered dielectric properties of the target inclusion, both sets of images are very similar and they closely match to the exact values. However, the correct size of the spherical inclusion is recovered more precisely in the 3D case. More specifically, since the radius of the inclusion was 1.5 cm and it was centered at z = 0 cm, the target should only appear in slices between z = –1.5 and z = 1.5 cm. While only a small trace of the inclusion shows on the subsequent 3D slices of z = –2 and z = 2 cm, a significantly higher mismatch appears on the corresponding 2D slices. This effect is even more prominent in 2D conductivity images, where higher property values arise on the first and last slices too. This confirms what Fang et al previously showed that due to the data-model mismatch, the 2D reconstructed microwave images may show traces of the object in planes that extend beyond its actual locations [67]. In addition to the size, dielectric properties of the target inclusion are also 42 slightly higher and closer to the exact values when the 3D model is used to reconstruct the images. Based on the observations above, the 3D FDTD model seems to be superior to its 2D counterpart. Nonetheless, the recovered 2D images are quite accurate, which makes our 2D microwave imaging a compelling alternative to 3D imaging when limited computational resources are available. 2.4.4. Data Selection for Reconstruction As described in Section 2.3.3, during the data selection process in the 3D GUI, several predefined data selection options, including full-data, in-plane, and multi-plane are available. In the current configuration of our imaging system, if data is acquired in NP planes, NP×16×15 2D inplane measurements are collected, whereas the equivalent number of collected in-plane and cross-plane measurements is NP×16×((NP-1)×8). In 2D, the number of measurements increases linearly with NP, however, in 3D the increase is exponential. If each antenna had the flexibility to move independent of other antennas, (i.e. if instead of only two independent moving plates, they were mounted on 16 independent moving plates) the number of 3D measured data would be about 1.8 times larger (NP×16×NP×15). The increased number of 3D measurements amplifies the computational complexity of the forward problem, which is the most time-consuming part of the reconstruction algorithm. However, using cross-plane data may potentially enable us to better characterize smaller targets that can fall between two consecutive in-planes (2D slices). Nonetheless, cross-plane data for receivers close to the transmitter might not be very useful, especially when separated by multiple levels. In order to study the effects of different data selections on the 3D reconstructed dielectric property distributions, images of the synthetic data described in Section 2.4 were reconstructed 43 using full-data and multi-plane data with two consecutive planes. Figure 2.19 shows 7 slices of the corresponding reconstructed images. (a) full-data (b) multi-plane data with two consecutive planes Figure 2.19. 3D reconstructed permittivity (top row) and conductivity (bottom row) slices (in zdirection) of the simulation experiment with -100 dBm added noise at 1300 MHz, using (a) the full-data (all 7 planes of data, including both in-planes and cross-planes), and (b) multi-plane data with two consecutive planes In terms of the size and location of the target, both sets of images are identical to Figure 2.18(b), where only in-plane data was used for reconstruction. However, the recovered properties 44 of the inclusion are slightly more accurate as the number of measured data increases from inplane to full-data set selection. More specifically, when full-data set is used for the reconstruction, the elevated inclusion region in the middle slice (z = 0 cm) in Figure 2.19(a) reaches its maximum value in both permittivity and conductivity images. As less data is progressively used in multi-plane and in-plane data selections, the corresponding inclusion values decrease accordingly. Nonetheless, the increased amount of data used in full-data set selection adds significantly to the computational complexity of the forward problem, which might not be fully compensated for by the subtle differences in the reconstructed property profiles. Moreover, the amount of background artifacts increases as more data is selected for the reconstruction. This effect is even more prominent in the conductivity images. Therefore, using multi-plane data seems to be an effective selection to tradeoff between the accuracy of the recovered properties and computational demand of the reconstruction process. In addition, when limited data or scarce computational resources are available, using only in-plane data can be sufficient to obtain reliable 3D dielectric property distributions. 2.4.5. Number of Iterations Due to the iterative nature of the reconstruction procedure, a stopping criterion is needed to terminate the algorithm. In our 2D reconstruction algorithm, convergence is typically achieved within less than 15 iterations. However, we allow all reconstructions to execute for 20 iterations [80] as a simple way to ensure convergence in each case and verify that the relative error is below a threshold of 0.2 (in most cases). In order to establish an analogous norm to terminate the 3D reconstruction algorithm, the relative error values of the simulation experiment in the previous section are plotted as a function of the number of iterations in Figure 2.20. 45 Figure 2.20. 3D Relative errors as a function of the number of iterations from a simulation experiment The relative errors reach the threshold of 0.2 after 17 iterations, and they do not decrease significantly after about iteration number 20. In order to verify this observation with the reconstructed property values, the relative permittivity and conductivity difference values with respect to the 20th iteration were calculated at each node using the following formula: D20,N = V20 − VN (2.9) where V20 and VN are the reconstructed property values at iteration number 20 and N, respectively. Figure 2.21 shows the relative permittivity (top) and conductivity (bottom) difference values for iterations 40, 60, 80, and 100 (i.e. D20,40, D20,60, D20,80, and D20,100) at all reconstructed nodes. For statistical purposes, the mean of the difference values, along with ± 2 standard deviations are also plotted in Figure 2.21. 46 (a) D20,40 (b) D20,60 47 (c) D20,80 (d) D20,100 Figure 2.21. Relative permittivity (top) and conductivity (bottom) difference values for iterations (a) 40, (b) 60, (c) 80, and (d) 100 (i.e. for D20,40, D20,60, D20,80, and D20,100, respectively), along with the mean of the difference values ± 2 standard deviations (i.e. 95% of all difference data) Overall, all difference values are very close to zero, which indicates that there are minor variations in the reconstructed property values after 20 iterations. More specifically, 95% of all difference data points fall in a range of 1 unit for relative permittivity and 0.05 S/m for conductivity profiles at all subsequent iterations, suggesting that similar to the 2D algorithm, the stopping criterion of 20 iterations is generally a reliable means of ensuring convergence in the 3D reconstruction algorithm. 48 2.5. Results – Phantom Experiments In order to study the behavior of our new 3D reconstruction process with experimental data, our clinical MIST system was used to perform several phantom experiments with different geometries and setups. 2.5.1. Two Cylindrical Inclusions In the first phantom experiment, two cylindrical inclusions (4.2 cm radius circular-based and 1.0 cm side-length square-based cylinders) were filled with 90:10 and 50:50 mixtures of glycerin:water, respectively. As illustrated in Figure 2.22, the square-based cylinder was tilted toward the right side of the imaging domain, and both inclusions were immersed in a coupling medium composed of 80:20 glycerin:water mixture. (a) (b) Figure 2.22. Setup for the phantom experiment with two cylindrical inclusions: (a) The squarebased inclusion tilted toward the right side of the imaging domain, (b) Both inclusions immersed in a coupling medium and surrounded by the antenna array 49 The corresponding permittivity and conductivity values of the phantom were measured using a dielectric probe (Agilent, Santa Clara CA, 85070E Dielectric Probe Kit – High temperature probe,) and they are reported in Table 2.3. Table 2.3. Independently measured dielectric properties of the phantom experiment with two cylindrical inclusions Frequency (MHz) εr, BK σBK εr, CircInc σCircInc εr, SqInc σSqInc 1100 24.8 1.12 9.7 0.78 53.7 1.15 1300 21.0 1.2 8.3 0.83 51.2 1.48 The phantom was then imaged in our clinical system (MIST) with the antenna array set to collect 3D data in 5 vertical steps of length 1 cm. 3D images were reconstructed at 1100 MHz on a cylindrical mesh composed of 3132 uniformly distributed nodes and 12707 tetrahedral elements. The reconstructed permittivity profile using multi-plane data (two consecutive planes), along with the iso-surfaced inclusion regions are shown in Figure 2.23. It should be noted that due to the small contrast between conductivity values of the inclusions and the background medium, the corresponding reconstructed conductivity image is not shown here. 50 Figure 2.23. 1100 MHz reconstructed permittivity values of the 3D phantom experiment with two cylindrical inclusions. The iso-surface thresholds are εr,CircInc = 10.0 and εr,SqInc = 33.0. Both circular-based (with lower permittivity values than those of the background) and square-based (with higher permittivity values than those of the background) cylindrical inclusions are detected successfully. Moreover, the location of the inclusions and the tilting of the square-based cylinder are reconstructed correctly. However, in terms of the property values, the recovered permittivity values of the square-based cylindrical inclusion are underestimated, which might be due to the small size of the cylinder and the large contrast between the permittivity values of the inclusion and those of the surrounding background medium. In order to assess the effects of different data selection schemes on the recovered property profiles, four data selection scenarios were used for 3D reconstructing the phantom experiment: (a) full-data set (i.e. all acquired data), (b) 5 multi-plane data sets with two consecutive planes, (c) 5 in-plane data sets (one reconstruction using 5 in-planes), and (d) 5 individual in-plane data sets (5 separate reconstructions, each using only one in-plane data set). In addition, for comparison purposes, the 5 individual in-plane data sets were used to reconstruct images in 2D (i.e. 5 separate 2D reconstructions at the matching vertical positions). Figure 2.24(a) – (d) show the extracted 2D slices of the 3D reconstructed images at 1300 MHz corresponding to the data 51 selection scenarios (a) – (d) described above. The equivalent 2D reconstructed slices are shown in Figure 2.24(e). (a) 3D – full-data (b) 3D – 5 multi-plane 52 (c) 3D – 5 in-plane (d) 3D – 5 individual in-plane (e) 2D – 5 individual in-plane 53 Figure 2.24. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) slices of the phantom experiment with two cylindrical inclusions using 3D reconstruction algorithm with (a) full-data set, (b) 5 multi-plane data sets with two consecutive planes, (c) 5 in-plane data sets, (d) 5 individual in-plane data sets, and (e) using 2D reconstruction algorithm with 5 individual inplane data sets In terms of the location and property values of the phantom, the reconstructed 3D images (a) – (d) are not significantly different, and they are all superior to the equivalent 2D images in (e). Overall, the more data that is used for the reconstruction, the more features show on the reconstructed images. The additional features do not all necessarily add to the accuracy of the recovered property values. In fact, the increased level of artifacts in (a), for example, may be due to the larger amount of measured data used in the full-data scheme. In Figure 2.24(c) and (d), the main difference between the selected data is that in (c) all 5 in-plane data sets are used at once to reconstruct a single 3D image on a cylindrical mesh; whereas in (d), each of the 5 in-plane data sets is used separately to reconstruct five 3D images on the same cylindrical mesh. Therefore, the amount of measured data used in (d) is five times smaller than that used in (c), which may explain the lower detail-level observed in Figure 2.24(d). 2.5.2. Breast-shaped Phantom Experiment with an Spherical Inclusion In order to perform a more anatomically-shaped phantom experiment, a spherical target inclusion of 2.0 cm radius made from a saline gel was suspended in a prototyped plastic breast model, as shown in Figure 2.25(a). The breast model was filled with a 88:12 glycerin:water mixture, and it was then submerged in the imaging tank filled with a matching liquid composed 54 of a 80:20 mixture of glycerin:water, as illustrated in Figure 2.25(b). The independently measured dielectric property values of the coupling medium, breast model, and the target inclusion at 1300 MHz are reported in Table 2.4. (a) (b) Figure 2.25. Breast-shaped phantom experiment setup: (a) The spherical saline gel target inclusion suspended in the plastic breast model, (b) Rapid–prototyped plastic breast model submerged in the imaging tank Table 2.4. Independently measured dielectric properties of the background coupling medium (BK), breast model (BR), and target inclusion (Inc) at 1300 MHz Frequency (MHz) εr,BK σBK εr,BR σBR εr,Inc σInc 1300 22.3 1.32 15.1 1.01 67.8 1.80 In the present phantom experiment, 3D microwave data was acquired in multiple planes. More specifically, the antenna array transmitted and received the signal at 9 equally–spaced (1 cm) vertical positions, and as a result, 9 × 9 × 240 measurements were collected. For comparison purposes, selective data was used to reconstruct images both in 2D and 3D. In the 2D 55 reconstruction, only in-plane data (with all 16 antennas at the same height collecting a total of 9 × 240 measurements) was used, whereas in the 3D reconstruction, several data selection scenarios including full-data, in-plane, and multi-plane (two consecutive planes) were considered. Figure 2.26 shows the 1300 MHz reconstructed images using the multi-plane data set on a 3D mesh composed of 4849 uniformly distributed nodes and 21381 tetrahedral elements. For illustration purposes, the 3D iso–surface of the breast and target inclusion are displayed in Figure 2.26. (a) (b) Figure 2.26. 3D reconstructed (a) permittivity and (b) conductivity images, using the multi-plane data set at 1300 MHz. The iso-surface thresholds are εr,BR = 15.0, σBR = 1.0 S/m, and εr,Inc = 65.0 and σInc = 2.3 S/m. In terms of the location, the target inclusion as well as the breast region is accurately detected in both permittivity and conductivity images. However, the recovered properties of the target inclusion are overestimated. More specifically, the recovered permittivity and conductivity of the target at the center of the inclusion are 91.5 and 2.46 S/m, respectively, which indicates 56 about 35% overestimation of dielectric properties. Moreover, the size of the reconstructed inclusion permittivity appears to be smaller than that of the conductivity. In order to compare the 3D reconstructed images with those reconstructed in 2D, 9 slices from the 3D volume in Figure 2.26 were extracted and are shown in Figure 2.27. Similarly, Figure 2.28 shows the 2D reconstructed dielectric properties at 9 equally–spaced (1 cm) vertical positions. In both figures, the images in the top row correspond to permittivity, whereas those in the bottom row represent the recovered conductivity values. 57 Figure 2.27. Extracted slices from the 3D reconstructed images in Figure 2.26:: Permittivity (top row) and conductivity (bottom row) 58 Figure 2.28. 2D reconstructed econstructed dielectric properties: Permittivity (top row) and conductivity (bottom row) Similar to the 3D reconstructed images iin Figure 2.27,, the breast as well as the target inclusion is successfully detected in the 2D images in Figure 2.28. In addition, as the reconstructed planes move toward the nipple (i.e. z decreases), the outline of the breast model becomes smaller and finally vanishes at about z = –2 cm in both sets of images. Nonetheless, the object seems to appear at different heights in the permittivity and conductivity images (at z = 4 cm in εr and at z = 3 cm in σ). In addition, the level of artifacts, especially in the background 59 region, is significantly increased in the 2D reconstructed images. This effect is even more prominent in the corresponding conductivity profiles. In terms of the location of the target, since the radius of the spherical inclusion was 2.0 cm, it should only appear in the first three reconstructed images from the top of the breast model (z = 4, 3, and 2 cm). While the 3D reconstructed images capture the correct location of the inclusion, the 2D reconstructed conductivity images show traces of the target at lower planes (as low as z = 1 cm), which indicates that the 3D reconstructed conductivity images are superior to those reconstructed in 2D, and confirms that due to the data-model mismatch, the 2D reconstructed microwave images may show traces of the object in planes that extend beyond its actual locations. Similar to the analysis performed in Section 2.5.1, three different data selection schemes including (a) full-data set (i.e. all acquired data), (b) 5 multi-plane data sets with two consecutive planes, and (c) 5 in-plane data sets were used for reconstructing the breast-shaped phantom experiment at 1300 MHz. Figure 2.29(a) – (c) show the extracted 2D slices of the 3D reconstructed images (permittivity on top and conductivity on bottom) at z = 3 cm, corresponding to the data selection scenarios (a) – (c), respectively. The equivalent 2D reconstructed slices are shown in Figure 2.29(d). Likewise, Figure 2.30 shows the counterpart reconstructed images at z = 0 cm, where only the distal parts of the breast are captured. 60 (a) (b) (c) (d) Figure 2.29. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) slices of the breast-shaped shaped phantom experiment at z = 3 cm using the 3D reconstruction algorithm with (a) full-data set, (b) 9 multi-plane data sets with two consecutive planes, (c) 5 in-plane plane data sets, and (d) using the 2D reconstruction algorithm with 9 individual in-plane data sets (a) (b) (c) (d) Figure 2.30. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) slices of the breast-shaped shaped phantom experiment at z = 0 cm using the 3D reconstruction algorithm with (a) 61 full-data set, (b) 9 multi-plane data sets with two consecutive planes, (c) 5 in-plane data sets, and (d) using the 2D reconstruction algorithm with 9 individual in-plane data sets Overall, the 3D reconstructed images are superior to those reconstructed in 2D, largely because the outline of the breast is more visible and the target inclusion is characterized more accurately. The property profiles in multi-plane and in-plane cases, (b) and (c), are very similar and they both present less background artifacts than the full-data case in (a). Despite the fact that the target inclusion is correctly detected in both permittivity and conductivity images in Figure 2.29, the size mismatch between the reconstructed inclusion permittivity and conductivity profiles can be seen in all cases, including the 2D reconstructions. This mismatch is even more prominent in (a) where full-data set was used for reconstruction, which may be due to the larger amount of measured data used in the full-data case. 62 2.6. Results – Initial Clinical Data Until very recently, all of our clinical breast data had been collected and reconstructed only in the 2D format (in-plane 2D slices). Since we have developed a viable and user-friendly interface to select and reconstruct 3D measurement data, we are planning to extend our data acquisition from 2D (only in-plane) to 3D (both in-plane and cross-plane). In Section 2.6.1, we will evaluate the performance of the reconstruction algorithms by comparing the 2D and preliminary 3D reconstructed breast images of a normal subject. In Section 2.6.2, we will present the 3D results for the initial clinical bone imaging, using a simple adaptation of our existing microwave imaging system. 2.6.1. Breast Imaging In this study, a normal subject (patient# 551, a volunteer with no known breast abnormalities) of 45 years old was imaged in our clinical system (MIST) using the 3D data acquisition, as shown in Figure 2.31. Based on the patient’s size and age, a matching liquid composed of a 80:20 mixture of glycerin:water was used, and six planes of microwave data (both in-plane and cross-plane, separated by 1 cm) were collected. The images were reconstructed both in 2D and 3D on triangular- and tetrahedral-element meshes, respectively. The 2D circular reconstruction mesh was comprised of 559 nodes and 1044 triangular elements with a 6.9 cm radius, whereas the 3D cylindrical reconstruction mesh had a radius of 6.9 cm, a height of 6 cm, and consisted of 6095 nodes and 29796 tetrahedral elements, respectively. 63 Figure 2.31. 3D microwave imaging of a normal subject in the clinical MIST system For comparison purposes, selective data was used to reconstruct images both in 2D and 3D. In the 2D reconstruction, only in-plane data (with all 16 antennas at the same height collecting a total of 6 × 240 measurements) was used, whereas in the 3D reconstruction, several data selection scenarios including full-data, multi-plane (four and two consecutive planes), and inplane were considered. Figure 2.32 shows 6 slices from the 3D reconstructed volume at 1100 MHz using the full-data set. Similarly, Figure 2.33 shows the 2D reconstructed dielectric properties at 6 equally–spaced (1 cm) vertical positions. In both figures, the images in the top row correspond to permittivity, whereas those in the bottom row represent the recovered conductivity values. 64 Figure 2.32. Extracted slices from the 1100 MHz 3D reconstructed images using the full-data set: Permittivity (top row) and conductivity (bottom row) Figure 2.33. 1100 MHz 2D reconstructed econstructed dielectric properties: Permittivity (top row) and conductivity (bottom row) Overall, the permittivity images appear to be superior to the conductivity counterparts, counterpar and despite the low contrasts between the dielectric properties of the coupling medium and the heterogeneously dense breast tissue, the outline of the breast is roughly detected in both 3D and 2D reconstructed images in Figure 2.32 and Figure 2.33.. In addition, as the reconstructed planes 65 move toward the nipple (i.e. z decreases), the contour of the breast becomes smaller and finally vanishes at about z = –1.5 cm in both sets of images. Based on the patient’s age and her breast density type, the elevated dielectric properties (especially the permittivity) at the right center of the reconstructed images closer to the chest wall may correspond to the fibroglandular tissue. While the 2D and 3D reconstructed permittivity images are comparable, the conductivity counterparts are more different. In fact, the higher property values of the fibroglandular region only appear in the first three planes (i.e. z = 2.5 to 0.5 cm) in the reconstructed 3D images, as well as in the 2D reconstructed permittivity profiles. However, this region arises in the lower planes in the 2D reconstructed conductivity images, which confirms one more time that due to the greater mismatch between the actual and 2D modeled magnitudes of the electric field, the 2D reconstructed conductivity images contain more artifacts, and may show traces of the tissue (in this case fibroglandular) in planes that extend beyond its actual location. In order to study the effects of different data selection schemes on the 3D recovered dielectric properties, four different cases were used for reconstructing the patient data described above: (a) the full-data set (i.e. all acquired data), (b) 6 multi-plane data sets with four consecutive planes (two planes above and two planes below), (c) 6 multi-plane data sets with two consecutive planes (one plane above and one plane below), and (d) 6 in-plane data sets. Figure 2.34(a) – (d) show the extracted 2D slices of the 3D reconstructed images (permittivity on top and conductivity on bottom) at z = 1.5 cm, corresponding to the data selection scenarios (a) – (d), respectively. The equivalent 2D reconstructed slices (using only one in-plane data set at z = 1.5 cm) are shown in Figure 2.34(e). It should be noted that for comparison purposes, the range of the recovered permittivity values in Figure 2.34 are different from those in Figure 2.33. 66 (a) (b) (c) (d) (e) Figure 2.34. 1100 MHz reconstructed pe permittivity rmittivity (top) and conductivity (bottom) slices of the breast-shaped shaped phantom experiment at z = 1.5 cm using the 3D reconstruction algorithm with (a) the full-data (FD) set, (b) 6 multi multi-plane data sets with four consecutive planes (4CS) (two planes above and two planes below), (c) 6 multi-plane data sets with two consecutive planes (2CS) (one plane above and one plane below), (d) 6 in-plane (IP) data sets,, and (e) using the 2D reconstruction algorithm with 6 individual in-plane data sets. At a first glance, all 3D reconstructed images in Figure 2.34 seem very similar and there are minor differences between them and their 2D counterparts counterparts. The outline off the breast appears to be more distinguishable in the 3D reconstructed permittivity images, compared to that reconstructed in 2D. In general, the 2D reconstructed images look sharper and appear to have higher quality compared to the 2D slices of the 3D re reconstructed images.. This is due to the fact that the 2D slices of the 3D reconstructed images transect the 3D mesh tetrahedral elements at oblique angles, and as a result, the breast structures in Figure 2.34(a) – (d) appear more blurred compared to those in Figure 2.34 34(e). 67 In order to evaluate the differences between the 3D reconstructed images using the data selection cases in Figure 2.34,, the recovered dielectric properties of the multi-plane plane (b and c) and in-plane (d) cases were subtracted from the full-data case in (a). The difference images are plotted in Figure 2.35 (a) FD – 4CS (b) FD – 2CS 68 (c) FD – IP Figure 2.35. Relative permittivity (top) and conductivity (bottom) difference images with respect to the recovered dielectric property values using the full-data (FD) set: (a) full-data data – multi-plane data sets with four consecutive planes (i.e. FD – 4CS) (b) full-data – multi-plane lane data sets with two consecutive planes (i.e. FD – 2CS), and (c) full-data – in-plane data sets (i.e. FD – IP) The amount of measured data decreases from the full-data set to multi-plane plane and then to inplane data selections.. Accordingly, the variations in the difference images increase from case (a) to (c) in Figure 2.35.. However, most of the variations appear on the exterior part of the reconstructed images and thee data selection scheme seems to have no major effects on the recovered dielectric properties of the breast. In order to quantify the variations, the relative permittivity (top) and conductivity (bottom) difference values in Figure 2.35 were extracted and they are plotted in Figure 2.36 36 as a function of the reconstruction nodes. In addition, for 69 statistical analysis purposes, the mean of the difference values, along with ± 2 standard deviations are also plotted in Figure 2.36. (a) FD – 4CS (b) FD – 2CS 70 (c) FD – IP Figure 2.36. Relative permittivity (top) and conductivity (bottom) difference values in Figure 2.35 as a function of the reconstruction nodes, along with the mean of the difference values ± 2 standard deviations (i.e. 95% of all difference data): (a) full-data – multi-plane data sets with four consecutive planes (i.e. FD – 4CS) (b) full-data – multi-plane data sets with two consecutive planes (i.e. FD – 2CS), and (c) full-data – in-plane data sets (i.e. FD – IP) Overall, all difference values are relatively close to zero, which indicates that there are minor variations in the 3D reconstructed property values using different data selection schemes. In case (c) where the largest deviations appear, 95% of all difference data points fall in a range of ± 2.5 units for the relative permittivity and ± 0.1 S/m for the conductivity profiles. These 71 variations may increase if less regularization is used during the reconstruction. For the present result, the Tikhonov regularization coefficient value of 50 was used. 2.6.2. Bone Imaging Microwave imaging can be used for detecting dielectric property changes in bone tissue with age or due to injury–related bone loss [81, 82]. Microwave tomography for heel imaging is particularly interesting, mainly because there is a significant dielectric property contrast between the bone and surrounding soft tissue. Moreover, like breast tissue, the heel is relatively easy accessible and it can fit into our microwave antenna array for imaging. In this section, we evaluate the feasibility of microwave heel imaging using a simple adaptation of our existing breast imaging system, as shown in Figure 2.37. (a) (b) Figure 2.37. Heel imaging in the microwave imaging system: (a) The heel extended through an aperture in the mounting plate, (b) The heel surrounded by the antenna array This was used to acquire four planes of microwave data (both in-plane and cross-plane, separated by 1 cm) of a normal heel (pt# 901, R) at 1300 MHz. The images were reconstructed 72 both in 2D and 3D on triangular triangular- and tetrahedral-element meshes, respectively.. The 2D circular reconstruction mesh was comprised of 559 nodes and 1044 triangular elements with a 6.9 cm radius, whereas the 3D cylindrical reconstruction mesh had a radiu radius of 6.9 cm, cm height of 5 cm, and consisted of 3132 nodes and 12707 tetrahedral ele elements respectively, as shown in Figure 2.38(a) and (b). (b) (a) Figure 2.38. 2D and 3D reconstruction meshes: (a) 6.9 cm radius circular mesh: 559 nodes and 1044 triangular elements, (b) 6.9 cm radius, 5 cm height cylindrical mesh: 3132 nodes and 12707 tetrahedral elements Figure 2.39 shows several views of the 3D reconstructed permittivity image of the heel with an iso-surface surface threshold of εr = 26, which captures the bone/soft tissue and the heel/coupling liquid interfaces.. The light blue shaded area surrounding the heel corresponds to the cylindrical reconstruction zone zone. The he underside of the surface exterior appears to conform nicely to the overall heel geometry. Moreover, several important features including, the large high-water-content content zones corresponding to the soft tissue nearer to the front of the foot, foot the 73 centrally located vertical calcane calcaneous us bone, and the narrowing of the heel towards the posterior of the foot, can be observed in the 3D reconstructed images in Figure 2.39. Figure 2.39. Several everal views of the 3D reconstructed permittivity image of the heel at 1300 MHz, with an iso-surface threshold of εr = 26 Figure 2.40(a) (a) and (b) show the corresponding 2D and horizontal slices of the 3D reconstructed permittivity (top) and conductivity (bottom) images, respectively. The four image pairs in Figure 2.40(a) (a) correspond only to the in-plane data acquired from the ankle (P1) down towards the heel apex (P4), whereas those in Figure 2.40(b) correspond to both in-plane and cross-plane data (two consecutive ive planes selection). In both cases, as the images progress towards the heel apex (P4), the outline of the heel shrinks until it completely disappears in the lower planes. Moreover, the calcaneous bone located at the center of the heel has lower properti properties es compared to the surrounding soft tissue, which 74 shows that the recovered property distributions are representative of the actual dielectric properties of the tissues. (a) 2D (b) 3D Figure 2.40. 1300 MHz,, (a) 2D and (b) slices of the 3D reconstructed permittivity (top) and conductivity (bottom) images at 4 planes of data separated by 1 cm,, starting from the ankle (P1) and moving down towards the heel apex (P4). 75 In terms of the size, location, feature shapes, and overall property distributions, the 3D reconstructed images are similar to their 2D counterparts. However, since the 2D slices of the 3D reconstructed images transect the 3D mesh tetrahedral elements at oblique angles, the heel structures in Figure 2.40(b) appear more blurred compared to those in Figure 2.40(a). Moreover, the calcaneous appears in the first two planes of the permittivity and the first three planes of the conductivity images in Figure 2.40(a), whereas in Figure 2.40(b), the bone is only visible in the first plane of the permittivity and the first two planes of the conductivity images. Based on the relative vertical position of the antennas to the heel during the exam, the recovered location of the calcaneous bone seems to be more accurate in the 3D reconstructed images. This observation confirms again that due to the data-model mismatch, the 2D reconstructed microwave images may show traces of the object in planes that extend beyond its actual locations. 76 3. Incorporation of a Priori Structural Information into the Microwave Image Reconstruction Algorithm 3.1. Introduction In prior studies, we have shown that MI has substantial potential for diagnosing breast lesions greater than 1cm in size [9, 17] – a result that was driven largely by the excellent specificity of the technique (rather than sensitivity as in the case of MR). We have also shown that MI can be a useful tool to monitor breast cancer treatment response to neoadjuvant chemotherapy [18, 83, 84]. However, MI suffers from a relatively poor spatial resolution compared to other conventional imaging modalities such as MRI or X-ray CT. Moreover, due to the use of the regularization method to stabilize the MI image reconstruction process, detecting boundaries of the regions of interest in the recovered microwave properties is a difficult task, and as a result, the diagnostic potential of MI is often limited. In response to these shortcomings, we have developed a method to combine the functional information available through MI (i.e. property contrast) with the high spatial resolution anatomical information available from other imaging modalities, such as MRI, to recover more accurate dielectric properties and potentially better characterize smaller masses [72]. This approach has been developed successfully in a related imaging system that utilizes a combination of NIR tomography and MRI methods for breast cancer detection [85-89]. A mathematical framework of microwave imaging reconstruction algorithm is presented in Section 3.2. Several regularization techniques, including the ones that incorporate the spatial priors, are also discussed in that section. In order to combine the functional information of the 77 MI with the spatial information from the MRI, we have developed a prototype MR-MI system, which is described in details in Section 3.3. Finally, the results of the 2D and 3D simulations, phantom experiments, and initial clinical data are presented in Sections 3.4 and 3.5, respectively. 78 3.2. A Mathematical Framework of Microwave Image Reconstruction Algorithm 3.2.1. Overview The reconstruction process in microwave imaging is based on determining the distribution of the constitutive parameters within tissue where the dielectric properties are embedded in the squared complex-valued wave number, which can be written as k 2 ( r ) = ω 2µ0ε ( r ) − jωµ 0σ ( r ) (3.1) where r is the position vector in the imaging domain, ω is the angular frequency, j is the imaginary unit, µ0 is the free-space permeability, ε is the permittivity, and σ is the conductivity [90]. The reconstruction algorithm mainly consists of solving two problems: the forward problem and the inverse (or optimization) problem. The forward problem involves computing the output (scattered field) from known inputs (microwave excitation) and system properties (dielectric property distribution of the tissue being imaged). In our algorithm, the calculation of the forward solution is based on Maxwell’s equations and it is computed using a finite difference time domain (FDTD) algorithm [79]. In the FDTD method, a frequency domain forward field response is produced for each transmitter and the individual field values are extracted at each receiver location [66]. The length of the vector k2 is N, the number of reconstruction parameters. The inverse problem, on the other hand, estimates the properties of an unknown volume (dielectric properties of the tissue) from a known input (microwave excitation) and measured field values. Since the inverse electromagnetic problem is non-linear, iterative Gauss-Newton schemes are well suited for the application [68, 69]. In our case, the solution of the inverse 79 problem is based on an improved Gauss-Newton (GN) iterative approach, called LevenbergMarquardt (LM), with a variance stabilizing transformation in which, the measured electric field vector Em is matched iteratively with the computed electric field vector Ec(k2) calculated using the forward model for a given distribution of the constitutive parameters stored in the vector k 2 [65, 91]. Although regularized Gauss-Newton methods can be susceptible to convergence to local minima [92, 93], we have found that log transformation serves to mitigate these effects [65, 94]. Indeed, recent studies of convergence indicate that our Gauss-Newton algorithm reaches the same image (i.e., solution) for widely different initial estimates suggesting that it is not easily trapped by local minima [95]. Completely different approaches to inverse problem solutions such as genetic algorithms and stochastic processes are also possible and have been applied successfully to MI as well [96-98]. In addition, due to the ill-posed nature of the inverse problem, some form of regularization method must be utilized in order to impose additional constraints [69, 70]. Regularization is often accomplished by introducing a priori information about the tissue being imaged, which may be necessary to ensure convergence of the reconstruction algorithm to the correct electromagnetic property distribution [71]. In the field of microwave imaging, a number of studies have investigated the incorporation of different types of priors ranging from the internal and/or external shape of the body to information about tissue dielectric properties including their upper and lower bounds [99-108]. For example, Crocco et al assumed that the object was homogeneous and that the permittivity value of the target was known [109]. In their shape reconstruction simulation study, El-Shenawee et al, used estimates of inclusion properties and location, and exact knowledge of the number of targets to aid convergence to a viable solution [110]. 80 In our algorithm, the Tikhonov regularization [69] is used to stabilize the reconstruction procedure, albeit with added smoothing. The objective function is 2 2 2 2 2 2 Ω = Γ m − Γ c (k 2 ) + Φ m − Φc (k 2 ) + λ L(k 2 − k02 ) m c m (3.2) c where Γ and Γ are the log magnitudes and Φ and Φ are the phases of the measured and computed field values, respectively [65, 66, 91, 94], λ is the weighting coefficient, also known as the Tikhonov regularization parameter, and L is a positive definite, dimensionless regularization 2 matrix. k02 is a prior estimate of k , and the two-norm V2 of a vector of length M (in this case M is the number of measurements) is the square root of the sum of the squares of the complex modulus of its elements: M v 2 = ∑v 2 (3.3) i i=1 In our previous and current studies, the choice of λ is derived empirically. After some manipulation, equation 3.2 can be solved for the iterative property update, ∆kη2 ∆ kη = J J + λ L L J T 2 T T −1 Γ m − Γc (kη2 ) T 2 2 − − L L ( k k ) m η 0 c 2 Φ − Φ (kη ) (3.4) where J is the Jacobian matrix, which has dimensions 2M × 2N, and it consists of derivatives of the log magnitude and phases of the computed field values with respect to the property values at each of the N reconstruction parameter mesh nodes. kη2 is the vector k2 at iteration η and is updated as ∆ kη2 = kη2+1 − kη2 (3.5) In addition, we use a dual-mesh approach where the forward solution is computed on a uniform FDTD grid (rectangular lattice in 2D and rectangular cuboid in 3D), while the 81 electromagnetic property parameters are reconstructed on a triangular (2D) or tetrahedral (3D) element mesh, placed concentrically within the antenna array [67]. In our original reconstruction algorithm, the regularization matrix L in equation 3.4 was set to the identity matrix, which applied the same weight to the values at all nodes within the imaging domain. In addition, for the right hand term in equation 3.4, k02 was set to kη2 as in the case for the Levenberg-Marquardt (LM) algorithm [68], leading to a simplified version of the update equation Γ m − Γ c ( kη2 ) ∆ kη = J J + λ I J m c 2 Φ − Φ ( kη ) 2 T −1 T (3.6) Throughout this chapter, The LM method, also referred as the “no priors”, corresponds to the update equation 3.6. 3.2.2. Soft prior Regularization In order to incorporate prior spatial information of the tissue (or phantom) being imaged into the microwave imaging, we have modified our reconstruction algorithm similarly to that reported in [85, 87-89, 111, 112]. In the new algorithm, called soft prior regularization, the spatial prior is considered “soft” because it does not force the property estimates inside an identified region to be constant. Instead, the known boundary data is used to adjust the regularization to smooth the property estimates within pre-identified regions, while limiting the smoothing across the boundaries (to preserve property changes at the interface with other tissues/regions). The basic idea behind the soft prior regularization is to more heavily weight the uniformity within regions that are assumed to have the same or similar dielectric properties. In addition, when two different regions share the same boundary, the smoothing across their common interface is penalized [111]. In the current implementation, we incorporate a priori 82 information about the anatomical structure of the tissue through the regularization matrix L in equation 3.4. According to known information about the tissue’s structure (derived from any high spatial resolution source, such as MR images), each node in the reconstruction mesh is labeled with an associated region number. Given two nodes, i and j, in the reconstruction mesh, with their associated matching regions, Ri and Rj, the corresponding entry in the L matrix is defined as: lij = 0 if Ri ≠ R j −1 N Ri if Ri = R j 1 if i = j (3.7) where is the number of nodes in region Ri – 1. Figure 3.1 illustrates how the regularization matrix L is constructed in a simplified case, where there are just two regions, R1 and R2, containing N1 + 1 and N2 + 1 nodes, respectively. Figure 3.1. Soft prior regularization matrix L with two distinct regions 83 T Based on the above construction of L, L L in equation 3.4 is an approximation of a second order Laplacian smoothing operator inside each region, which limits the smoothing across the boundary of distinct regions [113, 114]. Since the structure of the tissue being imaged does not change during the iterative image reconstruction algorithm, both the regularization matrix L and T Laplacian smoothing operator L L, can be calculated once and stored at the beginning of the reconstruction process. In this way, redundant calculations are avoided and the algorithm becomes more efficient. Effective encoding of spatial priors in microwave image reconstruction requires software tools for region identification and segmentation, in addition to mesh generation. We have used both commercial software (e.g. Mimics, SolidWorks) as well as our own algorithms to postprocess the structural information. In the current study, we apply a collection of user-selected options, such as region growing methods in Mimics, to segment regions of interest in the MR images. When prior spatial structures are available, a more specific and customized reconstruction parameter mesh is desired. In fact, sufficient nodes need to be deployed on the boundary of adjacent regions, so that the interfaces in the reconstructed image are smooth. We will use methods developed in-house for surface (for 2D) and volume (for 3D) mesh generation, which control the nodal sampling density to capture the irregular object geometry. Finally, the subsequent user-defined mesh needs to be co-registered with the acquired microwave data – This step usually involves shifting, rotating, and/or flipping the mesh. Due to the different coordinate systems, co-registration of the microwave and MR images can be a very challenging task. Therefore, we utilize additional markers during both imaging sessions (i.e. MI and MRI) to ensure accuracy in the co-registration process. 84 Further on in this chapter, we evaluate the performance of the MI with the soft prior regularization relative to no priors in a number of simulations, phantom experiments, and some preliminary patient data. The corresponding 2D and 3D results are presented in Sections 3.4 and 3.5, respectively. 3.2.3. Hard priors – Parameter Reduction Unlike the soft prior regularization technique, the hard priors approach uses a priori constraints through parameter reduction. The basic idea behind the hard priors is to reduce the size of the reconstruction parameters by assuming that there are homogeneous dielectric property distributions within each of the pre-identified regions [88]. As a result, given structural information with R distinct regions, single values for permittivity and conductivity are reconstructed within each region. The region identification process is identical to that in the soft priors. However, the parameter reduction is implemented during the property update in equation 3.6, by defining a new Jacobian matrix: Jɶ = JKɶ (3.8) , and K is the a priori parameter reduction matrix defined as: where 1 kij = 0 if node(i ) ∈ R j (3.9) if node(i ) ∉ R j where node(i) is the ith node in the parameter reconstruction mesh, and Rj is jth region identified from the prior structural information. The dimensions of K are the number of nodes by the number of regions (N×R), and since the dimensions of the original Jacobian matrix J are two times the number of measurements by two times the number of reconstruction nodes (2M × 2N), has the dimensions of 2M × R. In the new sensitivity matrix, all the columns corresponding to 85 the same region are summed, and as a result, the number of reconstruction parameters drops from N (number of nodes in the reconstruction mesh) to R (number of pre-identified regions of interest). Therefore, no matter how many nodes are placed in each region, the number of unknowns is set to the number of regions in the prior structural data R. As long as R < M, no regularization is required during the image reconstruction, and consequently equation 3.6 can be simplified to Γ m − Γ c ( kη2 ) T ɶ −1 ɶ T ɶ ∆ kη = J J J m c 2 Φ − Φ ( kη ) 2 (3.10) We compare the performance of the hard priors with the soft priors as well as no priors in a simulation study case in Section 3.5.1.7. 3.2.4. Error Analysis Throughout the remainder of this chapter, we will use an error analysis to compare quantitatively the soft prior and no prior results. Since the true values of the dielectric property distributions are known in simulations and phantom experiments, we can calculate the relative error between the true/independently-measured property distributions and the reconstructed values as: N Errw = ∑ wn n=1 recon exact V(n) −V(n) (3.11) exact V(n) where V(n)recon is the reconstructed dielectric property value (either permittivity or conductivity) at node n in the reconstruction mesh, whereas V(n)exact is the true value of the selected dielectric property at that location. To account for the fact that nodes in the reconstruction meshes may not 86 be uniformly distributed, a weighting factor is added at each location computed as wn = An /A, where An is the area of the elements surrounding node n and A is the total imaging area. 87 3.3. Integration of Microwave Tomographic Imaging System with MRI 3.3.1. Motivation Unlike microwave tomographic imaging, MRI has an excellent spatial resolution (down to 0.5 mm) and it is broadly used in clinical practice. In fact, MRI is the most sensitive medical exam for breast cancer surveillance, especially in the dense breast [115]. The high spatial resolution of internal anatomical structures that MRI provides is based on the natural tissue contrast in terms of the relaxation times [116]. Breast MRI has been used successfully to monitor treatment response to neoadjuvant chemotherapy and is deployed routinely for cancer staging and surgical treatment planning [117-119]. Despite the high cost of the imaging system and the relative expense of an exam, the number of clinical breast MR imaging studies is likely to continue to grow in the future, mainly due to the extraordinarily detailed 3D anatomical information that can be obtained [120, 121]. Nonetheless, breast MRI has substantial limitations, including a high false positive rate and the fact that it cannot provide information about a tissue’s physiological state [121]. As a result, MRI has been combined with other imaging modalities to improve the specificity of diagnosis and therapy [122-125]. For instance, in the combined NIR-MRI imaging, MRI contributes excellent spatial resolution, while the NIR tomography provides the complementary functional information about the tissue [126]. Similar to the NIR-MRI concept mentioned above, we proposed to combine the high spatial resolution of MRI with the high specificity of the microwave dielectric properties. One of the most apparent challenges for combining microwave tomography with MRI is the physical conflict between the small MR bore size and our requirement to utilize a dielectric coupling bath. Based on our experience from earlier studies, we have been able to strategically reduce the size 88 of the imaging tank enough so that it can accommodate many women even within an MR bore. In addition, the non-ferrous monopole antennas and their unique low profile design make them ideal for operating in the MR environment while inducing only minimal and acceptable MR image disruptions. Along with the microwave antennas affecting the MR signals, the RF gradients used in the MR gradient sequences produce high power signal bursts that can saturate and desensitize receivers and subsequently disrupt the desired electromagnetic signals. In the following four Sections, 3.3.2 – 3.3.5, we will discuss these issues in details and present strategies to overcome barriers in developing a multi-modality MI-MR imaging system with minimal unwanted effects. 3.3.2. Integration of MI System into MRI As with all research performed in an MR environment, the electronics are generally positioned outside of the examination area. In our case, this is achieved by using long, flexible coaxial cables which connect our imaging tank in the bore to a set of 16 SMA flange connectors positioned on a bulkhead panel between the exam room and the external control room. Figure 3.2(a) shows our clinical microwave imaging system located just outside and below the connector panel with short coaxial cables running from the connectors down to the transceiver modules inside the bed. 89 (a) (b) Figure 3.2. MI-MR imaging: (a) the portable data acquisition system in the control room with the auxiliary cables attached to the bulkhead connectors in the wall, (b) the imaging tank outside of the MR bore showing the semi-rigid cables extending from the tank While the extended cable lengths (7.9 dB return path attenuation at 1.5 GHz for the 6.5 m cables – LMR200 from Times Microwave Systems, Wallingford, CT) adds significant signal loss, it does allow us to utilize the MR–based tank while locating all electronics out of the MR room. Since signal losses in the tank increase with frequency, for future development, we intend to implement lower loss cables (LMR600 – 2.0 dB loss at 1.5 GHz for the return path loss over 6.5 m long cables), which will ultimately allow us to operate at higher frequencies. An overriding principle for the MI-MRI multi–modality integration is to minimize the amount of metal near the MR bore, especially removing anything that may be ferromagnetic. As can be seen in Figure 3.2(b), the cables leading directly away from the illumination tank are 3.58 mm diameter, semi–rigid coaxial cables made from a copper outer conductor and a silver–plated, copper center conductor. Once these have extended at least 1 meter from the tank, the cables are connected to the flexible coaxial lines described above, which are more convenient for the long 90 span across the room to the bulkhead interface connectors. The microwave antenna array configuration in this new system is identical to that of our original clinical imaging system, and each antenna is capable of operating as both transmitter and receiver in the frequency range of 700 MHz to 3 GHz. However, in our actual prototype system, the antennas are not movable and do not have the capability of collecting data at multiple planes, and as a result, only 240 measurements (16 transmitters × 15 receivers) can be collected for each image. 3.3.3. MR Image Artifact Reduction Strategies The most fundamental challenge in combining microwave imaging with an MR system is the need to bring metallic objects into the bore. Ferromagnetic materials are clearly prohibited, but many other types of metals such as titanium, aluminum, copper and silver are safe to position within the magnet. However, these metals can still disrupt the images. Figure 3.3 shows a 19 cm height, 13 cm radius microwave imaging tank used for a phantom experiment in the MR. The phantom was comprised of a 5.4 cm radius cylindrical container (84:16 glycerin:water mixture) with 2.1 cm and 1.0 cm radius cylindrical inclusions filled with 72% and 55% glycerin mixtures, respectively. The associated copper–clad, semi–rigid coaxial feed lines and the stainless steel flange connectors were placed on the under surface of the tank (the coaxial cables consisted of a 1.16 mm diameter center conductor, a 3.05 mm insulator diameter and 3.58 mm outer conductor diameter). 91 (a) (b) Figure 3.3. MI tall tank with a large cylinder breast phantom (green) and smaller inclusions surrounded by the monopole antenna array: (a) side view and (d) top-down view Figure 3.4 shows the corresponding T2–weighted MR images of the phantom experiment. Both Figure 3.4(a) and (b) are for silver–plated copper center conductors. The first imaging plane in (a) transects the monopole antennas, while the one in (b) intersects the coaxial feedlines 3 cm below the antenna/coaxial line interface plane. For the imaging plane in the active part of the antennas in Figure 3.4(a), there are minor artifacts surrounding each antenna, but these perturbations do not impact the images of the phantom. In fact, for co-registration purposes, the positions of the antennas can be accurately marked as a reference with respect to the microwave imaging array geometry. At the plane transecting the feedline in (b), the artifacts surrounding the antenna array have increased and there are minor perturbations to the shapes of the interrogated phantom and inclusions. For Figure 3.4(c), each monopole antenna and coaxial line used a steel center conductor, and the phantom only had the 2 cm diameter inclusion. The MR image in (c) is completely disrupted due to the presence of ferromagnetic materials in the magnet. This sequence of images emphasizes that the ferromagnetic materials cannot be used in the combined 92 MI-MRI system. In contrast, it also demonstrates that other metals can still disrupt the image quality, and that reducing the metal footprint clearly reduces artifacts. In this regard, the monopoles are optimal in that they employ the least amount of metal possible of any antenna. (a) (b) (c) Figure 3.4. T2–weighted MR images of the phantom experiment: (a) at the center plane of the active part of the antennas and (b) at a plane 3 cm below the interface between the semi-rigid coaxial line and the antenna. For (a) and (b) the phantom included a 5.4 cm radius outer cylinder and the 2.1 and 1.0 cm radius inclusions, and the coaxial line center conductors were silverplated copper. (c) Image at the center plane of the active part of the antennas for the 5.4 cm radius outer cylinder and just the 1.0 cm radius inclusion. In this case the coaxial cable center conductor was steel. 3.3.4. Multipath Signal Assessment Multipath signals are problematic for any communication or radar system. They are especially troublesome for near field systems where there is less opportunity for natural attenuation from simple physical isolation of the transmitters and receivers [127, 128]. One of the principle ways we have overcome this issue with our system has been the adoption of a lossy coupling bath. 93 In earlier work, we determined that the tank side walls could be placed as close as 7 cm from the antennas without signal reflections from the wall interfering with desired transmitting signals [129]. In addition, we found that we could operate with the tips of the monopole antennas as close as 2 mm from the top air/liquid interface, while only inducing minor artifacts in the microwave images [29]. However, our earliest systems outside of the MR used relatively deep tanks (usually 24 to 30 cm deep), which worked well and did not appear to create unwanted multipath signals. As we proceeded to translate this technology to integration with the MR, space was a premium and we initially reduced the internal tank height towards 8.0 cm with shorter, straight antenna feedlines shown in Figure 3.5(a). However, the reduced-height configuration induced multiple unwanted signals to our microwave data acquisition system. After performing a thorough analysis, we concluded that the offending multipath signals were the result of surface waves propagating along the outside coaxial feedlines and along the tank base/liquid interface [130]. Given that the majority of the previous loss derived primarily from the external coaxial feed path, we determined that it was essential to maintain a minimum feed length of 9.0 cm for sufficient attenuation. Utilizing straight feedlines, this would have required a minimum internal tank height of 11.7 cm (including 2.5 cm for the length of the monopole antenna and 2 mm between the antenna tips and liquid/air interface). However, by precision bending the semi-rigid coaxial lines into a serpentine shape as shown in Figure 3.5(b), we were able to retain the required shorter, 8.0 cm height MR-based tank. 94 (a) (b) Figure 3.5. Monopole antenna array: (a) in the MI shortened tank with a breast phantom (green), and straight 5.0 cm long feedlines (short), (b) with 9.0 cm long serpentine feedlines (tank wall removed) 3.3.5. MR RF Gradient Signal Disruption of Microwave Signals The RF gradient pulses from the MR coils are primarily at the Larmor frequency (128 MHz for the 3T system) and radiate freely into the MR bore. However, through our own internal measurements using the same configuration as in Section 3.3.3, except for replacing the custom data acquisition system with an Agilent E4402B Spectrum Analyzer (SA), we detected quite high signal levels for the Larmor signal and its harmonics. For these measurements, the signals were measured after being coupled to the monopole antennas submerged in the 80:20 glycerin:water coupling bath which then transit through the coaxial cables and through the bulkhead to the SA. In this case the MR-based signals had to propagate through the bath and couple to the antennas – the coupling bath attenuation decreases rapidly at lower frequencies, especially as low as the Larmor fundamental frequency. In addition, the coupling efficiency to 95 the antenna is suboptimal, especially given that the polarization of the signal is uncertain. However, this is the same configuration for when the data acquisition system is attached to the channel cables and is representative of the experimental conditions for both phantom tests and patient exams. The power levels measured for the Larmor fundamental and the first six harmonics (128, 255, 383, 510, 638, 766, and 893 MHz, respectively) were nominally 0, – 40, –50, – 67, –72, – 77, and –79 dBm, respectively. Because of the intermittent nature of the pulses, these power levels were estimates based on visual observations of the SA during continuous operation of the MR system. The signal strength roll off demonstrates the normal decrease expected for increasing harmonic numbers. In addition, the wave attenuation increases substantially in the glycerin bath with frequency, which also adds to the signal degradation. Other losses included the coaxial cable attenuation as mentioned in Section 3.3.3, which decreased somewhat at lower frequencies. Using our custom data acquisition system, Figure 3.6(a) shows plots of the raw measured signal amplitudes (uncalibrated) at 1300 MHz for 15 receivers collecting data from all 16 transmitters in a homogeneous bath due to the 1 mW (0 dBm) transmitted signals. Most of measured data exhibit the typical pattern for received signals, i.e. a symmetric parabolic plot. However, the power levels for the four corrupted sets (transmitters 2, 5, 13, and 16) are typically 20 to 40 dB lower than for the other receivers. Because the MR RF gradients are intermittent and not synchronized with the microwave data acquisition, the distortions did not affect all data sets. Had the signal distortion been from simple constructive or destructive addition of one of the Larmor harmonics with the desired signals, we would have expected positive and negative distortions around the predicted pattern with the disruptions visible mostly at the lower signal 96 levels. Instead, the distortions are most likely the result of amplifier de-sentization, whereby an unwanted signal is processed simultaneously with the desired one and has a substantially higher amplitude, such that it effectively turns off the small signal amplification process for all other frequencies. (a) No low pass filters (b) With two low pass filters Figure 3.6. Plots of the 1300 MHz amplitudes as a function of receiver number for all 16 transmitting antennas when the imaging tank is operating in the MR bore with MR data being 97 acquired simultaneously for the cases where (a) no low pass filters, and (b) two low pass filters installed, respectively. Given that the Larmor fundamental is roughly 0 dBm after being received by an antenna, after the first stage of 20 dB receiver amplification it would have easily saturated that amplifier and even more so in the subsequent cascaded element. It was more than apparent that this was the offending mechanism. After installing two VHF-740 Mini-Circuits high pass filters (combined attenuation of 122 dB at 128 MHz), the unwanted distortions completely disappeared, as shown in Figure 3.6(b). 98 3.4. Results – 2D In this section, we evaluate the performance of the MI with the soft prior regularization relative to no priors in 2D simulations, phantom experiments, and some preliminary patient data where the boundaries of different regions of interest were known. 3.4.1. Simulation Experiments 2D simulated measurement data was generated by the hybrid boundary element/finite element approach [63] for different target shapes. The images were reconstructed in 2D using the update equation (3.6) with no spatial information and the update equation (3.4) with the soft prior regularization algorithm described in Section 3.2.2. 3.4.1.1. Mesh Resolution and Noise Level Figure 3.7 illustrates a circular inclusion centered at (x, y) = (0, –3 cm) with a radius of 1.4 cm and dielectric properties of εr,Inc = 51.16 and σInc = 1.44 S/m embedded in a background medium with dielectric properties of εr,BK = 15.60 and σBK = 0.90 S/m. Figure 3.8 shows different parameter reconstruction meshes used for this experiment. Noise ranging from -110 dBm to -80 dBm was added to the synthetic measurements. Figure 3.9 shows reconstructed images at 1300 MHz (noise level of -100 dBm) using the no priors (on 473, 961, and 915 node meshes) and soft prior (on the 915 node mesh) regularizations, respectively. 99 Figure 3.7. Schematic of the imaging domains evaluated. The background diameter was 14 cm (antennas are positioned on a 15.2 cm diameter). The circular inclusion’s center with radius r = 1.40 cm was located at (0, –3 cm). (a) (b) (c) Figure 3.8. Reconstruction parameter meshes (a) uniformly distributed 473 node mesh, (b) uniformly distributed 961 node mesh, and (c) 915 node mesh with preferential node deployment in the inclusion in Figure 3.7 100 (a) (b) (c) (d) Figure 3.9. Simulated 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images (a) without priors on thee 473 node mesh, (b) without priors on the 961 node mesh, (c) without priors on the 915 node mesh with preferential node deployment in the inclusion, and (d) with the soft prior regularization rization on the 915 node mesh mesh. Without spatial priors, the reconstructed images are very similar, which indicates that increasing the number of nodes or changing their distribution to be preferentially greater within the inclusion in the reconstruction mesh does not improve the quality of the recovered images. Both regularizations ions recovered the inclusion in the permittivity and conductivity images, but those from the soft prior technique are more accurate in terms of the location of the inclusion and its property values. In addition, the recovered background permittivity and co conductivity nductivity values are more uniform in the soft prior case. Figure 3.9 confirms that these improvements are not due to the reconstruction mesh, but to the regularizat regularization matrix, which incorporates structural priors into the reconstruction procedure. Therefore, from this point on, only the 473 node mesh will be 101 used for the 2D reconstructions with no priors. The no priors weighted permittivity and conductivity errors (described in 3.2.4) were 0.172 and 0.128, respectively, while those with the soft prior regularization were reduced by 83% to 0.028 for εr and by 87% to 0.016 for σ. In order to compare the reconstructed and the true dielectric properties, vertical transects of the permittivity and conductivity profiles along the y-axis are shown in Figure 3.10 for added noise levels of -110, -100, -90, and -80 dBm. (a) –110 dBm (b) –100 dBm 102 (c) –90 dBm (d) –80 dBm Figure 3.10. Comparison of the 1300 MHz reconstructed permittivity (top) and conductivity (bottom) values using the no priors and soft prior regularization with different levels of added noise: (a) –110, (b) –100, (c) –90, (d) –80 dBm As expected, artifacts increased in both the permittivity and conductivity images without priors as the noise level rose, especially in the -80 dBm conductivity images where the fluctuations are significant. The soft prior regularization tolerates the added noise significantly better with relatively minor decreases in the recovered inclusion permittivity and only slightly greater reductions in conductivity at the highest noise level. Using the soft priors, the reconstructed permittivity values were underestimated (~10-15%) in the inclusion at all noise levels. Notwithstanding, the method successfully characterized the inclusion given the large 103 property contrast with the background. In addition, even when considerable noise was added to the measured data (-90 dBm), the algorithm with the soft prior regularization recovered the inclusion conductivity very accurately, and started to underestimate (~15%) its property values only at higher noise levels (-80 dBm). Based on the results shown above, -100 dBm noise is the level where we can see a degradation effect on the reconstructed images. In practice, the noise floor of our system is closer to -135 dBm, and as a result, we generally have a very large signal to noise ratio (SNR) even for the weakest signals transmitted from the opposing antennas [31]. 3.4.1.2. Estimation of the Parameter Coupling In order to study the effects of the no prior and the soft prior regularization on reconstruction parameter coupling, two simulation experiments with the same setup and the same background medium (εr, BK = 15.60 and σBK = 0.90 S/m) as the previous experiment (Figure 3.7) were performed. In the first case, no permittivity contrast existed in the inclusion region (dielectric properties were εr, Inc = 15.60 and σInc = 1.44 S/m). The second case had no conductivity contrast in the inclusion (dielectric properties were εr, Inc = 51.16 and σInc = 0.90 S/m). Figure 3.11 shows transects of the 1300 MHz reconstructed permittivity (top) and conductivity (bottom) profiles along the y-axis for the first and second experiments, using the 473 and 915 node meshes, respectively (with -100 dBm of added noise). 104 (a) (b) Figure 3.11. Comparison of 1300 MHz reconstructed permittivity (top row) and conductivity (bottom row) profiles for the no priors and the soft prior regularization: (a) no permittivity contrast (left column), (b) no conductivity contrast (right column). While both methods (no priors and the soft priors) handle the no-permittivity contrast case (a) effectively, the soft prior regularization is clearly superior when no conductivity contrast exists in the inclusion (b). The soft prior regularization profile almost overlaps the exact solution, whereas the no priors curve has significant differences with the exact solution. Quantitatively, the permittivity – conductivity parameter coupling is only about 5% with the soft prior regularization in both cases in Figure 3.11. 105 3.4.1.3. Arbitrarily-Shaped Target An arbitrarily-shaped target inclusion, as shown in Figure 3.12(a), with dielectric properties of εr, Inc = 40.0 and properties of σInc = 1.30 S/m was embedded in a coupling medium with dielectric εr, BK = 15.60 and σBK = 0.90 S/m. Based on the exact location of the target, the customized mesh in Figure 3.12(b) was generated and used for the soft prior regularization algorithm. (a) (b) Figure 3.12. (a) Schematic of the imaging domains evaluated. The arbitrarily shaped inclusion was located in the upper part of the imaging domain. The reconstructed values of the simulation experiment with arbitrarily shaped inclusion were extracted at 30 points evenly distributed along the line x = –2 cm. (b) Customized soft prior reconstruction mesh comprised of 1725 nodes and 3248 triangular elements. 106 Figure 3.13 shows the 1300 MHz reconstructed images with -100 dBm added noise both with and without priors, on the 1725 ((Figure 3.12(b))) and 473 node meshes, respectively. (a) (b) Figure 3.13. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images from a phantom tom experiment with arbitrarily arbitrarily-shaped target inclusion for (a) no priorss on the 473 node mesh, and (b) soft priors on the customized 1725 node mesh in Figure 3.12(b) While the no prior permittivit permittivity image shows the target at roughly the correct location, the recovered target is considerably shifted towards the top part of the counterpart conductivity image in Figure 3.13(a), (a), and the level of artifacts is significantly higher. In contrast, the reconstructed soft prior images not only recover the complex shape of the target and its exact location, but they also appear to recover the dielectric property distributio distributions ns more accurately. accurately In order to confirm this observation, transects of the reconstructed onstructed permittivity and conductivity profiles in Figure 3.13 were extracted along the line x = -2 cm,, and they are plotted along with the exact solution in Figure 3.14(a) (a) and (b), respectively respectively. 107 (a) (b) Figure 3.14. Comparison of the 1300 MHz reconstructed (a) permittivity and (b) conductivity profiles (along the line x = –2 cm cm) in Figure 3.13 phantom tom experiment with arbitrarily-shaped arbitrarily target inclusion for the no priors and soft prior regularization The transected profiles in Figure 3.14 confirm that the soft prior recovered property distributions are more accurate. In fact, tthe no priors weighted permittivity and conductivity errors (described in 3.2.4) were 0. 0.398 and 0.411,, respectively, while those with the soft prior regularization were ere reduced by 94% to 0.022 for εr and by 87% to 0.051 for σ. 3.4.1.4. Multiple Inclusions and Number of R Regions In order to study the effects of the number of inclusions, as well as the number of nodes and regions in the reconstruction mesh on the recovered property values when the soft prior regularization is used, three simulation experiments with one and three inclusions were performed. In all cases, -100 dBm noise was added to data. a) Single Inclusion: In the first experiment, a single inclusion of 1.5 cm radius centered at (x, y) = (3, 0.0 cm) with dielectric properties of εr, I1 = 40.0 and σI1 = 2.0 S/m was placed in a 7 108 cm radius circular region (mimicking an 80:20 mixture of glycerin:water background medium) located at the center of the imaging domain (x = y = 0) with dielectric properties of εr, BK = 22.4 and σBK = 1.23 S/m. Two customized reconstruction meshes with 402 and 1203 nodes, as shown in Figure 3.15(a) and (b), respectively, were then used to reconstruct the dielectric property values of the inclusion and the background region at 1300 MHz. (a) (b) (c) (d) (e) (f) Figure 3.15. Customized soft prior reconstruction meshes used in 1 and 3–inclusion simulation experiments: (a) 402 node mesh with one inclusion and two distinct regions (RBK and RI1), (b) 1203 node mesh with one inclusion and two distinct regions (RBK and RI1), (c) 407 node mesh with three inclusions and two distinct regions (RBK and RI1), (d) 1203 node mesh with three 109 inclusions and two distinct regions (RBK and RI1), (e) 407 node mesh with three inclusions and four distinct regions (RBK, RI1, RI2, and RI3), and (f) 1203 node mesh with three inclusions and four distinct regions (RBK, RInc1, RInc2, and RInc3) Since the recovered properties within each region were very similar, the maximum reconstructed values for both meshes along with the exact properties in each region were selected and are displayed in Figure 3.16.. (a) (b) Figure 3.16. Soft prior reconstructed (a) permittivity and (b) conductivity values using two meshes with different number of nodes, along with the exact properties of background and inclusion at 1300 MHz Both permittivity and conductivity values of the background medium are reconstructed reco very closely to the exact properties. Moreover, the inclusion is successfully characterized with both meshes. However, the recovered property values of the target-inclusion inclusion, especially the permittivity values, are closer to exact exact, when the finer mesh (1203 nodes) is used for reconstruction. b) Multiple Inclusions with Identical Dielectric Properties: Similar to the first experiment, in the second experiment a circular background region with a radius of 7 cm was 110 placed at the center of the imaging domain. In this case, three target inclusions (I1, I2, and I3) with the same radius of 1.5 cm were centered at (x, y)I1 = (3, 0 cm), (x, y)I2 = (–2, 3 cm), and (x, y)I3 = (–2, –3 cm), respectively. The dielectric properties of all inclusions were set to εr, I1 = εr, I2 = εr, I3 = 40.0 and σI1 = σI2 = σI3 = 2.0 S/m; however, two different cases were considered for the reconstruction: b1) All reconstruction nodes in the three inclusions were assigned the same region number (i.e. RI1 = RI2 = RI3), and as a result, only two distinct regions were considered: one for the background (RBK) and one for all inclusions (RI). In this case, two customized meshes with different numbers of nodes (407 and 1203 nodes) were used for the reconstruction, as shown in Figure 3.15(c) and (d), respectively. b2) The reconstruction nodes in each inclusion were assigned a different region number, assuming different inclusions and permitting their property distributions to be recovered independent of each other. In this case, four distinct regions were considered for the reconstruction: one for the background (RBK) and three for inclusions (RI1, RI2, and RI3). Similar to (b2), two customized meshes with different number of nodes (407 and 1203 nodes) but with four distinct regions were used for the reconstruction, as shown in Figure 3.15(e) and (f), respectively. Figure 3.17 shows the maximum soft prior reconstructed (a) permittivity and (b) conductivity values of each region, along with their exact property values at 1300 MHz, using the customized reconstruction meshes in Figure 3.15(c) – (f) 111 (a) (b) Figure 3.17. Maximum soft oft prior recovered (a) permittivity and (b) conductivity values along with the exact properties of the background and target inclusions at 1300 MHz using the customized reconstruction struction meshes in Figure 3.15(c) – (f) Similar to the results of the first experiment with only one inclusion, both permittivity and conductivity values of the background region are recovered very accurately, and all inclusions are successfully characterized.. When two-region meshes (Figure 3.15(c) (c) and (d)) (d) are used, the recovered values of all targets are the same and they are all very close to the exact property distributions. However, when four four-region meshes (Figure 3.15(e) and (f))) are used, some variations appear between the recove recovered properties of inclusions, though all inclusions were 112 assigned the same property values. In these cases, the deviation from the exact solution is more prominent in the recovered permittivity values. A complete summary of the relative dielectric property errors using different reconstruction meshes with two and four distinct regions is presented in Table 3.1. The negative sign indicates that the reconstructed property value is underestimated (i.e. less than the exact property value). Table 3.1. Relative (a) permittivity and (b) conductivity errors for the multiple inclusions simulation experiment using different reconstruction meshes with two and four distinct regions Rel. Err. Rel. Err. Rel. Err. Rel. Err. εr, BK εr, I1 εr, I2 εr, I3 Mesh 407 – 2r Mesh 407 – 4r Mesh 1203 – 2r -0.4% -0.4% 0.0% -11.5% -10.8% -7.5% -11.5% -11.3% -7.5% -11.5% -12.0% -7.5% Mesh 1203 – 4r 0.0% -7.5% -8.0% -8.2% Rel. Err. Rel. Err. Rel. Err. Rel. Err. σ BK σ I1 σ I2 σ I3 0.0% 0.0% 0.8% 0.0% -2.0% -5.0% 1.0% -2.5% -2.0% -4.0% 1.0% -0.5% -2.0% 2.0% 1.0% 5.0% Reconstruction Mesh (a) Reconstruction Mesh Mesh 407 – 2r Mesh 407 – 4r Mesh 1203 – 2r Mesh 1203 – 4r (b) When a finer soft prior mesh (1203 node) is used, the reconstructed properties are overall closer to the exact solution, and both relative permittivity and conductivity errors decrease in all regions. Moreover, restricting variations between different regions that have the same properties (i.e. using meshes with two distinct regions instead of four) can improve the accuracy of the soft prior reconstructed dielectric properties. 113 c) Multiple Inclusions with Different Dielectric Properties: This experiment was identical to the second experiment in (b), except for the property values of the target inclusions. In this case, dielectric properties of S/m, and εr, I1 = 27.0 and σI1 = 1.45 S/m, εr, I2 = 35.0 and σI2 = 1.66 εr, I3 = 40.0 and σI3 = 2.0 S/m were assigned to I1, I2, and I3, respectively. Since dielectric properties of the inclusions were different, only customized meshes with four distinct regions (BK, I1, I2, and I3) shown in Figure 3.15(e) and (f) were used for the soft prior reconstructions. Figure 3.18 shows the maximum soft prior recovered dielectric property values of all regions along with the exact property values at 1300 MHz. 114 (a) (b) Figure 3.18. Maximum soft oft prior recovered (a) permittivity and (b) conductivity values along with the exact properties of the background and target inclusions at 1300 MHz using the customized reconstruction meshes in Figure 3.15(e) and (f) Similar to the results of the first and second experiments, both permittivity and conductivity values of the background region are recovered very accurately, and all target inclusions are successfully characterized characterized. Moreover, the correct pattern of property increase from I1 to I3 can be seen in the soft prior reconstructed permittivity and conductivity profiles. A complete error analysis of the recovered soft prior property values is presented in Table 3.2. The negative egative sign indicates that the reconstructed pproperty roperty value is underestimated (i.e. less than the exact property value). 115 Table 3.2. Relative (a) permittivity and (b) conductivity errors for the multiple inclusion simulation experiment with different property contrasts using the soft prior reconstruction meshes with four distinct regions in Figure 3.15(e) and (f) Reconstruction Mesh Mesh 407 – 4r Mesh 1203 – 4r Rel. Err. Rel. Err. Rel. Err. Rel. Err. εr, BK εr, I1 εr, I2 εr, I3 -0.9% -0.4% -3.7% -3.0% -8.3% -5.4% -12.5% -7.8% (a) Rel. Err. Rel. Err. σ BK σ I1 σ I2 σ I3 Mesh 407 – 4r 1.6% 0.0% -3.6% -4.0% Mesh 1203 – 4r 1.6% 0.0% -1.8% -1.0% Reconstruction Mesh Rel. Err. Rel. Err. (b) Similar to the results of cases (a) and (b), when a finer soft prior mesh (1203 node) is used, the relative permittivity and conductivity errors decrease. This effect is more prominent on larger relative errors as the dielectric property contrast between the background medium and target inclusions increases from I1 to I3. In the first and second experiments in cases (a) and (b), the true dielectric properties of I1 centered at (x, y)I1 = (– 3, 0 cm) were identical and equal to εr, I1 = 40.0 and σ I1 = 2.0 S/m. Therefore, in order to study the effects of additional regions of interest (in this case I2 and I3) on the recovered dielectric properties of I1, the maximum recovered soft prior property values of I1 in the first and second simulation experiments (on 1203 node mesh with one and three target inclusions shown in Figure 3.15(b) and (f), respectively) are plotted in Figure 3.19. 116 (a) (b) Figure 3.19. Comparison between the soft prior reconstructed (a) permittivity and (b) conductivity values of I1 with and without additional targets (three inclusions and one inclusion cases, respectively) The difference in recovered dielectric property values in Figure 3.19 is minimal (2.75% (2.75 in permittivity and 2.5% in conductivity conductivity), which suggests that additional targets do not have a major effect on the 2D reconstructed soft prior dielectric property profiles. 3.4.1.5. Sensitivity to a False Region egion of Interest The behavior of the soft prior regularization in regions which are identified prior to the MI property estimation, but do not actually have contrast, is of interest because of the potential for an approach to identifying false inclusions based solely on the presumed structural information. In order to study how the dielectric property contrast between the background and a target zone would affect the sensitivity of the soft prior regularization algorithm to a false region of interest, three simulation experiments were performed. Dielectric properties of the background coupling medium and the actual target inclusion in th thee experiments were set to those of the 86:14, 40:60, 55:45, and 70:30 mixtures of glycerin and water, respectively. The corresponding property values are notes in Table 3.3. 117 Table 3.3. Dielectric properties of different mixtures of glycerin and water used for the simulation experiments with a false region of interest Glycerin:water Mixture (case) Relative Permittivity Conductivity (S/m) 40:60 Inc – Case (a) 55:45 Inc – Case (b) 70:30 Inc – Case (c) 86:14 BK – (all cases) 62.6 51.2 40 15.6 1.12 1.44 1.64 0.90 Figure 3.20 shows the customized reconstruction mesh with a false inclusion region (in red) of 1.4 cm radius centered at (0, 3 cm) along with the actual circular target zone (in green) centered at (0, – 3 cm) with a radius of 1.4 cm. 118 Figure 3.20. Customized reconstruction mesh comprised of 1196 nodes and 2215 triangular elements, accounting for two inclusion regions: a false circular region of interest (in red) with radius of 1.4 cm centered at (0, 3 cm), and the actual circular target zone (in green) with radius of 1.4 cm centered at (0, -3 cm) In all three cases ((a), (b), and (c) in Table 3.3, -100 dBm of noise was added to the synthetic measurements and images were reconstructed at 1300 MHz using the soft prior regularization algorithm. The corresponding transect plots of the reconstructed permittivity (top) and conductivity (bottom) along the y-axis are shown in Figure 3.21. 119 (a) (b) (c) 120 Figure 3.21. Transect plots of the 1300 MHz reconstructed soft prior permittivity (top) and conductivity (bottom) profiles along the y-axis for three simulation cases with a false region of interest. Dielectric properties of the actual inclusion were set to those of the (a) 40:60, (b) 55:45 and (c) 70:30 mixtures of glycerin and water, reported in Table 3.3 The false region appears as a weak decrease in both permittivity and conductivity images when the contrast between the true target-inclusion and the background medium (86:14 mixture of glycerin and water) is higher, as in case (a) 40:60 glycerin:water inclusion. The presence of the false region declines as the concentration of glycerin and water in the actual target-inclusion gets closer to that of the background medium (86% glycerin). In order to quantify the effect of inclusion property contrast on the sensitivity to the false region of interest, the relative errors between the maximum soft prior reconstructed properties and the actual values are plotted for all three cases in Figure 3.22. 8% % False Inclusion Permittivity Error 6% 4% % False Inclusion Conductivity Error 2% 0% 40% glyrecin 55% glycerin 70% glycerin Figure 3.22. Relative errors between the maximum soft prior recovered properties and the actual values in Figure 3.21 In both permittivity and conductivity images, the relative error in the false inclusion region is always under 7%. In fact, this error decreases as the concentration of glycerin and water in the actual target inclusion gets closer to that of the background medium. In general, conductivity 121 images are less affected by the presence of the false prior structural information. In particular, when dielectric properties of the inclusion are set to those of 55% and 70% glycerin, the relative errors of the reconstructed conductivity values in the false inclusion region with respect to the true background conductivity values are very close to zero. These results show that the soft prior regularization algorithm can effectively handle misleading situations where an incorrect region of interest is identified. The reconstructed dielectric property distributions in an erroneously identified region of interest are clearly different from those in the true region of interest, and they tend to be very close to those of the background medium. 3.4.1.6. Simulation Based On the MR Images of a Normal Patient In order to study the behavior of the soft prior regularization approach in a more realistic scenario where prior clinical structural information of a subject is available, an MR image of a patient’s breast, as illustrated in Figure 3.23(a), was used to create several breast-shaped simulation experiments. The same MR image was also used to generate a two-region conformal soft prior mesh in Figure 3.23(b), composed of 1465 nodes and 2706 triangular elements. 122 (a) (b) Figure 3.23. (a) Patient’s breast MR image (b) The customized soft prior mesh composed of 1465 nodes and 2706 triangular elements containing only internal structure of the breast – adipose (orange) and fibroglandular (blue) tissue, in (a). Reconstructed dielectric property values were extracted along the dark blue line across the breast and fibroglandular tissues. Four sets of dielectric property values were assigned to the identified region of interest (fibroglandular region), resulting in four different breast simulation experiments. These values, along with the dielectric properties of the breast and the background medium at 1300 MHz are presented in Table 3.4. 123 Table 3.4. Dielectric property values of the background medium, breast, and fibroglandular region, in four breast simulation experiments, cases (a) – (d), at 1300 MHz εr σ Background Medium All cases Breast All Cases Region of interest Case (a) Region of interest Case (b) Region of interest Case (c) Region of interest Case (d) 15.3 0.98 11.0 0.90 27.0 1.14 35.0 1.45 40.0 1.66 48.0 2.00 Figure 3.24(a) – (d) show the 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images of the breast-shaped simulation experiments with -100 dBm added noise using the soft priors (left) and no priors (right) for cases (a) – (d) (as described in Table 3.4), respectively. The soft prior images were reconstructed on the two-region conformal mesh in Figure 3.23(b), whereas the no prior images were reconstructed on a uniform 1193 node mesh within the imaging domain (containing the breast and part of the background medium around it.) (a) (b) 124 (c) (d) Figure 3.24. Simulated 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images with -100 dBm added noise, using the soft priors (left) and no priors (right) for cases (a) – (d) in Table 3.4 Comparing the two regularization techniques in terms of finding the correct shape and location of the region of interest, the soft prior regularization technique is significantly superior to the no priors, mainly due to the fact that the prior structural information is used in the soft prior algorithm. In the no prior images on the right, the shape of the region of interest is different in the reconstructed permittivity and conductivity images. More specifically, the ‘C’ shape appeared in the conductivity images (especially when the dielectric property contrast is higher, as in case (c) and case (d)) and is closer to the exact shape of the fibroglandular region. In accordance with the true property contrast increase from cases (a) to (d), both the soft prior and no prior recovered property distributions increase from cases (a) to (d). In order to analyze the results more quantitatively, transects of the reconstructed permittivity (top) and conductivity 125 (bottom) profiles in Figure 2.24 were extracted along an arbitrary line (across the breast and fibroglandular region, as shown in Figure 3.23(b)) and are plotted in Figure 3.25. (a) (b) 126 (c) (d) Figure 3.25. Comparison of the 1300 MHz recovered permittivity (top) and conductivity (bottom) values along an arbitrary line (illustrated in Figure 3.23(b)) in Figure 2.24: (a) – (d) correspond to different target values in case (a) – (d) in Table 3.4 While the no prior images show the fibroglandular region at a roughly correct location, both permittivity and conductivity components of the soft prior images contain the complex shape and exact location of this region. In terms of the recovered dielectric property distributions, the soft prior permittivity values nearly overlap the exact solution in all 4 cases (a) – (b), while the conductivity counterparts are underestimated in the breast region and overestimated in the fibroglandular region. The soft prior conductivity mismatch in the fibroglandular region becomes more prominent as the dielectric property contrast between regions increases from case (a) toward case (d). In particular, the relative soft prior conductivity 127 error in the fibroglandular region increases from 4% in case (a) by nearly a factor of 8 to 31% in case (d). Gauss–Newton image reconstruction in microwave imaging is formulated in terms of the complex wave number squared (k2), from which the relative permittivity and conductivity values 2 are extracted. The magnitude of the average real and imaginary components of k can be significantly out of balance depending on the operating frequency and material characteristics of the object being imaged, and as a result, the reconstructed property distribution values can be imbalanced toward one parameter (relative permittivity or conductivity) or the other. Therefore, in order to balance the property recovery, a pre-scaling procedure is required at the property update stage of the reconstruction [131]. Since the soft prior recovered conductivity values in Figure 3.25 are significantly overestimated, the soft prior images were re-reconstructed with exactly the same reconstruction parameters, but with a 10 times smaller scaling factor applied to the conductivity values at the nodes inside the fibroglandular region. Figure 3.26 shows the transects of the new reconstructed permittivity (top) and conductivity (bottom) profiles with the original and new scaling factor along the same arbitrary line in Figure 3.25 using the soft prior regularization. 128 (a) (b) (c) (d) 129 Figure 3.26. Comparison of the 1300 MHz soft prior reconstructed permittivity (top) and conductivity (bottom) values with the standard (Stnd w) and new (w*0.1) scaling factor along the same arbitrary line shown in Figure 3.23(b): (a) – (d) correspond to different target values in case (a) – case (d) in Table 3.4. Since the scaling factor applied to the permittivity values was kept the same, the reconstructed permittivity values in Figure 3.25 and Figure 3.26 are very similar. On the contrary, due to the using of a smaller scaling factor, the conductivity mismatch in Figure 3.26 is significantly reduced. In fact, the relative soft prior conductivity errors in the fibroglandular region in cases (a), (b), (c), and (d) decrease to 4%, 6%, 9%, and 15%, respectively. A complete summary of the relative dielectric property errors in all four cases with two different scaling factors is shown in Table 3.5. The negative sign indicates that the reconstructed property value is underestimated (i.e. less than the exact property values). Table 3.5. Soft prior relative dielectric property errors in the fibroglandular region, breast-shaped simulation experiments with two different scaling factors Case (a) Case (b) Case (c) Case (d) Relative Error εr, Default Scaling Factor 1% -0.3% -0.5% -1% Relative Error εr, Smaller Scaling Factor 1% -1% -2% -3% Relative Error σ, Default Scaling Factor 4% 17% 27% 31% Relative Error σ, Smaller Scaling factor 4% 6% 9% 15% 130 3.4.2. Phantom Experiments In order to study the behavior of the soft prior regularization with experimental data, several phantom experiments with different geometries and setups are performed at multiple frequencies. The results are presented in the following sections. 3.4.2.1. No priors Versus the Soft prior Regularization The geometry used in this study was the same as the simulations with the circular inclusion described in 3.4.1.1. Specifically, a 1.4 cm radius, thin-walled plastic cylinder filled with a mixture of 55% glycerin and 45% water was offset by 3 cm along the y-axis in the background medium (86:14 glycerin:water mixture). Figure 3.27 shows the 1300 MHz reconstructed images using the no priors and soft prior regularization. 131 (a) (b) Figure 3.27. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images from a phantom experiment with a single circular inclusion for (a) no priors on the 473 node mesh and (b) soft priors on the 915 node mesh in Figure 3.8 In both instances, the inclusion is evident. Several background artifacts appear in the no prior images, which are more prominent in the conductivity parameter. Incorporating the spatial priors substantially improves the quality of both the permittivity and conductivity images. Weighted permittivity and conductivity errors decrease from 0.329 and 0.302 to 0.015 and 0.045, respectively, when the spatial structure of the phantom is incorporated through the soft prior regularization. Figure 3.28 shows the extracted reconstructed dielectric properties along the yaxis from the two approaches relative to the exact values. 132 Figure 3.28. Comparison of the 1300 MHz reconstructed permittivity (left) and conductivity (right) profiles extracted along the y-axis in Figure 3.26 phantom experiment using the no priors and soft prior regularizations. Clearly, the soft prior regularization dramatically reduces the spatial oscillations within the background, but it also recovers the dielectric properties in the inclusion more accurately. Improvements are even more prominent in the conductivity images in which case the recovered inclusion values are over-estimated and displaced (toward the boundary) without priors. These artifacts are eliminated when the spatial structure of the phantom is incorporated through the soft prior regularization. 3.4.2.2. Choice of the Soft Prior Coefficient In all previous results, the soft prior coefficient, λ in equation 3.4, was set to unity as a default, but in this section, a more detailed study of the effects of λ as a function of frequency is presented. Data from the same experimental setup described in 3.4.2.1 was used and images were reconstructed at 900, 1100, 1300, 1500, 1700, 1900, 2100, 2300, and 2500 MHz. The 133 independently measured dielectric properties of the coupling medium and the inclusion are reported in Table 3.6 as a function of frequency. Table 3.6. Independently measured dielectric properties of the background medium (BK) and target inclusion (Inc) over the range of frequencies evaluated Frequency (MHz) εr, BK σBK εr, Inc σInc 900 1100 1300 1500 1700 1900 19.3 17.3 15.6 14.4 13.6 12.8 0.72 0.82 0.90 0.97 1.03 1.10 56.1 53.5 51.2 48.9 46.4 44.5 0.80 1.12 1.44 1.77 2.09 2.42 2100 2300 2500 12.1 11.5 11.1 1.14 1.20 1.24 42.3 40.4 38.8 2.74 3.04 3.35 Based on testing over a wide range of values in simulation and phantom experiments, a spectrum of the soft prior coefficients (λ = 0.01, 0.1, 1, 10, and 100) was used for the reconstruction. Lower values, such as 0.001 for λ, allowed the solution to diverge in some cases, while higher values tended to suppress the recovered inclusion properties. Transects along the y-axis of the reconstructed images for εr and σ at 900, 1100, 1300, 1500, 1700, 1900, 2100, 2300, and 2500 MHz for the range of soft prior coefficients are shown in Figure 3.29. 134 (a) 900 MHz (b) 1100 MHz (c) 1300 MHz 135 (d) 1500 MHz (e) 1700 MHz (f) 1900 MHz 136 (g) 2100 MHz (h) 2300 MHz (i) 2500 MHz 137 Figure 3.29. Comparison of the reconstructed permittivity (left), and conductivity (right) values from a phantom experiment using different spatial prior coefficient λ values at (a) 900, (b) 1100, (c) 1300, (d) 1500, (e) 1700, (f) 1900, (g) 2100, (h) 2300, and (i) 2500 MHz In general, at lower frequencies up to 1700 MHz (e), no significant difference occurs between the reconstructed values when λ = 0.1, λ = 1.0 and λ = 10. When λ = 0.01, the recovered permittivity values in the inclusion region are close to the exact values, while the conductivity counterparts are slightly underestimated. In addition, the recovered properties of the background region begin to deviate from the true levels as the frequency increases. At lower frequencies, the λ = 100 reconstructions appear to estimate the inclusion conductivity values closer to the true levels than the corresponding permittivity values, which are noticeably underestimated across all frequencies. This increasing mismatch may be due to the excessive smoothing of the soft prior coefficient. In order to analyze the results quantitatively, the weighted permittivity and conductivity errors associated with each image were computed and are plotted as a function of the soft prior coefficient λ in Figure 3.30. 138 (a) (b) Figure 3.30. Weighted (a) εr and (b) σ errors for a phantom experiment over a range of frequencies from 900 to 2500 MHz using six different soft prior coefficients: λ = 0.01, 0.1, 1, 10, 100, and 1000 139 The weighted property errors for λ = 100 and λ = 1000 are elevated at all frequencies, when compared to the λ = 1 and λ = 10 cases. Similar behavior is observed for λ = 0.01 and λ = 0.1, but reducing the soft prior weighting coefficient λ appears to be more sensitive at higher frequencies, specifically over 1700 MHz. These results suggest that values in the range of λ=1.0 and λ=10 appear to be appropriate soft prior weighting coefficients over our reconstruction frequency range (usually from 1100 to 1700 MHz). However, the stability of the weighting parameter λ investigated here is based on the present set of experiments, and may need to be adjusted in other configurations. 3.4.2.3. Size of the Target Three thin-walled plastic cylinders of 0.65, 1.4, and 2.1 cm radii filled with a mixture of 50% glycerin and 50% water were centered at (x, y)Inc = (0, 3) cm in a background medium comprised of 80:20 glycerin:water mixture. Microwave data was then collected at a frequency range from 900 to 1700 MHz in steps of 200 MHz using our experimental imaging Metrikos System (Metrikos, Inc., Hopkinton, MA), and the images were reconstructed using the no priors on the 473 node mesh in Figure 3.8(a), and the soft prior regularization on the customized reconstruction meshes illustrated in Figure 3.31. 140 (a) (b) (c) Figure 3.31. Customized soft prior reconstruction meshes used for the phantom experiments with a circular inclusion of radius: (a) 0.65 cm – 1329 nodes and 2530 triangular elements, (b) 1.4 cm – 1037 nodes and 1943 triangular elements, and (c) 2.1 cm – 813 nodes and 1498 triangular elements Figure 3.32 shows transects of the reconstructed permittivity (left) and conductivity (right) images along the y-axis at a frequency range from 900 to 1700 MHz in steps of 200 MHz, using the no priors (LM) and soft prior (SP) along with the exact properties of the inclusion. (a) 900 MHz 141 (b) 1100 MHz (c) 1300 MHz (d) 1500 MHz 142 (e) 1700 MHz Figure 3.32. Comparison of the reconstructed permittivity (left), and conductivity (right) values from three phantom experiments with different sized inclusions (of radius 0.65, 1.4, and 2.1 cm) using the no prior (LM) and soft prior (SP) regularizations at (a) 900, (b) 1100, (c) 1300, (d) 1500, and (e) 1700 MHz Clearly, the soft prior regularization not only reduces the spatial oscillations within the background in all cases and at all frequencies, but also recovers the dielectric properties in the smallest inclusion (of 0.65 cm radius) much more accurately. The soft prior improvements are even more prominent in the conductivity images at higher frequencies in which case the recovered inclusion values are over-estimated, and in some cases the inclusion is displaced (toward the boundary) without priors. These artifacts are eliminated when the spatial structure of the phantom is incorporated through the soft prior regularization. Except for the last frequency (1700 MHz), the recovered soft prior permittivity values of the mid-sized inclusion (of 1.4 cm radius) are generally closer to the exact solution, while the most accurate reconstructed soft prior conductivity values appear in the smallest inclusion. 143 3.4.2.4. Shape of the Target Three different base shaped thin-walled plastic tubes including a 1.3 cm diameter circleshaped, a 1.35 cm side length square-shaped and a 1.35 cm side length diamond-shaped cylinders were filled with a mixture of 50% glycerin and 50% water, and then were centered at (x, y)Inc = (0, 3) cm in a background medium comprised of 80:20 glycerin:water mixture. Similar to the previous section, microwave data was collected at a frequency range from 900 to 1700 MHz in steps of 200 MHz using the Metrikos System, and the images were reconstructed with no priors on the 473 node mesh in Figure 3.8(a) and with the soft prior regularization on the customized reconstruction meshes illustrated in Figure 3.33. (a) (b) (c) Figure 3.33. Customized soft prior reconstruction meshes used for the phantom experiments with: (a) circular inclusion – 1329 nodes and 2530 triangular elements, (b) square-shaped inclusion – 1096 nodes and 2096 triangular elements, and (c) diamond-shaped inclusion – 1040 nodes and 1984 triangular elements 144 Figure 3.34 shows transects of the reconstructed permittivity (left) and conductivity (right) images along the y-axis at a frequency range from 900 to 1700 MHz in steps of 200 MHz, using the no priors (LM) and soft prior (SP) along with the exact properties of the inclusion. (a) 900 MHz (b) 1100 MHz 145 (c) 1300 MHz (d) 1500 MHz (e) 1700 MHz 146 Figure 3.34. Comparison of the reconstructed permittivity (left), and conductivity (right) values from a phantom experiment with three different shaped but similar sized inclusions (1.3 cm diameter circular, 1.3 cm side length square, and 1.3 cm side length diamond-shaped cylinders) using the no prior (LM) and soft prior (SP) regularizations at (a) 900, (b) 1100, (c) 1300, (d) 1500, and (e) 1700 MHz Similar to the previous phantom experiment results presented in 3.4.2.3, when the soft prior regularization is used, the spatial oscillations within the background are significantly damped in all cases and at all frequencies. Overall, the no prior (LM) reconstructed dielectric properties are very similar in all three inclusions, which indicates that the Levenberg-Marquardt regularization is not very sensitive to the shape of the region of interest. Alternatively, the soft prior regularization appears to respond quite differently to different inclusion shapes. In particular, over all frequencies the recovered dielectric properties in the circular inclusion are more accurate than those in the square and diamond-shaped inclusions in which case the recovered soft prior permittivity values are underestimated, whereas their conductivity counterparts are overestimated. This behavior may be due to the sharp corners in the square and diamond-shaped cylinders. 3.4.2.5. Dielectric Property Contrast of the Target In order to study the effects of dielectric property contrast between the target and background medium on the soft prior reconstructed property distributions, five phantom experiments with a relatively small target-inclusion were performed. The schematic setups in all cases were alike: A circular thin-walled plastic tube of 0.65 cm radius filled with different mixtures of glycerin and water (50:50, 60:40, 70:30, 90:10, and 100:0 glycerin:water solutions) 147 was centered at (x, y)Inc = (0, 3) cm in a background medium comprised of 80:20 glycerin:water mixture. Microwave data was collected at a frequency range from 900 to 1700 MHz in steps of 200 MHz using the Metrikos System, and the soft prior images were reconstructed on the customized 1329 node mesh, illustrated in Figure 3.33(a). The maximum recovered soft prior permittivity (left) and conductivity (right) of the five target-inclusions with different contrast-levels against the independently measured properties (i.e. the exact values) are plotted in Figure 3.35. For analytical purposes, the linear least-square fits, along with the linear regression slopes and their corresponding R2 values were calculated and they are shown in Figure 3.35. (a) 900 MHz 148 (b) 1100 MHz (c) 1300 MHz (d) 1500 MHz 149 (e) 1700 MHz Figure 3.35. Comparison of the soft prior reconstructed permittivity (left) and conductivity (right) values of five different contrast-level target-inclusions (50:50, 60:40, 70:30, 90:10 and 100:0 mixture of glycerin:water) at (a) 900, (b) 1100, (c) 1300, (d) 1500, and (d) 1700 MHz, along with linear least-square fits In Figure 3.35, concentration of glycerin (%g) in the target-inclusion is labeled next to each data point. Two metrics are considered in this study: the slope of the linear least-square fit shows the pattern of the soft prior reconstructed properties at different contrast-levels, whereas the R2 value is a measure of goodness-of-fit. In an ideal case, the linear fit would have both slope and R2 value of unity – i.e. all data points would fall on the line y = x. Data points above y = x indicate that the recovered dielectric properties are overestimated (reconstructed value > exact value), whereas those underneath represent underestimated property values. In general, the permittivity R2 values across all frequencies and the conductivity R2 values at 1100, 1300 and 1500 MHz are very close to unity (> 0.93), which confirms a strong linear correlation between the soft prior reconstructed property values of different contrast-level inclusions. In addition, the slopes of the linear fits at 1100, 1300 and 1500 MHz are also close to unity (1), indicating the 150 correct pattern of the soft prior reconstructed properties at different contrast-levels. Overall, the soft prior recovered properties of the inclusions are more degraded at the lowest and the highest frequencies (900 and 1700 MHz). In fact, the mismatch between the reconstructed and actual conductivity values is significantly more prominent. Nonetheless, it should be noted that the range of actual conductivity values of different mixture of glycerin and water used in this experiment is constricted at 900 MHz (varying only from 0.33 to 0.88 S/m). 3.4.2.6. Number of the Reconstruction Nodes In order to study how the soft prior recovered dielectric properties are affected by the number of nodes in the customized reconstruction mesh, measurement data from the singletarget phantom experiment described in 3.4.2.1 was used to reconstruct soft prior images at 1300 MHz, using a set of customized meshes with different number of nodes in each region. In this analysis, three cases were considered: (a) increasing number of nodes in the target-inclusion region (NInc), while keeping number of nodes in the background region (NBK) in the same bulk part, (b) increasing number of nodes in the background region (NBK), while number of nodes in the target-inclusion region (NInc) stays in the same range, and (c) increasing number of nodes in both regions (NBK and NInc) proportionally. Table 3.7 contains all details about the corresponding soft prior meshes in each case, including the total number of nodes (NTotal), number of nodes on the perimeter of the background (NΩBK) and target inclusion (NΩInc) regions, maximum area of triangular elements (maxA) in each region, and number of nodes in each region (N). 151 Table 3.7. Properties of the soft prior customized meshes used for the phantom experiment with a single target inclusion of 1.4 cm radius Case (a) more nodes in the target-inclusion region NTotal 255 NΩBK 30 maxABK 0.5 NΩInc 6 maxAInc 3 NBK 247 NInc 8 251 269 288 337 30 30 30 30 0.5 0.5 0.5 0.5 10 15 20 25 0.5 0.3 0.1 0.05 235 243 233 233 16 26 55 104 754 30 0.5 25 0.01 267 487 1304 1955 5224 30 0.5 25 0.005 315 30 0.5 25 0.003 327 30 0.5 25 0.001 410 Case (b) more nodes in the background region 989 1628 4814 NTotal 405 1159 2303 NΩBK 30 40 50 maxABK 0.3 0.1 0.05 NΩInc 10 10 15 maxAInc 0.5 0.5 0.3 NBK 389 1137 2267 11327 60 NTotal 251 NΩBK 30 maxABK 0.5 NΩInc 10 maxAInc 0.5 NBK 235 NInc 16 414 1204 2366 11686 30 40 50 60 0.3 0.1 0.05 0.01 15 25 30 30 0.3 0.1 0.05 0.01 388 1141 2260 11167 26 63 106 519 0.01 15 0.3 11238 Case (c) more nodes in both regions NInc 16 22 36 89 For each of the three cases above, the maximum soft prior recovered property values in the inclusion region were extracted, and they are plotted as a function of the number of reconstruction nodes (in logarithmic scale) in Figure 3.36. It should be noted that since the 152 recovered dielectric properties of the background medium nearly overlap the exact properties in all cases, only reconstructed values of the target are presented here. (a) NInc (b) NBK 153 (c) NTotal Figure 3.36. Comparison of the soft prior reconstructed permittivity (left), and conductivity (right) values of a target inclusion using several soft prior customized meshes with different node densities: case (a) number of nodes in the inclusion region (NInc) changes, case (b) number of nodes in the background region (NBK) changes, and case (c) number of nodes in both regions (NTotal) changes proportionally In all three cases, the soft prior recovered inclusion conductivity values (plots on the right) appear to be less sensitive to the number of reconstruction nodes. In fact, all reconstructed conductivity values are within roughly 6% of the exact values. Nonetheless, the optimal inclusion conductivity value is reached when approximately 25-100 nodes are placed in the inclusion region. That is when the 269 node mesh is used in case (a) (i.e. 3rd point from the left in Figure 3.36(a)), the 1159 node mesh is used in case (b) (i.e. 2nd point from the left in Figure 3.36(b)), and the 414 node mesh is used in case (c) (i.e. 2nd point from the left in Figure 3.36(c)). On the other hand, changing the number of reconstruction nodes has a more significant effect on the recovered inclusion permittivity values. More specifically, as the number of nodes increases, generally the soft prior permittivity values of the target inclusion are recovered more 154 accurately in all three cases (as much as %20). This observation suggests that the denser reconstruction mesh can improve the accuracy of the soft prior reconstructed permittivity values. However, no obvious reason is found for the slight decrease in the permittivity values of the 4th data point in Figure 3.36(a). Intuitively, when additional nodes are inserted in the non-target region (here, the background) as in case (b), the recovered dielectric properties of the target inclusion should not change. However, the reconstructed inclusion permittivity values in Figure 3.36(b) are improved as the number of nodes in the background region increases. This may be due to the fact that the number of nodes in the inclusion region also increases (although at a much lower rate). In effect, as shown in Table 3.7, the number of inclusion nodes increases from 16 to 89 as the number of background nodes rises from 389 to 11238. It would be desirable to keep the number of inclusion nodes constant, but in order to comply with certain mesh quality criteria, additional nodes ought to be inserted in the inclusion region as the number of background nodes was increased. 3.4.2.7. Sensitivity to a False Region of Interest In order to study the sensitivity of the soft prior algorithm to false regions of interest in measurement data, we used the data generated in the phantom experiment in 3.4.2.1, where only a single inclusion at the lower location (0, –3 cm) actually existed. The customized reconstruction mesh in Figure 3.20 containing a false inclusion region of 1.4 cm radius centered at (0, 3 cm), and the actual target zone was used to reconstruct the soft prior images at 1300 MHz. Figure 3.37 shows the 1300 MHz reconstructed images using the soft prior regularization coefficients of λ = 0.01, λ = 0.1, λ = 1.0 and λ = 10, respectively. The false region appears as a 155 weak increase (~ 3 to 6%) in the permittivity images, but with a more prominent decrease (~ – 4% to – 20%) in the conductivity images. (a) λ = 0.01 (b) λ = 0.1 (c) λ = 1.0 (d) λ = 10 Figure 3.37. 1300 MHz soft prior reconstructed permittivity (top) and conductivity (bottom) images of a phantom experiment using the reconstruction mesh with a false region of interest in Figure 3.20 Transect plots along the y-axis in the reconstructed permittivity and conductivity images are presented in Figure 3.38. 156 Figure 3.38. Comparison of the 1300 MHz soft prior reconstructed permittivity (left), and conductivity (right) profiles from a phantom experiment using the reconstruction mesh with a false region of interest in Figure 3.20 Consistent with the images in Figure 3.37, the permittivity values within the false inclusion region in Figure 3.38 are close to the true background, while the counterpart conductivity profiles are more noticeably affected by the misleading information and they exhibit lower properties than those of the background medium, especially when larger values of λ are used. Since the contribution from the soft prior regularization is reduced when λ = 0.01, the dielectric properties are not significantly influenced by the presence of the false inclusion region (~ 3% and ~ – 4% error in the permittivity and conductivity images, respectively). Nonetheless, for smaller values of λ, more artifacts are observed in the background region. However, as the soft prior weighting coefficient λ increases to 1.0 and 10, the level of the error increases in the false region of interest (~ 6% and ~ – 20% in the permittivity and conductivity images, respectively). A summary of the weighted soft prior property errors is presented in Table 3.8. 157 Table 3.8. Weighted soft prior εr and σ errors for a phantom experiment using a reconstruction mesh with a false region of interest for different soft prior coefficients: λ = 0.01, 0.1, and 1.0 at 1300 MHz Frequency Errw, ε r (MHz) λ =0.01 1300 0.232 Errw, σ λ =0.01 Errw, ε r λ =0.1 Errw, σ λ =0.1 Errw, ε r λ =1 Errw, σ λ =1 Errw, ε r λ =10 Errw, σ λ =10 0.246 0.143 0.159 0.131 0.162 0.127 0.157 For λ = 1, Errw, ε r = 0.131 and Errw, σ = 0.162, which are larger relative to the case when the exact spatial structure of the phantom was used in Section 3.4.2.1 ( Errw, ε r = 0.015 and Errw, σ = 0.045), but significantly lower than those obtained without priors ( Errw, ε r = 0.329 and Errw, σ = 0.302). 3.4.2.8. Sensitivity to Imperfect Spatial Priors Obtaining perfect a priori structural information of the tissue being imaged may not be feasible in practice, and therefore, evaluating the sensitivity of the soft prior technique to imperfect priors is essential. In this section, two types of priors’ imperfections, namely the size and location of the target, are studied: a) Imperfect prior size of the target inclusion: Table 3.9 presents some attributes of a set of different inclusion-sized soft prior meshes that was used to reconstruct the same data generated in the phantom experiment described in 3.4.2.1, with a target inclusion of 1.4 cm radius centered at (0, –3 cm). 158 Table 3.9. Characteristics of the soft prior reconstruction meshes used for a phantom experiment with imperfect prior size of the target inclusion Inclusion Location Actual Inclusion Prior Inclusion Number of (In the actual Total Size Size (In the Nodes in the experiment and in Number of (In the experiment) mesh) Inclusion the mesh) Nodes Radius (cm) Radius (cm) Region Center (cm) (0, –3 cm) (0, –3 cm) (0, –3 cm) (0, –3 cm) (0, –3 cm) (0, –3 cm) (0, –3 cm) 1.4 1.4 1.4 1.4 1.4 1.4 1.4 0.8 1.0 1.2 1.4 1.6 1.8 2.0 873 1054 969 1036 978 954 874 217 288 288 318 318 326 318 Figure 3.39 shows the 1100 MHz reconstructed images using different soft prior meshes in Table 3.9. (a) 0.8 (b) 1.0 (c) 1.2 (d) 1.4 (e) 1.6 (f) 1.8 (g) 2.0 cm Figure 3.39. 1100 MHz reconstructed permittivity (top) and conductivity (bottom) images of a phantom experiment using different prior inclusion-sized meshes: (a) 0.8 cm, (b) 1.0 cm, (c) 1.2 cm, (d) 1.4 cm (the true size), (e) 1.6 cm, (f) 1.8 cm, and (g) 2.0 cm were considered as the radius 159 of the inclusion in the reconstruction mesh. The actual size (radius of 1.4 cm) and location (0, – 3 cm) of the inclusion is outlined by a circle on the reconstructed images. As the inclusion radius in the reconstruction mesh increases from Figure 3.39(a) – (g) the recovered inclusion permittivity values decrease, while the counterpart conductivity values increase. In order to compare these values with the actual property distributions, the reconstructed permittivity and conductivity values of the background and target inclusion in Figure 3.39 were extracted and they are plotted in Figure 3.40 as a function of the inclusion radius in the reconstruction mesh. Figure 3.40. Comparison of the extracted 1100 MHz reconstructed permittivity (left), and conductivity (right) values in Figure 3.39 along with the exact property distributions for a phantom experiment using different inclusion–sized soft prior meshes Consistent with the reconstructed images in Figure 3.39, as the radius of the inclusion in the soft prior mesh increases from 0.8 to 2 cm, the inclusion recovered permittivity values in Figure 3.40 decrease, while the counterpart conductivity values increase. In addition, the recovered and exact background values overlap very closely. As expected, when the exact size of 160 the inclusion is used in the reconstruction mesh, the most accurate inclusion permittivity value is obtained. Using prior structural information with smaller–sized inclusions misleads the soft prior algorithm to overestimate the inclusion permittivity values, while using larger–sized inclusions results in underestimation of properties. However, imperfect prior information of the size of the inclusion has a very different effect on the recovered conductivity values. As the inclusion radius in the soft prior mesh increases, the recovered conductivity values approach the exact solution, and as a result, the most accurate values are not obtained with the exact–sized target inclusion. At 1100 MHz, there is about 30% change in recovered permittivity and 27% change in recovered conductivity values when approximately 78% change is applied to the area of the target region (changing the inclusion radius from 1.2 to 1.6 cm). In order to analyze the sensitivity to imperfect size of the target at other operating frequencies, the reconstructed permittivity and conductivity values of the background and target– inclusion at a frequency range from 1300 to 1700 MHz in steps of 200 MHz were extracted and plotted versus the inclusion radius in the reconstruction mesh in Figure 3.41. (a) 1300 MHz 161 (b) 1500 MHz (c) 1700 MHz Figure 3.41. Comparison of the extracted (a) 1300, (b) 1500, and (c) 1700 MHz reconstructed permittivity (left), and conductivity (right) profiles along with the exact property distributions for a phantom experiment using different inclusion–sized meshes Overall, at higher frequencies in Figure 3.41 the same pattern as the 1100 MHz reconstructed property values in Figure 3.40 is seen: As the inclusion radius in the reconstruction mesh increases from 0.8 to 2 cm, the recovered inclusion permittivity values decrease, while the 162 counterpart conductivity values increase. In terms of the relative property variations, when approximately 78% change is applied to the area of the target region, the recovered permittivity values at 1300, 1500 and 1700 MHz vary about 28%, 30%, and 28% respectively, while the counterpart conductivity changes decrease from 13% at 1300 MHz to only 3% at 1700 MHz. b) Imperfect prior location of the target inclusion: Table 3.10 presents some attributes of a set of different inclusion–located soft prior meshes that was used to reconstruct the same data generated in the phantom experiment described in 3.4.2.1, with a target inclusion of 1.4 cm radius centered at (0, – 3 cm). Table 3.10. Characteristics of the soft prior customized meshes used for a phantom experiment with imperfect prior location of the target inclusion Inclusion Size Inclusion Location Prior Inclusion Number of (In the actual Total (In the Location Nodes in the experiment and in Number of experiment) (In the mesh) Inclusion the mesh) Nodes Center (x, y) cm Center (x, y) cm Region Radius (cm) 1.4 1.4 (0, –3) (0, –3) (0, –4.5) (0, –4.2) 882 895 288 288 1.4 1.4 1.4 1.4 1.4 (0, –3) (0, –3) (0, –3) (0, –3) (0, –3) (0, –3.9) (0, –3.6) (0, –3.3) (0, –3) (0, –2.7) 897 906 914 912 911 288 288 288 288 288 1.4 1.4 1.4 (0, –3) (0, –3) (0, –3) (0, –2.4) (0, –2.1) (0, –1.8) 909 903 904 288 288 288 1.4 1.4 1.4 1.4 1.4 (0, –3) (0, –3) (0, –3) (0, –3) (0, –3) (0, –1.5) (0, –1.2) (0, –0.9) (0, –0.6) (0, –0.3) 894 913 924 899 920 288 288 288 280 288 163 1.4 1.4 (0, –3) (0, –3) (0, 0) (0, 0.3) 899 916 288 288 1.4 1.4 1.4 1.4 (0, –3) (0, –3) (0, –3) (0, –3) (0, 0.6) (0, 0.9) (0, 1.2) (0, 1.5) 899 909 916 869 288 288 288 280 Figure 3.42 shows the 1100 MHz reconstructed images using some of the soft prior meshes in Table 3.10. (a) -3.9 (b) -3.0 (c) -2.1 (d) -1.2 (e) -0.3 (f) 0.6 (g) 1.5 cm Figure 3.42. 1100 MHz reconstructed permittivity (top) and conductivity (bottom) images of a phantom experiment using different prior inclusion–located meshes: x = 0 and y–coordinates of (a) – 3.9, (b) – 3.0 (the true location), (c) – 2.1, (d) – 1.2, (e) – 0.3, (f) 0.6, and (g) 1.5 cm were considered as the center of the inclusion in the reconstruction mesh. The actual size (radius of 1.4 cm) and location (0, – 3 cm) of the inclusion is outlined by a circle on the reconstructed images. As the inclusion in the reconstruction mesh moves away from its actual location in the phantom experiment setup (black circle in Figure 3.42), the recovered permittivity profiles decrease toward the background, while the corresponding conductivity profiles slightly increase and then decrease lower than the background conductivity values. 164 In order to analyze the results more quantitatively, the reconstructed permittivity and conductivity values of the background and target inclusion in Figure 3.42 were extracted along the y–axis and they are plotted in Figure 3.43(a) versus the inclusion center in the reconstruction mesh. The same set of property values reconstructed at 1300, 1500, and 1700 MHz are presented in Figure 3.43(b), (c), and (d), respectively. (a) 1100 MHz (b) 1300 MHz 165 (c) 1500 MHz (d) 1700 MHz Figure 3.43. Comparison of the extracted (a) 1100, (b) 1300, (c) 1500, and (d) 1700 MHz reconstructed permittivity (left), and conductivity (right) values along with the exact property distributions for a phantom experiment using different inclusion–located meshes Overall, when the prior and actual target locations match, the reconstructed inclusion permittivity values approach the exact profiles. However, the most accurate values are not recovered when the mismatch is zero, but when the inclusion is offset along the y–axis about 0.3 166 cm. In fact, this can be due to the slight misplacement of the actual inclusion tube during the phantom experiment. When the partial overlap between the actual and prior target area decreases, the recovered inclusion permittivity profiles approach the background values. The imperfect prior location of the target seems to have a dissimilar effect on the recovered conductivity profiles. In general, as the center of the inclusion in the reconstruction mesh shifts down from its actual location, the reconstructed inclusion conductivity values are overestimated. Similar behavior is seen when the prior center of the inclusion moves up from y = – 3 cm to about y = – 1.5 cm. The conductivity profiles then start to decrease significantly to the extent that they even go below the background values. In terms of the recovered background properties, the most accurate values at all frequencies are obtained around the exact location of the inclusion. As the mismatch between the prior and actual target locations increases – i.e. the y–coordinate of the inclusion center in the reconstruction mesh moves away from – 3 cm, both reconstructed background permittivity and conductivity values increase, which may be an indication of averaged property distribution between the background and actual target regions. The results presented in this section verified that the soft prior algorithm is not overly sensitive to the exact size and location of the target region and it can tolerate priors’ imperfection to some extent. 3.4.2.9. Breast–Shaped Phantom Experiments In order to study the behavior of the soft prior regularization technique in more practical circumstances, several phantom experiments with different sized target inclusions were performed. In order to conduct real breast–shaped experiments, MR scans of a real breast of cup 167 size B (425 ml) was used to create a breast mold, as shown in Figure 3.44(a). This mold was then used to make the rapid–prototyped plastic breast model in Figure 3.44(b). (a) (b) Figure 3.44. Real breast–shaped model: (a) breast mold, (b) rapid–prototyped breast model The plastic used in the model was chosen to be thin enough so that the microwave signals would not be significantly distorted as they propagated through it. Two cylinders – of 0.65 and 1.0 cm radius – were inserted into the breast model as tumor inclusions. The dielectric properties of the breast were chosen to mimic those of the scattered breast by using a 88:12 mixture of glycerin and water, whereas the inclusions were set to be similar to a tumor with a mixture of 50:50 glycerin:water. The breast phantom was also submerged in a coupling liquid composed of a 80:20 mixture of glycerin and water. The 1300 MHz relative permittivity and conductivity values of the actual phantom components are presented in Table 3.11. 168 Table 3.11. Independently measured dielectric properties of the background medium (BK), breast (BR), and inclusions (Inc) at 1300 MHz Frequency (MHz) εr, BK σBK εr, BR σBR εr, Inc σInc 1300 15.3 0.98 13.0 0.83 58.6 1.02 The required a priori structural information of the phantom setup for the soft prior regularization was obtained before microwave data acquisition by using our 3D optical scanner shown in Figure 3.45 (VIUscan color laser scanner, CreaForm Inc. Quebec, Canada). Figure 3.45. 3D optical scanning of the model First, the breast model as well as the two tumor inclusions was scanned by the laser scanner. Next, the resulting 3D surface mesh, illustrated in Figure 3.46(a), was imported into the SolidWorks (SolidWorks Corp., Concord, MA) software package, where the solid object in Figure 3.46(b) was formed and 2D slices were extracted. 169 (a) (b) Figure 3.46. Exported 3D model (a) 3D surface mesh, (b) Formation of a solid object and extraction of 2D slices in the SolidWorks software package Finally, different 2D customized meshes containing the exact location of the tumor inclusions with respect to the breast region were obtained and used in the reconstruction algorithm with soft prior regularization. The soft prior meshes used in the present experiment are illustrated in Figure 3.47. 170 (a) (b) Figure 3.47. 2D customized soft prior meshes used in a breast–shaped phantom experiment with circular target inclusions of (a) 0.65 cm and (b) 1.0 cm radius, as well as a polyline along which the actual and reconstructed dielectric property distributions were extracted. The mesh in (a) is comprised of 1778 nodes and 3468 triangular elements, whereas the one in (b) contains 1788 nodes and 3488 triangular elements. Figure 3.48 shows transects of the 1300 MHz reconstructed permittivity (top) and conductivity (bottom) profiles along the arbitrary polyline across the imaging domain in Figure 3.47, using the no priors (on 473 node mesh) and soft prior regularization (on the 1778 and 1788 node meshes). 171 (a) (b) Figure 3.48. Comparison of the 1300 MHz reconstructed permittivity (left) and conductivity (right) profiles along an arbitrary polyline for the breast–shaped phantom experiments with (a) small inclusion of 0.65 cm radius and (b) medium–sized inclusion of 1 cm radius, using the no priors and soft prior regularization The tumor inclusion is detected both with and without using prior structure of the phantom setup in (a) and (b). Nonetheless, the soft prior images are clearly superior in terms of the recovered dielectric properties of the different regions. The reconstructed soft prior permittivity values of the background medium and the breast nearly overlap the exact solution in Figure 3.48(a) and (b). However, the corresponding tumor inclusion values are not as precise as expected, possibly due to the significant contrast between the breast and tumor inclusion. The 172 mismatch between the recovered soft prior properties and the exact permittivity values of the inclusion deteriorates as its radius decreases from 1 to 0.65 cm. The recovered soft prior conductivity profiles are overall more accurate than the permittivity counterparts; however, similar to the permittivity profile, as the size of the target inclusion becomes smaller, the mismatch between the reconstructed soft prior properties and the exact conductivity values increases. 3.4.2.10. Phantom Experiment in the MR Bore Accurate a priori structural information of the object being imaged can be obtained from different sources. In this section we show the results of combining microwave data with actual MR images obtained simultaneously in the MR bore, as described in 3.3.3. In this experiment, the MR system was acquiring data during the microwave acquisition. The tank was filled with an 80% glycerin solution with dielectric properties of εr, BK = 20.9 and σBK = 1.22 S/m, respectively. The imaging target was a 5.4 cm radius tube filled with an 84% solution of glycerin with 2.1 cm and 1.0 cm radius offset cylindrical inclusions having 72% and 55% glycerin:water mixtures. The properties of the three regions were: 1.02 S/m, εr, Br = 16.1 and σBr = εr, Fg = 30.2 and σFg = 1.51 S/m, and εr, Tu = 51.2 and σTu = 1.44 S/m, respectively. While the three regions are nominally meant to mimic the properties of a fatty breast (Br), fibroglandular tissue (Fg) and tumor (Tu), in actuality they only achieve this for the permittivity. Because of the relaxation behavior of the glycerin:water solutions, their conductivity values in the 1300 MHz region are only fractionally different. However, for the purposes of testing the robustness of the hardware and software under challenging, high contrast situations, this 173 configuration is more than adequate. Figure 3.49 shows the corresponding MR T2–weighted image of the horizontal plane through the target and the center of the active part of the antennas. Figure 3.49. T2–weighted MR images of the phantom experiment at the center plane of the active part of the antennas The location of the target and associated inclusions with respect to the antenna positions is clearly evident. Figure 3.50(a) shows the no prior reconstruction mesh comprised of 559 uniformly distributed nodes and 1044 triangular elements, while Figure 3.50(b) shows the soft prior reconstruction mesh containing 1034 nodes and 1821 elements discretized from the MR image. (Note that the MR image and mesh are horizontal mirror images of each other since the MR system produces images viewed from top down while the microwave system is configured for en face view images – i.e. from underneath – during breast patient exams.) 174 (a) (b) Figure 3.50. 2D reconstruction meshes used for the phantom experiment in MR: (a) No priors, with 559 uniformly distributed nodes and 1044 triangular elements, and (b) soft priors, with 1034 nodes and 1821 triangular elements Figure 3.51(a) and (b) show the recovered 1300 MHz permittivity (top) and conductivity (bottom) images of the target using the no priors and the soft prior regularization exploiting the MR–based prior information. 175 (a) (b) Figure 3.51. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images of the multi–region region phantom case using the (a) no priors (circles indicate the exact ct object locations) loca and (b) soft prior regularization. In Figure 3.51(a), the breast region and both target inclusions are visible in the permittivity and conductivity images, but with a fair number of artifacts within the breast phantom and the surrounding bath. Circles have been drawn over the conventional images to indicate the exact locations of the objects. Overall, tthe soft prior images in Figure 3.51(b) (b) appear to recover the property distributions more accurately accurately,, especially the permittivity profile of the tumor inclusion. inclusion In order to evaluate whether using tthe hard priors would improve the characterization of the tumor inclusion, the MR–based based prior information was used to reduce the size of the 176 reconstruction parameters, as described in 3.2.3. Then, the images were reconstructed using the hard priors. The actual and maximum recovered property values of the fibroglandular and tumor inclusion regions both with and without priors are reported in Table 3.12(a). While the reconstructed properties with no priors are significantly different from the actual values, the soft prior and hard prior values are generally closer to the exact properties. In order to quantify the improvement, the corresponding relative dielectric property errors were calculated and they are presented in Table 3.12(b). The negative sign indicates that the reconstructed property value is underestimated (i.e. less than the exact property values). Table 3.12. (a) The actual and averaged recovered property values of the fibroglandular and tumor inclusion regions in the phantom experiment in the MR bore using the no priors, soft priors, and hard priors, (b) The corresponding relative dielectric property errors εr, Fg σFg εr, Tu σTu Actual values 30.2 1.51 51.2 1.44 Avg. Reconstructed – No priors 25.1 1.28 29.1 1.25 Avg. Reconstructed – Soft priors 28.6 1.41 47.3 1.27 Avg. Reconstructed – Hard priors 28.7 1.41 50.7 1.12 (a) εr, Fg σFg εr, Tu σTu % Relative Error – No priors -16.9% -15.2% -43.2% -13.2% % Relative Error – Soft priors -5.3% -6.6% -7.6% -11.8% % Relative Error – Hard priors -5.0% -6.6% -1.0% -22.2% (b) 177 Both soft and hard priors appear to have a similar impact on the recovered property values of the fibroglandular region (~70% and ~55% less errors in the reconstructed permittivity and conductivity profiles, respectively). However, they improve the reconstructed properties of the target inclusion differently. Specifically, the hard priors significantly improve the permittivity profile of the tumor region (decreasing the error to only 1%), while they have a negative effect on the counterpart conductivity values (the error increases to 22%). The soft priors, on the other hand, have a more stable effect on both permittivity and conductivity images, improving them by 42% and 10%, respectively. 3.4.3. Bone Imaging: Monitoring Changes in the Bone Dielectric Property Distributions Bone quality assessment is one of the areas where microwave imaging can be used to detect and monitor dielectric property changes. Bone dielectric properties have been studied extensively up to 5 MHz [132, 133]. These investigations have demonstrated strong correlations between the electrical and clinically important mechanical properties of bone. For instance, a recent study by Peyman et al showed that dramatic changes in bone dielectric properties occur in vivo with age (not typical in other tissues) and the age–related dynamics of bone physiology may be evident through interrogation of dielectric properties as a means of bone health monitoring [134]. In the current study, we will investigate the correlation between dielectric properties in the microwave frequency range and the bone volume fraction, which is an established measure of overall bone health. The results are supported by a series of ex–vivo bone imaging, as well as some initial clinical patient data from a recent pilot study for microwave imaging of the heel using a simple adaptation of our existing breast imaging system. 178 3.4.3.1. Ex–vivo Bone Imaging We performed a series of ex ex–vivo imaging on several de–marrowed marrowed porcine trabecular bone specimens. Figure 3.52 illustrates a flow diagram used for testing and processing the bone samples. Figure 3.52. Experiment flow diagram for evaluating the correlation between bone mineral and dielectric properties of porcine rcine trabecular bone specimens In each cycle, first a desiccator was used to remove the water and moist moisture ure from the bone pores. Then, a micro–CT CT scan was used to calculate the bone volume fraction (BVF) of each sample. Next, the specimens were immersed in 18 mm diameter test tubes containing a 0.9% saline solution to remove any air traces during the microwave measurements. As shown in Figure 3.53,, the test tubes were then placed at a position 3 cm offset from the center of the microwave imaging tank filled with a matching liquid comprised of an 80:20 mixture of glycerin:water solution.. 2D microwave data was acquired at a frequency range from 500 to 2500 MHz at a step size of 200 MHz.. Finally, the bone samples were de de–calcified calcified by being put into an 179 acid treatment procedure. At the end of each cycle, all samples were rinsed in tap water for about 10 minutes. Figure 3.53. Microwave data acquisition of the bone sample Since prior structural information of the bone samples (i.e. size and location of the test tubes) was available from the X–ray micro–CT, the soft prior regularization technique was used to obtain higher quality images and recovering more accurate dielectric property distributions. Figure 3.54 shows the customized soft prior reconstruction mesh used for recovering the dielectric properties of all bone samples. The mesh consists of 1107 nodes and 2092 triangular elements. 180 Figure 3.54. Soft prior microwave image reconstruction mesh comprised of 1107 nodes and 2092 triangular elements Figure 3.55 shows the representative 1100 MHz soft prior permittivity (top) and conductivity (bottom) images of the (a) first, (b) second, and (c) final (5th) microwave scans. (a) (b) (c) Figure 3.55. The 1100 MHz reconstructed permittivity (top) and conductivity (bottom) images of the (a) first, (b) second, and (c) final (5th) microwave scans 181 It is interesting to note that while the dry diameters for all samples decreased about 12.2% due to the acid treatments, when en the samples were immersed in the saline, they expanded slightly to snuggly fit up against the test tube walls. The properties of the recovered targets demonstrated a clear increasing trend from the first to the last imaging session. In order to evaluate the correlation between microwave dielectric properties and BVF, the peak recovered permittivity and conductivity values of the saline–saturated bone specimens in the test tube were extracted and are plotted as a function of bone volumee fraction in Figure 3.56(a) and (b) respectively. (a) (b) Figure 3.56. Relationship between the dielectric properties of saline–saturated saturated bones at 1100 MHz and BVF of dry bone specimens In both Figure 3.56(a) and (b) (b),, as the BVF increases, the dielectric property values decrease and the corresponding R2 correlation coefficients are 0.351 and 0.544, respectively. In terms of testing for statistical significance, we applied the Pearson correlation test [135],, which produced coefficient values of – 0.592 for permitti permittivity and -0.738 0.738 for conductivity. For the case of N = 30, absolute values of the Pearson coefficient greater than 0.373 are statistically significant to a p– p 182 value of 0.05 or better. As reported in Table 3.13, the results from a similar analysis performed on the 900 and 1300 MHz recovered dielectric properties were also statistically significant. Table 3.13. Pearson coefficient for the recovered dielectric properties at various frequencies Frequency (MHz) Permittivity Conductivity Pearson r–value Pearson r–value 900 1100 – 0.526 – 0.592 – 0.679 – 0.738 1300 – 0.603 – 0.753 These results suggest a significant correlation between the tissue dielectric properties (especially conductivity) and the bone volume fraction at our microwave frequency range, which agree with the strong relationship between bone content and the dielectric properties at lower frequencies, previously found by Sierpowska et al [132]. 3.4.3.2. Initial Clinical Heel Imaging We recruited patients who had required external support of one leg for approximately 4 – 6 weeks. Significant trabecular bone loss can be expected in the heel of a non–weight bearing leg, especially for the immediate period after cast or support removal. Studies have shown that trabecular-rich bones such as those in the ankle and calcaneus undergo dramatic changes in mineralization levels – decreases as much as 35% – when patients reduce weight-bearing for six weeks or more [136-138]. Figure 3.57 shows the CT cross sections through the calcaneous of both heels of a patient (#902) who wore a boot on one leg for at least 6 weeks (while using crutches) and recovered for several weeks (2 – 4) without support prior to imaging. According to the radiologist, the nominal Hounsfield units for the affected leg were 36.5 while those for the 183 normal weight–bearing leg were 49.2, indicating a 25.8% drop in bone density during the non– or limited weight–bearing time. (a) (b) Figure 3.57. X–ray CT images of the patient’s heels for the (a) normal weight–bearing and (b) injured leg (White zones within the calcaneous bones are the radiologists’ labels.) Since prior anatomical structure of the heels was available from X–ray CT images, customized reconstruction meshes were created (Figure 3.58) and used in the soft prior regularization technique. 184 (a) (b) Figure 3.58. 2D customized reconstruction meshes of the (a) normal weight–bearing and (b) affected leg comprised of 1196 nodes and 2264 triangular elements, and 815 nodes and 1519 triangular elements, respectively Figure 3.59(a) and (b) show the 2D recovered soft prior permittivity (top) and conductivity (bottom) images for the normal (fully weight–bearing) and injured (non– or limited weight– bearing) heels, respectively, suggesting slightly elevated values for the affected calcaneous. In this case, the recovered dielectric properties in the normal calcaneous were εr = 13.5 and σ = 0.81 S/m; while the values in the injured calcaneous were εr = 16.7 and σ = 0.92 S/m. Since bone has significantly lower dielectric properties comparing to the surrounding soft tissue, major trabecular bone loss in the heel of a non–weight bearing leg (i.e. injured heel) can result in higher bone dielectric properties. The findings are consistent with expectations in that the permittivity and conductivity values were about 24% and 14% higher for the injured calcaneous, respectively. 185 (a) (b) Figure 3.59. 1300 MHz soft prior permittivity (top) and conductivity (bottom) images for the heel of a (a) normal and an (b) injured leg utilizing the CT images as prior information during the reconstruction process 3.4.4. Clinical Breast Microwave Imaging in MR In order to assess the performance of our combined MI-MRI system, we imaged a 43 year old woman who had heterogeneously dense breasts, a body mass index of 25 and did not have any known cancers. Figure 3.60(a) shows the patient lying prone with her left breast suspended through an aperture in the top of the illumination tank right before entering the MR bore. Figure 3.60(b) shows the patient inside the bore. 186 (a) (b) Figure 3.60. Photograph of a patient in the combined MI-MRI system (a) before entering the MR bore and (b) inside the bore For this pilot test case, we took care to choose a woman with fairly low body mass index while also having modestly sized breasts to allow for sufficient tissue to be pendant down to the plane of the antennas. In addition, while previously categorized as having dense breasts, we deliberately reviewed earlier normal patient MR images and selected this subject based on there being a fair amount of distinct fibroglandular regions which were useful in testing the soft prior regularization scheme. It should be noted that we took particular care to not inadvertently create any conductive loops (i.e. from the various feed cables under the patient) that might have caused localized heating. Figure 3.61(a) shows the anatomically coronal MR T1-weighted image at the plane of the active part of the monopole antennas. The antennas are clearly visible and can be accurately registered with respect to the breast tissue. As discussed in Section 3.3.3, there is no noticeable distortion in the breast portion of the images due to the proximity of metal objects. After 187 segmenting the adipose and fibroglandular tissues in the MR image, the data was analyzed and co-registered with the corresponding microwave images. Then, our in-house developed codes were used to create a customized reconstruction mesh for use in the soft prior regularization algorithm. The resulting mesh shown in Figure 3.61(a) was comprised of 2473 nodes and 4647 triangular elements. For the fibroglandular tissue, there are three discrete zones (one larger and two smaller), which we tied together as a single zone for the purposes of the soft prior reconstruction. (a) (b) Figure 3.61 (a) MR T1-weighted image of the patient’s left breast at the plane of the active part of the monopole antennas, and (b) the corresponding segmented soft prior mesh Figure 3.62 shows the 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images with the no priors and soft prior regularization, respectively. For the former, the breast perimeter is roughly detected in both permittivity and conductivity (largely based on the fact that the breast properties are generally lower than those of the background), but with 188 considerable artifacts in and outside of the breast. The soft prior images are essentially homogeneous for the two zone types in both property images. (a) (b) Figure 3.62. 1300 MHz reconstructed permittivity (top) and conductivity (bottom) images of the patient’s left breast utilizing the (a) no priors and (b) soft prior regularizations, respectively While we do not have exact validation of the tissue properties for this patient, the recovered values are reasonable compared with published ex-vivo values for the associated tissue types [139]. In this instance, the recovered adipose property values were εr, Adp = 8.8 and σAdp = 0.09 S/m, while those for the fibroglandular tissue were εr, Fg = 40.7 and σFg = 0.96 S/m. Comparable published values for these tissue types were: εr, Adp = 4.8 and σAdp = 0.05 S/m; εr, Fg = 39.4 and 189 σFg = 0.75 S/m, respectively (assuming roughly 50% adipose content in the fibroglandular range) [139]. The patient results from our combined MI-MRI are encouraging, because they demonstrate that it is possible to collect both MR and microwave data on a patient even in such a confined space along with the various integration challenges. These first images are promising in that they demonstrate that the scattered signals contain information about the target that indicates which subzones have lower and higher properties. This will be particularly valuable when imaging patients in later studies where this approach may add valuable diagnostic information concerning the specific nature of suspicious lesions. 190 3.5. Results – 3D As discussed in Section 3.4, incorporating structural information of the object being imaged into our 2D image reconstruction algorithm through the soft prior regularization approach significantly improved the quality of the recovered images, and as a result, a more accurate reconstruction of dielectric property distributions was obtained. In this section, we will extend the idea of using prior spatial information with microwave imaging to our 3D image reconstruction algorithm. In addition to a complete analysis of the behavior of the soft prior technique in 3D, we will perform a comprehensive comparison between the reconstructed 2D and 3D images using the prior spatial information. 3.5.1. Simulation Experiments In order to study the effects of different factors such as the soft prior coefficient, reconstruction mesh density, imperfect spatial priors, and size, shape, and number of target inclusions on the 3D soft prior technique, several 3D simulation experiments will be performed. The schematic setup of the simulation experiments studied in this section is the same as the one described in Section 2.4. 3.5.1.1. No Priors versus Soft Prior Regularization Figure 3.63(a) and (b) show the 1300 MHz 3D reconstructed images of the synthetic data with and without incorporating prior spatial information, respectively. For visualization purposes, horizontal slices normal to the XY plane at z = 0, along with iso-surfaces of the target inclusion are presented in Figure 3.63. 191 (a) (b) Figure 3.63. 1300 MHz 3D reconstructed permittivity (top) and conductivity (bottom) images of the simulation experiment with – 100 dBm added noise using (a) the soft prior regularization and (b) no prior spatial information While the target inclusion is successfully detected in both cases, the soft prior reconstructed dielectric properties are much closer to the actual property values in both permittivity and conductivity images. Moreover, the level of background artifacts is significantly reduced in the images with incorporated structural information. In order to compare quantitatively the no priors and soft prior images, the relative root mean square error (RRMSE) for the inclusion region is calculated as Vnrecon −Vnexact RRMSE = ∑ exact V n n=1 N 2 N (3.12) 192 where N is the total number of nodes in the inclusion region, Vn recon is the reconstructed dielectric property value (either permittivity or conductivity) at node n (in the reconstruction mesh), and Vn exact is the true value of the selected dielectric property at that location. The permittivity and conductivity RMS errors associated with the inclusion region for each image in Figure 3.63 are summarized in Table 3.14. Table 3.14. The permittivity and conductivity RMS errors associated with the inclusion region in Figure 3.63 Soft priors No priors Permittivity RRMS Error 0.8% 22% Conductivity RRMS Error 8% 18% The soft prior regularization significantly improves the recovered dielectric property distributions, especially the permittivity values. In fact the permittivity RRMS errors are about 27 times smaller when the soft prior regularization is used. This improvement is less significant in the counterpart conductivity profiles, as the soft prior RRMS errors are only about two times smaller than that of the no prior case. 3.5.1.2. Number of the Reconstruction Nodes In order to study how the number of nodes in the customized reconstruction mesh affects the soft prior recovered dielectric properties, a set of customized meshes with different node densities was used to reconstruct the synthetic data in Section 2.4 with the soft prior regularization at 1300 MHz. A list of these meshes (including the number of nodes and elements) 193 is reported in Table 3.15. In order to minimize the effects of other factors, all parameters except the reconstruction mesh were kept the same in this experiment. Table 3.15. Customized reconstruction meshes used for the simulation experiment with soft prior regularization Case Number of Nodes Number of Elements (a) (b) (c) 2180 3083 4509 11033 16257 21330 (d) 7972 40686 In all cases, both permittivity and conductivity profiles of the target inclusion were effectively characterized, however, they did not appear to be significantly different as the number of reconstruction nodes increased. In order to analyze this observation quantitatively, the reconstructed dielectric properties of the inclusion and the background region were extracted at the center of the inclusion (3, 0, 0 cm) and at the center of the imaging domain (0, 0, 0 cm), respectively. Figure 3.64 shows the extracted property values along with the exact solutions, as a function of the number of reconstruction nodes in Table 3.15. 194 Figure 3.64. Extracted 3D recovered permittivity (left) and conductivity (right) values along with the exact property values, as a function of the number of nodes in the soft prior reconstruction mesh Results in Figure 3.64 confirm that the recovered soft prior dielectric properties are not highly sensitive to the reconstruction mesh density, especially in the permittivity images, as the recovered properties overlap with the exact values in all four cases. Nonetheless, the recovered conductivity values improve marginally, as the reconstruction mesh density increases from 2000 to over 4000 nodes, suggesting that for a configuration of this size, it appears to be sufficient to use a mesh composed of approximately 3000 – 4000 nodes. 3.5.1.3. Choice of the Soft-Prior Coefficient In the previous 3D soft prior reconstructions, the soft prior coefficient, λ in equation 3.4, was set to 10 as a default. However, in this section a more detailed study of the effects of λ on the recovered dielectric properties is presented. A spectrum of the soft prior coefficients (λ = 0.01, 0.1, 1, 10, and 100) was used for reconstructing the synthetic data described in Section 2.4. Lower values, such as 0.001 for λ, allowed the solution to diverge in some cases, while higher 195 values tended to suppress the recovered inclusion properties. Figure 3.65 3 shows the corresponding soft prior reconstructed images at 1300 MHz. (a) λ = 0.01 (b) λ = 0.1 196 (c) λ = 1.0 (d) λ = 10 (e) λ = 100 197 Figure 3.65. 1300 MHz 3D reconstructed permittivity (top) and conductivity (bottom) profiles for the simulation experiment using different soft prior coefficients λ = (a) 0.01, (b) 0.1, (c) 1.0, (d) 10, and (e) 100 In order to analyze the results quantitatively, the reconstructed dielectric properties of the inclusion and the background region were extracted at the center of the inclusion (3, 0, 0 cm), and at the center of the imaging domain (0, 0, 0 cm), respectively. Figure 3.66 shows the extracted property values along with the exact solutions, as a function of the soft prior coefficient λ. Figure 3.66. Extracted 3D recovered permittivity (left) and conductivity (right) values along with the exact property distributions, as a function of the base – 10 logarithm of the soft prior coefficient λ As the soft prior coefficient decreases, less regularization is applied and as a result, the overall quality of the reconstructed images degrades. This effect is even more prominent in the conductivity images, where the recovered properties of the target inclusion are significantly lower than the exact solution, and also closer to those of the coupling medium. On the other end, higher regularization values may cause excessive smoothing effects on the reconstructed images, 198 and as a result, the recovered properties of the inclusion are underestimated. For the intermediate values, no observable difference can be seen between the reconstructed values when λ = 1.0 to λ = 100 are used. When λ = 10, both recovered permittivity and conductivity values in the inclusion region are the closest to the exact solution, and those of the background match the true levels. The corresponding permittivity and conductivity RRMS errors as a function of the soft prior regularization coefficient λ are plotted in Figure 3.67. Figure 3.67. Relative RMS errors of the recovered dielectric property distributions as a function of the base – 10 logarithm of the soft prior regularization coefficient λ The relative RMS error plot confirms that the recovered soft prior permittivity values are in general more accurate than the counterpart conductivity values. Moreover, the minimum error is achieved when the soft prior coefficient of λ = 10 is used. 199 3.5.1.4. Sensitivity to Imperfect Spatial Priors Obtaining perfect a priori structural information of the object being imaged may not be possible in practice, therefore, it is essential to evaluate the sensitivity of the soft prior technique to imperfect priors. In this section, we extend the analysis performed in Section 3.4.2.8 to the 3D soft prior algorithm, and study two types of priors’ imperfections (the size and location of the target inclusion) in the simulation experiment described in Section 2.4. a) Imperfect prior size of the target inclusion: Table 3.16 presents some attributes of a set of different inclusion-sized 3D soft prior meshes that was used to reconstruct the synthetic data with a spherical target inclusion of 1.5 cm radius centered at (3, 0, 0 cm). In all reconstruction meshes, the total number of nodes, as well as the number of nodes in the inclusion region, was kept within approximately the same order of magnitude. Table 3.16. Characteristics of the 3D soft prior customized meshes used for a simulation experiment with imperfect prior size of the target inclusion Inclusion Location (In the actual experiment and in the soft prior mesh) Center (cm) Actual Prior Inclusion Size Inclusion Size Total Number of (In the soft prior (In the number of nodes in the mesh) experiment) nodes inclusion Radius (cm) Radius (cm) (3, 0, 0 cm) (3, 0, 0 cm) (3, 0, 0 cm) 1.5 1.5 1.5 0.7 0.8 0.9 3255 3209 3431 252 284 333 (3, 0, 0 cm) (3, 0, 0 cm) (3, 0, 0 cm) 1.5 1.5 1.5 1.0 1.1 1.2 3078 3358 3257 332 149 170 (3, 0, 0 cm) (3, 0, 0 cm) (3, 0, 0 cm) (3, 0, 0 cm) 1.5 1.5 1.5 1.5 1.3 1.4 1.5 1.6 3323 3249 3083 3007 182 198 240 241 200 (3, 0, 0 cm) (3, 0, 0 cm) 1.5 1.5 1.7 1.8 3693 2808 276 339 (3, 0, 0 cm) (3, 0, 0 cm) (3, 0, 0 cm) (3, 0, 0 cm) 1.5 1.5 1.5 1.5 1.9 2.0 2.1 2.2 2932 2945 3059 3016 358 438 504 546 (3, 0, 0 cm) (3, 0, 0 cm) (3, 0, 0 cm) 1.5 1.5 1.5 2.3 2.4 2.5 3213 2898 2815 608 613 627 In order to compare the recovered dielectric properties with the exact solution, the reconstructed permittivity and conductivity values of the background and the target inclusion were extracted at the center of each region, and are plotted in Figure 3.68, as a function of the spherical inclusion radius in the reconstruction mesh (i.e. the radius of the inclusion in the prior information). Figure 3.68. Extracted soft prior 3D recovered permittivity (left) and conductivity (right) values along with the exact property distributions, using different inclusion–sized meshes as prior information 201 In Figure 3.68, the exact size of the target inclusion is marked with a green vertical line, which indicates that when the exact prior information is used, the most accurate property distributions are recovered. More specifically, as the radius of the inclusion in the reconstruction mesh increases, both reconstructed permittivity and conductivity values of the target inclusion gradually approach those of the background medium (Bk). On the other end, as the prior size of the inclusion decreases, the reconstructed conductivity values also approach those of the background medium, while the permittivity counterparts first increase and then decreased toward the background values. In all cases, the recovered and exact background values overlap very closely. In terms of the sensitivity of the algorithm to the imperfect prior size of the target, there is about a 31% change in the recovered permittivity and 15% change in the recovered conductivity values, when approximately 126% change is applied to the volume of the target region (changing the inclusion radius from 1.3 to 1.7 cm). Comparing these numbers with those reported in 2D soft prior case in Section 3.4.2.8(a), reveals that the 3D soft prior algorithm is less sensitive to the imperfect prior size of the target. b) Imperfect prior location of the target inclusion: Table 3.17 presents some attributes of a set of different inclusion–located soft prior meshes that was used to reconstruct the synthetic data with a spherical target inclusion of 1.5 cm radius centered at (3, 0, 0 cm). In all reconstruction meshes, the total number of nodes, as well as the number of nodes in the inclusion region, was kept approximately the same. 202 Table 3.17. Characteristics of the 3D soft prior customized meshes used for a simulation experiment with imperfect prior location of the target inclusion Inclusion Size Inclusion Location (In the actual (In the experiment and in experiment) the soft prior mesh) Center (cm) Radius (cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) 1.5 (3, 0, 0 cm) Prior Inclusion Location (In the soft prior mesh) Center (cm) (0, 0, 0 cm) (0.2, 0, 0 cm) (0.4, 0, 0 cm) (0.6, 0, 0 cm) (0.8, 0, 0 cm) (1.0, 0, 0 cm) (1.2, 0, 0 cm) (1.4, 0, 0 cm) (1.6, 0, 0 cm) (1.8, 0, 0 cm) (2.0, 0, 0 cm) (2.2, 0, 0 cm) (2.4, 0, 0 cm) (2.6, 0, 0 cm) (2.8, 0, 0 cm) (3.0, 0, 0 cm) (3.2, 0, 0 cm) (3.4, 0, 0 cm) (3.6, 0, 0 cm) (3.8, 0, 0 cm) (4.0, 0, 0 cm) (4.2, 0, 0 cm) (4.4, 0, 0 cm) (4.6, 0, 0 cm) (4.8, 0, 0 cm) (5.0, 0, 0 cm) (5.2, 0, 0 cm) (5.4, 0, 0 cm) 203 Number of Total Nodes in Number of the Nodes Inclusion Region 2925 226 3041 238 3053 233 3024 241 3036 243 3005 239 3016 239 3087 235 3060 235 3007 241 3120 233 3157 241 3118 233 3094 240 3068 237 3083 240 2936 230 3044 230 2925 231 3046 233 2954 227 3267 235 3274 263 3204 253 3134 267 3047 253 3033 263 2969 260 In order to compare the recovered dielectric properties with the exact solution, the reconstructed permittivity and conductivity values of the background and the target inclusion were extracted at the center of each region. They are plotted in Figure 3.69, as a function of the X-coordinates of the center of the prior spherical inclusion in the reconstruction mesh (i.e. the center of the inclusion in the prior information). Figure 3.69. Extracted soft prior 3D recovered permittivity (left) and conductivity (right) values along with the exact property distributions, using different inclusion–located meshes as prior information In Figure 3.69, the exact location of the target inclusion is marked with a green vertical line. Overall, the imperfect prior location of the target seems to have a similar effect (hill-shaped curve) on both recovered permittivity and conductivity profiles. When the prior and actual target locations match, the reconstructed permittivity and conductivity values of the inclusion approach their exact profiles. However, the most accurate values are not recovered when the mismatch is zero, but when the prior inclusion is slightly offset along the X–axis. Nonetheless, the offset directions for the permittivity and conductivity profiles are the opposite of each other, which makes the recovered values at the actual target location (green line) the most accurate 204 simultaneous solutions for both permittivity and conductivity. When the partial overlap between the actual and prior target volumes decreases, the recovered property distributions of the inclusion approach those of the background medium. In all cases, the recovered and exact background values overlap very closely and only minor variations occur when imperfect prior information is used. However, the most accurate values are obtained at the exact prior location of the inclusion. In terms of the sensitivity of the algorithm to the imperfect prior location of the target, there is only about a 4% change in the recovered permittivity and 15% change in the recovered conductivity values when the prior target region is shifted approximately 87% of its radius (i.e. shifting the X-coordinates of the center of the inclusion from 2.6 to 3.4 cm). The results presented in this section verify that the soft prior algorithm is not overly sensitive to the exact size and location of the target region, and that it can tolerate priors’ imperfection to some extent. 3.5.1.5. Sensitivity to a False Region of Interest In order to study the sensitivity of the 3D soft prior algorithm to false regions of interest, the synthetic data described in 2.4 was used. In the actual experiment, there was only one spherical inclusion of 1.5 cm radius centered at (3, 0, 0 cm). However, as a part of prior information, an additional false region of interest of the same size was identified in the soft prior reconstruction mesh. The center of the false inclusion was set to 5 random locations. Some attributes of these meshes are presented in Table 3.18. In all reconstruction meshes, the total number of nodes, as well as the number of nodes in the true and false inclusion regions, was kept approximately the same. 205 Table 3.18. Characteristics of the 3D soft prior customized meshes used for a simulation experiment with a false target inclusion (Bk = background, TI = True Inclusion, and FI = False Inclusion) Case Number of Elements Number of Nodes Number of Bk Nodes Number of TI Nodes Number of FI Nodes TI Center (cm) FI Center (cm) (a) (b) (c) 19496 19494 18834 3687 3687 3583 3293 3292 3173 204 205 204 190 190 206 (3, 0, 0) (3, 0, 0) (3, 0, 0) (-1, -1, 0) (-2, 2, 1.2) (1, -4.5, -1) (d) 18353 3490 3081 206 203 (3, 0, 0) (1, 4.5, -1.2) (e) 18480 3514 3117 206 191 (3, 0, 0) (-5, 0, -1) Figure 3.70(a) – (e) show the 1300 MHz 3D reconstructed images using the corresponding soft prior meshes with a false region of interest in Table 3.18. For illustration purposes, the reconstructed images in each case are sliced in such a way that both the true and false inclusion regions can be seen. Moreover, iso-surface thresholds of εr = 39 and σ = 1.85 S/m are applied to the reconstructed permittivity and conductivity images, respectively. 206 (a) FI centered at (-1, -1, 0) (b) FI centered at (-2, 2, 1.2) (c) FI centered at (1, -4.5, 4.5, -1) (d) FI centered at (1, 4.5, -1.2) 207 (e) FI centered at (-5, 0, -1) Figure 3.70. 1300 MHz 3D reconstructed images using the corresponding soft prior meshes with a false region of interest in Table 3.18 The false region of interest does not appear prominently in any of the images in Figure 3.70. It onlyy presents as a subtle variation from the background properties. In order to study the sensitivity of the 3D soft prior algorithm to a false region of interest, the recovered dielectric properties of all 5 cases were extracted at the center of the true and false inclusions, and they are presented in Table 3.19,, along with their corresponding error values. 208 Table 3.19. Recovered dielectric properties (εr and σ) of the true (TI) and false inclusion (FI) extracted at the center of each region Case εr, TI %Error σTI %Error (a) 39.40 2% 1.89 (b) (c) (d) (e) 39.44 39.51 39.32 39.53 1% 1% 2% 1% 1.90 1.89 1.89 1.89 εr, TI εr, FI %Error 5% 22.51 5% 5% 6% 6% 22.38 22.48 22.46 22.29 σTI σFI %Error 0.5% 1.28 4% 0.1% 0.4% 0.3% 0.5% 1.25 1.26 1.26 1.25 2% 3% 3% 1% εr, FI σFI From Table 3.19, for all cases the false region appears with a very small increase or decrease (less than 1%) in the permittivity images, but appears with a slightly larger increase (from 1% to 4%) in the conductivity images. In the true inclusion region, while the permittivity errors are less than 2%, the counterpart conductivity errors are somewhat larger, but still less than 6%. The results in this section show that the 3D soft prior algorithm is able to manage situations where a false region of interest has been identified in the prior structural information. In such situations the false region of interest only shows a weak variation of the background properties. 3.5.1.6. Multiple Inclusions and Number of Regions In order to study the effects of the number of inclusions as well as the number of nodes and regions in the reconstruction mesh on the recovered property values, when the 3D soft prior algorithm is used, the simulation experiment described in Section 2.4 was repeated, but this time there were three spherical inclusions of 1.5 cm radius. The target inclusions I1, I2, and I3 were centered at (x, y, z)I1 = (3, 0, 0 cm), (x, y, z)I2 = (– 2, 3, 0 cm), and (x, y, z)I3 = (– 2, – 3, 0 cm), 209 respectively. Similar to the study performed on the 2D soft prior algorithm in Section 3.4.1.4, the following cases were considered: a) Multiple Inclusions with Identical Dielectric Properties: In the synthetic data, the dielectric properties of all inclusions were set to εr, I1 = εr, I2 = εr, I3 = 40.0 and σI1 = σI2 = σI3 = 2.0 S/m. The 3D soft prior images were then reconstructed using four different customized meshes. The first two meshes were composed of 5545 nodes and 27389 tetrahedral elements, whereas the other two meshes were composed of 2459 nodes and 12656 tetrahedral elements. In the first and third mesh, all inclusion regions were assigned the same region number (i.e. RI1 = RI2 = RI3); therefore, only two distinct regions were considered: one for the background (Bk) and one for all inclusions (I). In the second and forth meshes, each inclusion region was assigned a different region number (i.e. RI1, RI2, and RI3 had independent properties); therefore, four distinct regions were considered in this case: one for the background (Bk) and three for each inclusion (I1, I2, and I3). A summary of the soft prior reconstruction meshes used for this simulation experiment is presented in Table 3.20. Table 3.20. Soft prior reconstruction meshes used for the simulation experiment with three spherical target inclusions Mesh5545_2r Mesh5545_4r Mesh2459_2r Number of Nodes 5545 5545 2459 Number of Elements 27389 27389 12656 Number of Distinct Regions 2 4 2 Mesh2459_4r 2459 12656 4 3D Mesh ID 210 Figure 3.71 shows the 3D soft prior reconstructed dielectric properties extracted at the center of each region, along with their exact property values values, using four customized reconstruction meshes in Table 33.20. (a) Permittivity (b) Conductivity Figure 3.71. 3D soft prior recovered (a) permittivity and (b) conductivity values extracted at the center of each region, along with the exact properties properties, using customized reconstruction meshes in Table 3.20 Overall, all recovered property distributions are fairly close to the exact solution. When 2region meshes (Mesh5545_2r and Mesh2459_2r)) are used, the recovered values of all targets are 211 the same. On the other hand, when 4-region meshes (Mesh5545_4r and Mesh2459_4r) are used, some variations appear between the recovered properties of the inclusions, although they were assigned the same property values. In these cases, the deviation from the exact solution is more prominent in the recovered conductivity profiles. A complete summary of the relative dielectric property errors using different reconstruction meshes with two and four distinct regions is shown in Table 3.21. Table 3.21. Relative (a) permittivity and (b) conductivity errors for the multi-inclusion simulation experiment using different reconstruction meshes with two and four distinct regions. The negative sign indicates that the reconstructed property value is underestimated. (i.e. less than the exact property values) Rel. Err. Rel. Err. Rel. Err. Rel. Err. εr, BK εr, I1 εr, I2 εr, I3 Mesh5545_2r Mesh5545_4r Mesh2459_2r -0.6% -0.6% -0.9% 0.1% -0.2% -2.2% 0.1% -0.2% -2.2% 0.1% 0.2% -2.2% Mesh2459_4r -0.9% -1.7% -2.5% -3.0% Rel. Err. Rel. Err. Rel. Err. Rel. Err. σ BK σ I1 σ I2 σ I3 2.5% 2.5% 2.4% 2.4% -2.1% -4.0% -4.7% -6.3% -2.1% -2.2% -4.7% -4.8% -2.2% -0.8% -4.8% -4.4% Reconstruction Mesh (a) Reconstruction Mesh Mesh5545_2r Mesh5545_4r Mesh2459_2r Mesh2459_4r (b) When a finer soft prior mesh (5542 nodes) is used, the reconstructed properties are overall closer to the exact solution and both relative permittivity and conductivity errors decrease in 212 almost all regions. Moreover, it seems that restricting variations between different regions that have the same properties (i.e. using meshes with two distinct regions instead of four) can improve the accuracy of the soft prior reconstructed dielectric properties. b) Multiple Inclusions with Different Dielectric Properties: This case was identical to the previous experiment, except for the property values of the target inclusions. In this case, dielectric properties of εr, I1 = 27.0 and σI1 = 1.45 S/m, εr, I2 = 35.0 and σI2 = 1.66 S/m, and εr, I3 = 40.0 and σI3 = 2.0 S/m were assigned to I1, I2, and I3, respectively. Since the dielectric properties of the inclusions were different, only customized meshes with four distinct regions (Mesh5545_4r and Mesh2459_4r Mesh2459_4r) were used for the soft prior reconstruction. Figure 3.72 shows the 3D soft prior reconstructed dielectric properties extracted at the center of each region, region along with their exact property values. (a) 213 (b) Figure 3.72. Soft prior 3D recovered (a) permittivity and (b) conductivity values extracted at the center of each region, along with the exact properties properties, using 4-region customized reconstruction meshes in Table 3.20 Similar to the previous experiment’s results, both permittivity and conductivity values of the background region are recovered very accurately, and all inclusions are characterized successfully. Moreover, the correct pattern of property increase from I1 to I3 can be seen in the soft prior reconstructed permittivity and conductivity profiles. A complete error analysis of the reconstructed soft prior property values is reported in Table 3.22. The negative egative sign indicates that the reconstructed property roperty value is underestimated. (i.e. less than the exact property values) 214 Table 3.22. Relative (a) permittivity and (b) conductivity errors for the multi-inclusion simulation experiment with different property contrasts using two soft prior reconstruction meshes with 4 distinct regions Rel. Err. Rel. Err. Rel. Err. Rel. Err. εr, BK εr, I1 εr, I2 εr, I3 – 0.5% – 0.6% 0.6% 0.3% 0.1% – 1.6% – 0.1% – 3.4% Rel. Err. Rel. Err. Rel. Err. Rel. Err. σ BK σ I1 σ I2 σ I3 Mesh5545_4r 2.5% 0.2% – 0.7% – 1.2% Mesh2459_4r 2.4% – 1.0% – 2.4% – 4.6% Reconstruction Meshes Mesh5545_4r Mesh2459_4r (a) Reconstruction Mesh (b) In general, when a finer soft prior mesh (5545 nodes) is used, the relative permittivity and conductivity errors decrease. This effect is more prominent on larger errors as the dielectric property contrasts between the background medium and target inclusions increase from I1 to I3. In the simulation experiments described in Section 2.4, the target inclusion was centered at the same location as I1 – i.e. (x, y, z)I1 = (3, 0. 0 cm), and the dielectric properties were identical and equal to εr, I1 = 40.0 and σI1 = 2.0 S/m. Therefore, in order to study the effects of additional regions of interest (in this case I2 and I3) on the recovered dielectric properties of I1, the synthetic data in 2.4 was reconstructed, using the 5545 node mesh with only one inclusion region I1 (I2 and I3 regions in the reconstruction mesh were assigned the same region number as that of the background medium). The recovered soft prior property values of I1 along with the exact solution for the single and multi-inclusion cases are plotted in Figure 3.73. 215 (a) (b) Figure 3.73. Comparison between the soft prior reconstructed (a) permittivity and (b) conductivity values of I1 with and witho without additional targets (i.e. multi-inclusion inclusion and single inclusion cases, respectively) The difference in recovered dielectric property values in Figure 3.73 is minimal (0.5% (0.5 in permittivity and 3.5% in conductivity conductivity), which suggests that additional targets do not have a major effect on the 3D reconstructed soft prior dielectric property profiles. 3.5.1.7. Size and Shape of the Target In this section, we study the effects ooff the size and shape of the target inclusion on the recovered dielectric properties both with and without using prior structural information. We also compare the corresponding reconstructed images in 2D and 3D. Similar to the previous simulation experiments, synthetic 3D measurement data with –100 dBm of added noise was generated by the FDTD forward solver for different shaped and different sized target inclusions. More specifically, 48 monopole antennas were configured in three evenly–spaced circles with a 15.2 cm diameter and a 0.5 cm separation from each other. Five different sized spherical and cylindrical shaped inclusionss of 0.25, 0.5, 0.75, 1.0, and 1.25 cm radius, centered at (x, y, z)) = (3 (3, 0, 0 cm), with dielectric properties of εr, Inc = 40.0 and σInc = 216 2.0 S/m were embedded in a background medium with dielectric properties of εr, Bk = 22.4 and σBk = 1.23 S/m. For 3D reconstructions, the full-data selection containing both in-plane and cross-plane data was used, whereas for 2D reconstructions, only the middle in-plane data was selected. All no prior 2D and 3D images were reconstructed on a circular mesh composed of 559 uniformly distributed nodes and 1044 triangular elements, and on a cylindrical mesh comprised of 1298 uniformly distributed nodes and 5259 tetrahedral elements, respectively. On the other hand, the hard prior and soft prior images were reconstructed on customized 2D and 3D meshes summarized in Table 3.23. Table 3.23. Characteristics of the customized 2D and 3D meshes used for the simulation experiments with different sized and shaped target inclusions Inclusion Shape Cylinder Cylinder Inclusion Size 2D/3D Number of Number of Radius (cm) Mesh Nodes Elements 0.25 2D 803 1537 0.25 3D 3000 13863 Cylinder Cylinder Cylinder Cylinder Cylinder Cylinder 0.5 0.5 0.75 0.75 1.0 1.0 2D 3D 2D 3D 2D 3D 512 2160 421 1891 418 1726 962 9220 788 7850 780 7029 Cylinder Cylinder Sphere Sphere 1.25 1.25 0.25 0.25 2D 3D 2D 3D 429 1628 803 6382 800 6562 1537 35131 Sphere Sphere Sphere 0.5 0.5 0.75 2D 3D 2D 512 3493 421 962 18534 788 217 Sphere Sphere Sphere Sphere 0.75 1.0 1.0 1.25 3D 2D 3D 2D 1990 418 1556 429 9846 780 7428 800 Sphere 1.25 3D 1175 5305 For each case, the recovered dielectric property distributions were extracted at the center of the inclusion. Figure 3.74(a) and (b) show the reconstructed values of the cylindrical and spherical inclusions as a function of the inclusion size, using different reconstruction algorithms (i.e. 2D and 3D, with no priors, soft priors, and hard priors). 218 (a) (b) Figure 3.74. Recovered permittivity (top) and conductivity values at the center of the (a) cylindrical and (b) spherical inclusions as a function of the inclusion size, using different reconstruction approaches (2D and 3D, with no priors, soft priors, and hard priors) Overall, as the size of the inclusion decreases, the recovered dielectric properties become less accurate. In addition, the properties of the cylindrical inclusion are superior to those of the spherical inclusion, especially for the larger size targets. This improvement is even more prominent when prior structural information is incorporated into the reconstruction algorithm. 219 The no prior 3D reconstructed permittivity values (red dotted lines in the top parts of Figure 3.74) are, overall, similar to those reconstructed in 2D with and without prior information (all blue lines). However, the corresponding counterpart conductivity values (the bottom graphs in Figure 3.74) are significantly superior when the radius of the inclusion is greater than 0.5 cm, which confirms that the 3D reconstruction algorithm has a greater improvement on the recovered conductivity profiles. Comparing the two priors methods (i.e. hard priors and soft priors), it appears that the hard prior reconstructed properties are slightly superior to those reconstructed with the soft prior regularization. This effect is significantly more enhanced in the 3D reconstructed profiles (red dotted lines) of smaller sized spherical inclusions in Figure 3.74(b). This may be due to the fact that with the hard priors, the number of unknowns – i.e. the size of the reconstruction parameters – is drastically reduced from 418 – 803 in 2D and 1175 – 3000 in 3D to only two parameters for each of the two identified regions (BK and Target), and as a result, no regularization was required during the reconstruction. On the other hand, in the soft prior technique, the number of unknowns (i.e. the number of nodes in the customized reconstruction mesh) is significantly greater than the number of measurements, which requires additional regularizations to stabilize the image reconstruction. 3.5.2. Phantom Experiments In order to study the robustness of the 3D soft prior technique in more complex-shaped scenarios, and to evaluate the accuracy of the reconstructed dielectric property distributions in real measured data, two phantom experiments are performed. The MR images of the corresponding experiments are acquired separately, and the prior spatial information is used for the 3D soft prior algorithm. 220 3.5.2.1. Cylindrical and Cup-shaped Gelatin Inclusions In the first phantom experiment, a thin-walled plastic tube of 1.0 cm radius, filled with a 60:40 mixture of glycerin:water was inserted into a cup-shaped gelatin inclusion, as illustrated in Figure 3.75(a). The inclusions were then submerged in the imaging tank filled with an 86:14 mixture of glycerin:water, as illustrated in Figure 3.75(b). Figure 3.75(c) shows the top view of the phantom setup in the microwave imaging system. The independently measured dielectric properties of the coupling medium (Bk), cylindrical inclusion (CylInc), and the gelatin inclusion (GelInc) at 1300 MHz are reported in Table 3.24. (a) (b) (c) Figure 3.75. Complex-shaped phantom experiment with a cylindrical and a cup-shaped gelatin inclusion: (a) cylindrical inclusion (CylInc) inserted into the cup-shaped gelatin inclusion (GelInc), (b) both inclusions submerged in the imaging tank, and (c) top view of the phantom setup in the microwave imaging tank 221 Table 3.24. Independently measured dielectric properties of the background medium (Bk), cylindrical inclusion (CylInc), and cup-shaped gelatin inclusion (GelInc) at 1300 MHz Frequency (MHz) εr, Bk σBk 1300 10.8 0.75 εr, CylInc σCylInc εr, GelInc σGelInc 50.60 1.29 42.45 2.36 In the present phantom experiment, 3D microwave data was acquired at multiple planes. More specifically, the antenna array transmitted and received the signal at 8 equally–spaced (1 cm) vertical positions, and as a result, 8 × 8 × 240 measurements were collected. For the 3D reconstructions, the multi-plane (two consecutive planes) data selection was used. In order to obtain the corresponding structural information, the phantom inclusions were also imaged with the MRI. For co-registration purposes, the inclusions were placed in exactly the same location in an empty tank with an equal-sized antenna array, as shown in Figure 3.76(a) and (b). Figure 3.76(c) shows the top view of the phantom configuration in the MRI scanner. (a) (b) 222 (c) Figure 3.76. MRI Imaging of the phantom experiment with a cylindrical and a cup-shaped gelatin inclusions: (a – b) pictures of the phantom placed in an empty microwave imaging tank inside the MR bore, and (c) top view of the schematics of the phantom in the MRI scanner Figure 3.77(a) shows one of the corresponding MR images of the inclusions, which were post-processed and segmented to create the resultant 3D soft prior mesh in Figure 3.77(b). 223 (a) (b) Figure 3.77. (a) One of the corresponding MR images of the inclusions, (b) the corresponding 3D soft prior mesh composed of 3393 nodes and 17027 tetrahedral elements The customized soft prior mesh in Figure 3.77(b) is composed of 3393 nodes and 17027 tetrahedral elements. While this mesh was used to reconstruct microwave images with prior spatial information, the no prior images were reconstructed on a 6.9 cm radius, 10 cm height cylindrical mesh composed of 5301 uniformly distributed nodes and 23714 tetrahedral elements. Figure 3.78(a) and (b) show the 1300 MHz reconstructed permittivity (top) and conductivity (bottom) with and without prior spatial information, respectively. In order to view the recovered dielectric properties of the regions of interest, the reconstructed volumes in Figure 3.78 were sliced vertically through the inclusion regions. 224 (a) (b) Figure 3.78. A vertical slice through (a) soft prior and (b) no prior reconstructed permittivity (top) and conductivity (bottom) images of the phantom experiment with a cylindrical and a cupshaped gelatin inclusions at 1300 MHz When the soft prior regularization is used, the contrast between the inclusions is successfully detected in both permittivity and conductivity images. In fact, the higher permittivity values of the cylindrical inclusion versus the higher conductivity values of the gelatin inclusion are depicted in the recovered dielectric property profiles in Figure 3.78(a). However, when no prior structural information is used during the reconstruction, the recovered target inclusions are not distinguishable from each other. Moreover, while the cylindrical 225 inclusion appears indistinctly on the top-center of the reconstructed conductivity image, there are no indications of such region in the counterpart permittivity image in Figure 3.78(b). In order to analyze the recovered dielectric properties quantitatively, the soft prior and no prior RRMS errors of the recovered properties are calculated in both inclusion regions (i.e. Cylindrical – CylInc and Gelatin – GelInc inclusions), and the results are summarized in Table 3.25. Table 3.25. Soft-prior and no prior RRMS errors of the recovered properties in each region of the breast phantom experiment with two target inclusions Soft priors No priors εr,CylInc σCylInc εr,GelInc σGelInc εr,Total σTotal RRMS Error RRMS Error RRMS Error RRMS Error RRMS Error RRMS Error 0.07 0.55 0.41 0.22 0.15 0.41 0.05 0.37 0.14 0.43 0.15 0.35 The soft prior regularization significantly improves the recovered dielectric property distributions. In fact, the RRMS errors are overall about 3 to 8 times smaller when the soft prior regularization is used. This improvement is less significant in the conductivity profile of the cylindrical inclusion, as the soft prior RRMS error is only about 2 times smaller than that of the no prior case. 3.5.2.2. Breast Phantom with Two Target Inclusions In this experiment we assessed the 3D soft prior algorithm in a more realistic configuration, where two arbitrarily–shaped target inclusions were placed inside a prototyped plastic breast model, as shown in Figure 3.79(a). The breast model was filled with a 88:12 glycerin:water mixture (mimicking the dielectric properties of a scattered breast), and it was then submerged in 226 the imaging tank filled with a ma matching tching liquid composed of a 86:14 mixture of glycerin:water, as illustrated in Figure 3.79(b). Figure 3.79(c) shows the schematic configuration of the phantom. The independently measured dielectric properties of the coupling medium, breast model, and two target inclusions at 1100 MHz are reported in Table 3.26. (a) (b) (c) Figure 3.79. Breast–shaped shaped phantom experiment with two arbitrarily shaped target inclusions: (a) Two arbitrary-shaped gelatin target inclusion inclusions suspended in the plastic breast model, model (b) rapid– prototyped plastic breast model submerged in the imaging tank, and (c) schematic configuration of the breast phantom in the microwave im imaging tank. 227 Table 3.26. Independently measured dielectric properties of the background medium (Bk), breast model (Br), and the target inclusions (Inc1 and Inc2) at 1100 MHz Frequency (MHz) εr, Bk σBk εr, Br σBr 1100 12.2 0.71 14.2 0.77 εr, Inc1 σInc1 35.0 1.17 εr, Inc2 σInc2 35.0 1.17 In the present phantom experiment, 3D microwave data was acquired at multiple planes. More specifically, the antenna array transmitted and received the signal at 11 equally–spaced (1 cm) vertical positions, and as a result, 11 × 11 × 240 measurements were collected. For the 3D reconstructions, the multi-plane (two consecutive planes) data selection was used. In order to obtain the corresponding structural information, the breast model was also imaged with the MRI. For co-registration purposes, the breast model with the target inclusions was placed in exactly the same location in an empty tank with an equal-sized antenna array, as shown in Figure 3.80(a) and (b). Figure 3.80(c) shows the top view of the phantom configuration in the MRI scanner. For illustration purposes, an MRI coronal slice of the phantom containing both target inclusions is placed at the center of the imaging domain. (a) (b) 228 (c) Figure 3.80. MRI Imaging of the breast phantom experiment with two inclusions: (a – b) pictures of the phantom placed in an empty microwave imaging tank inside the MR bore, and (c) top view of the phantom configuration in the MRI scanner After the experiment was done, the MR images were post-processed and segmented. Figure 3.81(a) shows a stack of binary-segmented images used to create the corresponding 3D soft prior mesh in Figure 3.81(b). 229 (a) (b) Figure 3.81. (a) Stack of binary-segmented MRI images of the breast phantom experiment with two target inclusions, (b) the corresponding 3D soft prior mesh composed of 7540 nodes and 34591 tetrahedral elements The customized soft prior mesh in Figure 3.81(b) is composed of 7540 nodes and 34591 tetrahedral elements. While this mesh was used to reconstruct microwave images with prior spatial information, the no prior images were reconstructed on a 6.9 cm radius, 11 cm height cylindrical mesh composed of 5720 uniformly distributed nodes and 25800 tetrahedral elements. Figure 3.82(a) and (b) show the 1100 MHz reconstructed permittivity (top) and conductivity (bottom) with and without prior spatial information, respectively. In order to view the recovered dielectric properties of the regions of interest, the reconstructed volumes in Figure 3.82 were sliced vertically through the breast region. 230 (a) (b) Figure 3.82. A vertical slice through the (a) soft prior and (b) no prior reconstructed permittivity (top) and conductivity (bottom) images of the breast phantom experiment with two target inclusions at 1100 MHz When the soft prior regularization is used, the target inclusions are successfully characterized in both permittivity and conductivity images, although they appear weaker in the conductivity counterpart. When no prior spatial information is incorporated into the microwave reconstruction algorithm, the target inclusions are only successfully detected in the permittivity images. In the counterpart conductivity images, only the lower target (Inc2) appears to be detected. In addition, its size is significantly larger and its location is drastically shifted from the 231 actual location. In all cases, the recovered property distributions of the breast are closely matched to the exact values. In order to analyze the recovered dielectric properties of the phantom experiment quantitatively, the soft prior and no prior RRMS errors of the recovered properties were calculated in each region (i.e. Br, Inc1 and Inc2), and the results are summarized in Table 3.27. Table 3.27. Soft prior and no prior RRMS errors of the recovered properties in each region of the breast phantom experiment with two target inclusions εr, Br σBr εr, Inc1 σInc1 εr, Inc2 σInc2 εr, Total σTotal RRMS Error RRMS Error RRMS Error RRMS Error RRMS Error RRMS Error RRMS Error RRMS Error Soft priors 0.06 0.11 0.09 0.28 0.09 0.11 0.06 0.12 No priors 0.29 0.25 0.44 0.44 0.35 0.23 0.30 0.25 The soft prior regularization significantly improves the recovered dielectric property distributions, especially the permittivity values. In fact the permittivity RRMS errors are about 5 times smaller when the soft prior regularization is used. This improvement is less significant in the counterpart conductivity profiles, as the soft prior RRMS errors are only about 2 times smaller than that of the no prior case. In order to evaluate how the incorporated prior structural information improves the contrast between the target inclusions and the background breast region, the percentage of the contrast enhancement (CE) of the targets is calculated as: CEti = SP −100 × (CtNP − C t i i ) (3.13) CtNP i 232 where and are the no priors and soft priors contrasts of the ith target inclusion with respect to the true properties of the breast region, respectively, and they are calculated as NP ti C ∑ = N NP VnNP N n=1 ExactBR V and SP ti C ∑ = N SP V SP N n=1 n ExactBR (3.14) V where VnNP and VnSP are the no priors and soft priors reconstructed dielectric property values (either permittivity or conductivity) at node n (in the corresponding reconstruction meshes with ExactBR NNP and NSP nodes), respectively, whereas V is the true dielectric property value of the breast region. The computed contrast enhancements of the target inclusions with the incorporated prior structural information are presented in Table 3.28. Table 3.28. The computed contrast enhancements of the target inclusions with respect to the breast region, when the prior structural information of the phantom is used during the reconstruction % Contrast Enhancement (CE) εr σ Upper target (Inc1) Lower target (Inc2) 44.3 28.6 25.2 8.3 Positive values in Table 3.28 indicate that incorporating prior spatial information of the phantom enhances the contrast between the target inclusions and background breast region. Overall, this improvement is more prominent in the permittivity images, as the contrast enhancement is as high as 45% in Inc1. The corresponding contrast enhancements in the counterpart conductivity images are relatively smaller, especially in the lower inclusion (Inc2), 233 which shows that using prior structural information can improve the contrast between the lower target and its background by only about 8%. 234 4. 3D Microwave Image Reconstruction GUI In order to enable the user to interactively select the input data and the desired reconstruction parameters, we have developed a MATLAB Graphical User Interface (GUI) for 3D microwave image reconstruction. Prior to image reconstruction, the user selects the raw measurement and calibration (only the homogeneous bath) data acquired by our clinical microwave imaging system (MIST), as shown in Figure 4.1. This data can be in either 2D or 3D format. (a) (b) Figure 4.1. Input file selection: (a) file format and (b) raw measurement file selection windows 235 The information about the acquired data (i.e. patient ID, side, exam date, used frequencies, and antenna configuration including the number of planes for plates A and B in Figure 1.2(a) and the corresponding slice thicknesses) is automatically extracted from the header of the measured data file. Based on this information, the location of each antenna (either transmitting or receiving) in the 3D space is calculated and a plot containing the antenna configuration is displayed to the user. An example of such a plot is shown in Figure 4.2. Figure 4.2. Sample antenna configuration plot where 6 planes of data with a slice thickness of 1 cm have been acquired Once the measured and calibration data are subtracted from each other, the calibrated data is phased unwrapped and stored in two multi-dimensional matrices, one for the amplitude and one for the phase. At this stage, the data redundancy is handled as described in Section 2.3.2. Next, as illustrated in Figure 4.3, the user can select the output file format. It should be noted that the output file here is the input file for the reconstruction algorithm. That is to say, if the “2D” 236 output file format is selected, the calibrated measured data at each plane is stored in a file with the right format to be used in the 2D reconstruction algorithm. Similarly, the “3D” output file format selection corresponds to input data files for the 3D reconstruction algorithm. Figure 4.3. Output file format selection window Beside the file format, the antenna numbering system is different in the input files to the 2D and 3D reconstruction algorithms. More specifically, in the 2D case, the in-plane antennas on a circle are numbered sequentially (1 – 16), whereas in the 3D case, first the antennas on plate B (in Figure 1.2(a)) are numbered sequentially for each plane (plane 1: 1 – 8, plane 2: 9 – 16, etc.), and then, the same process is applied to those on plate A. The reason we have chosen this antenna numbering system is that, unlike the 2D case where all antennas on plates A and B move together, during the 3D data acquisition there could be a different number of planes with different slice thicknesses for plates A and B. Therefore, sequentially numbering the antennas based on the plates on which they are mounted is the most consistent way for the 3D antenna numbering system. If the “3D” output file format is selected, the user can decide what portion of the acquired data should be used for the reconstruction. As described in Section 2.3.3, this data selection is done through a matrix of checkboxes, where each checked box represents a flag corresponding to (TX, RX) pair (i, j), as shown in Figure 4.4. Diagonal and impossible pairs of the matrix are disabled so the user cannot switch them on. A few buttons, including SELECT ALL, CLEAR 237 ALL, SELECT IN–PLANE, and MULTI PLANES are also designed to help the user make selections faster and more efficiently. Based on the number of planes used to acquire measurements, the user can select a specific set of in-plane data (for example P1, P2, …) by checking/un–checking the corresponding checkboxes, in Figure 4.4. Figure 4.4. 3D data selection window The MULTI PLANES option allows the user to select multiple planes of data with a specific number of consecutive planes. For example, as shown in Figure 4.4, if 2 consecutive planes (2CS) are selected, for each transmitting antenna, all those in the same plane in addition to those on the other plate located one plane above and one plane below are automatically checked as receivers. Figure 4.5 shows an example of a multi-plane data selection with 2 consecutive planes. 238 Figure 4.5. Example of a multi-plane data selection with 2 consecutive planes Once the transmitter/receiver selection process is done, the user presses CONTINUE and proceeds to the next step, which is the 3D reconstruction input/parameter selection. A list of these inputs/reconstruction parameters along with a brief description of them is reported in Table 4.1. 239 Table 4.1. A list of inputs/3D reconstruction parameters along with a brief description of them Input/parameter ptID exam_date side Description Patient ID Exam (data acquisition) date ‘L’ or ‘R’ frequency Frequency list n_iteration Number of iterations fbgmedia meshID IsCoformalMesh IsProfiling IsParallel ParallelConfig IsSoftPrior SP_coef n_region IsHardPrior IsSimu SimuNoise PMLLayerNo Background medium property filename Mesh ID (the corresponding nodes and element files have the same ID) ‘y’ if a conformal mesh is used, ‘n’ otherwise: The mesh will be centered with respect to the center of the imaging domain (z = 0 is always the midplane ‘y’ if the user wants to use MATLAB profiling, ‘n’ otherwise ‘y’ if the user wants to use MATLAB parallel toolbox for the forward solution, ‘n’ otherwise If IsParallel = ‘y’, the user needs to select one of the preset parallel configurations (‘2nodes’, ‘4nodes’, etc.) ‘y’ if the user wants to use the soft prior regularization technique. In that case, the soft prior reconstruction mesh (with different regions) is required. ‘n’ otherwise If IsSoftPrior = ‘y’, the user needs to select a weighting factor as a soft prior coefficient. If IsSoftPrior = ‘y’, the user needs to specify how many different regions are in the soft prior reconstruction mesh ‘y’ if the user wants to use the hard prior technique. In that case, the soft prior reconstruction mesh (with different regions) is required. ‘n’ otherwise ‘y’ if reconstructing a simulation experiment, ‘n’ otherwise If IsSimu = ‘y’, the user needs to specify how much noise has been added to the synthetic data (the format is the |power level|). Example: 120, 100, or 0 for no added noise Number of PML layers (by default = 5) WaveLengthDiv Number of cells per wavelength (D in equation 2.8) TikhonovAlpha Tikhonov Regularization coefficient (α) Once the input parameters are selected, the 3D reconstruction algorithm begins to run in the background. The user can return later to access the 3D reconstructed images, which are conveniently in .vtk format. The open source scientific visualization software, ParaView 240 (Kitware, Inc., Clifton Park, New York), can be used to display the reconstructed images in 3D. In addition, the results are stored in a “Tecplot-ready” file for viewing or manipulation, using the Tecplot 360 (Tecplot, Inc., Bellevue, WA). 241 5. Summary and Conclusion This project focused on the computational aspect of tomographic microwave imaging for biomedical applications. Throughout this work, three main tasks were completed: 1) Our existing microwave 3D image reconstruction algorithm was optimized and redesigned to a faster, more practical, flexible, and also user–friendlier 3D image reconstruction module which is computationally feasible and balances the tradeoff between the accuracy and efficiency of the model. Through several simulation, phantom experiments, and clinical exams it was shown that using 3D microwave imaging can be beneficial in a number of ways including: a) The increased amount of data used during the reconstruction procedure can reveal more details about the internal structure of the imaged object/tissue, b) Despite the fact that our 2D reconstructed images have been proved to be very informative and clinically useful, the 3D reconstruction algorithm can help improve the recovered dielectric property distributions even further, and c) due to the greater data–model mismatch, the 2D reconstructed images (especially the conductivity) may show traces of the object in planes that extend beyond its actual location. Using a 3D model during the reconstruction process can improve the accuracy of the images in terms of both the location of the target and its recovered dielectric properties. In addition, it was shown that using different data selections (full-data, multi-plane, and inplane) in the 3D reconstruction algorithm can impact the quality of the reconstructed images. Specifically, the results suggested that using the multi-plane scheme with some consecutive planes of data enables us to utilize only a subset of data during the reconstruction and to get rid of the unwanted or potentially corrupted measurements. 242 The redesigned input–file system, as well as the optimization of the reconstruction algorithm has made it possible to develop a not only viable and fast, but also a user–friendly 3D reconstruction module, which is computationally feasible and balances the tradeoff between the accuracy of efficiency of the model. Taking advantage of the parallel computing and MATLAB executable files (MEX) has enabled us to reconstruct a 3D image in minutes, instead of previously days or even weeks. 2) We implemented a new approach combining the functional information of the MI with the high spatial resolution of other imaging modalities such as MRI and X–ray CT. This approach was based on developing new image reconstruction/regularization strategies (called soft and hard priors) that exploit structural information about the object being imaged through additional constraints to the microwave property contrast. We evaluated the performance of both 2D and 3D image reconstruction algorithms with the inclusion of prior structural information through a series of simulations, phantom experiments, and initial clinical data. Different factors such as: frequency, noise level, shape, size and number of regions of interest, imperfect or false priors, and contrast levels were studies in depth to assess the robustness of the new image technique. The results showed that incorporating prior high–resolution spatial information into the MI can significantly enhance the quality and accuracy of the recovered dielectric property distributions. In some cases, the weighted permittivity and conductivity errors were reduced by over 80% when the soft prior regularization was used. The hard priors had a similar effect, reducing the relative errors by up to 95%. In addition, it was shown that additional target regions do not have a major effect on the reconstructed soft prior dielectric property profiles. Moreover, the soft prior regularization 243 algorithm is not overly sensitive to the exact size and location of the target, and it can tolerate priors’ imperfection to some extend. For example, in 3D phantom experiments, changing the priors’ volume of the target region over 125% only produced 31% change in the recovered permittivity and 15% change in the counterpart conductivity values. Similarly, shifting the prior target region by over 85% only changed the recovered permittivity and conductivity values by 4% and 15%, respectively. Additionally, it was shown that the soft prior technique can effectively handle misleading situations where even an incorrect region of interest is identified. In those cases, the relative error in the false target region was always under 7%, and in fact, this error became even smaller (approaching zero) as the contrast between the actual target and the background medium was decreased. 3) The final outcome of this project was two graphical user interfaces (GUI) for the 2D and 3D FDTD microwave image reconstructions. These GUIs were designed to enable the user to control and adjust important parameters in the reconstruction procedure, and to visualize the final reconstructed images more conveniently. 244 REFERENCES [1] American Cancer Society., "Cancer Facts & Figures 2011," American Cancer Society, Inc, Atlanta, Ga.2011. [2] K. Berger, J. Bostwick, and G. E. Jones, A woman's decision: breast care, treatment & reconstruction, 4th ed. St. Louis, Mo.: Quality Medical Pub, 2011. [3] R. Smith-Bindman, P. Chu, D. L. Miglioretti, C. Quale, R. D. Rosenberg, G. Cutter, B. Geller, P. Bacchetti, E. A. Sickles, and K. Kerlikowske, "Physician predictors of mammographic accuracy," Journal of the National Cancer Institute, vol. 97, pp. 358-367, Mar 2 2005. [4] R. D. Rosenberg, W. C. Hunt, M. R. Williamson, F. D. Gilliland, P. W. Wiest, C. A. Kelsey, C. R. Key, and M. N. Linver, "Effects of age, breast density, ethnicity, and estrogen replacement therapy on screening mammographic sensitivity and cancer stage at diagnosis: Review of 183,134 screening mammograms in Albuquerque, New Mexico," Radiology, vol. 209, pp. 511-518, Nov 1998. [5] J. E. Joy, E. E. Penhoet, D. B. Petitti, National Cancer Policy Board (U.S.). Committee on New Approaches to Early Detection and Diagnosis of Breast Cancer., National Research Council (U.S.). Policy and Global Affairs., and National Research Council (U.S.). Board on Science Technology and Economic Policy., Saving women's lives : strategies for improving breast cancer detection and diagnosis. Washington, D.C.: National Academies Press, 2005. [6] M. F. Ernst and J. A. Roukema, "Diagnosis of non-palpable breast cancer: a review," Breast, vol. 11, pp. 13-22, Feb 2002. 245 [7] T. Hata, H. Takahashi, K. Watanabe, M. Takashi, K. Taguchi, T. Itoh, and S. Todo, "Magnetic resonance imaging for preoperative evaluation of breast cancer: A comparative study with mammography and ultrasonography," Journal of the American College of Surgeons, vol. 198, pp. 190-197, Feb 2004. [8] K. D. Paulsen, P. M. Meaney, and L. C. Gilman, Alternative breast imaging : four modelbased approaches. New York: Springer Science+Business Media, 2005. [9] P. M. Meaney, M. W. Fanning, T. Raynolds, C. J. Fox, Q. Q. Fang, C. A. Kogel, S. P. Poplack, and K. D. Paulsen, "Initial clinical experience with microwave breast imaging in women with normal mammography," Academic Radiology, vol. 14, pp. 207-218, Feb 2007. [10] N. K. Soni, A. Hartov, C. Kogel, S. P. Poplack, and K. D. Paulsen, "Multi-frequency electrical impedance tomography of the breast: new clinical results," Physiol Meas, vol. 25, pp. 301-14, Feb 2004. [11] S. Srinivasan, B. W. Pogue, B. Brooksby, S. D. Jiang, H. Dehghani, C. Kogel, W. A. Wells, S. P. Poplack, and K. D. Paulsen, "Near-infrared characterization of breast tumors in vivo using spectrally-constrained reconstruction," Technology in Cancer Research & Treatment, vol. 4, pp. 513-526, Oct 2005. [12] E. E. Van Houten, M. M. Doyley, F. E. Kennedy, J. B. Weaver, and K. D. Paulsen, "Initial in vivo experience with steady-state subzone-based MR elastography of the human breast," Journal of Magnetic Resonance Imaging, vol. 17, pp. 72-85, Jan 2003. [13] S. S. Chaudhary, R. K. Mishra, A. Swarup, and J. M. Thomas, "Dielectric properties of normal & malignant human breast tissues at radiowave & microwave frequencies," Indian Journal Biochemical Biophys, vol. 21, pp. 76-9, Feb 1984. 246 [14] W. T. Joines, Y. Zhang, C. Li, and R. L. Jirtle, "The measured electrical properties of normal and malignant human tissues from 50 to 900 MHz," Medical Physics, vol. 21, pp. 547-50, Apr 1994. [15] A. J. Surowiec, S. S. Stuchly, J. B. Barr, and A. Swarup, "Dielectric properties of breast carcinoma and the surrounding tissues," IEEE Transactions on Biomedical Engineering, vol. 35, pp. 257-63, Apr 1988. [16] M. Lazebnik, D. Popovic, L. McCartney, C. B. Watkins, M. J. Lindstrom, J. Harter, S. Sewall, T. Ogilvie, A. Magliocco, T. M. Breslin, W. Temple, D. Mew, J. H. Booske, M. Okoniewski, and S. C. Hagness, "A large-scale study of the ultrawideband microwave dielectric properties of normal, benign and malignant breast tissues obtained from cancer surgeries," Phys Med Biol, vol. 52, pp. 6093-115, Oct 21 2007. [17] S. P. Poplack, K. D. Paulsen, A. Hartov, P. M. Meaney, B. W. Pogue, T. D. Tosteson, M. R. Grove, S. K. Soho, and W. A. Wells, "Electromagnetic breast imaging: average tissue property values in women with negative clinical findings," Radiology, vol. 231, pp. 57180, May 2004. [18] S. P. Poplack, T. D. Tosteson, W. A. Wells, B. W. Pogue, P. M. Meaney, A. Hartov, C. A. Kogel, S. K. Soho, J. J. Gibson, and K. D. Paulsen, "Electromagnetic breast imaging: results of a pilot study in women with abnormal mammograms," Radiology, vol. 243, pp. 350-9, May 2007. [19] E. C. Fear, S. C. Hagness, P. M. Meaney, M. Okoniewski, and M. A. Stuchly. (2002) Enhancing breast tumor detection with near-field imaging. IEEE Microwave magazine. 4856. 247 [20] F. Bardati and S. Ludicello, "Modeling the visibility of breast malignancy by a microwave radiometer," IEEE Transactions on Biomedical Engineering, vol. 55, pp. 214-221, Jan 2008. [21] A. H. Barrett, P. C. Myers, and N. L. Sadowsky, "Detection of Breast-Cancer by Microwave Radiometry," Radio Science, vol. 12, pp. 167-171, 1977. [22] K. R. Foster and E. A. Cheever, "Microwave Radiometry in Biomedicine - a Reappraisal," Bioelectromagnetics, vol. 13, pp. 567-579, 1992. [23] G. Giaux, J. Delannoy, D. Delvalee, Y. Leroy, A. Mamouni, J. C. Vandevelde, and B. Bocquet, "Microwave Radiometry Imaging - Characterization of Breast-Tumors," Journal of Photographic Science, vol. 39, pp. 164-165, 1991. [24] H. Hinrikus, J. Riipulk, and H. Podra, "Design of microwave radiometer for early detection of cancer," Physica Medica, vol. 13, pp. 324-325, Dec 1997. [25] S. Iudicello and F. Bardati, "Functional Imaging of Compressed Breast by Microwave Radiometry," Applied Computational Electromagnetics Society Journal, vol. 24, pp. 64-71, Feb 2009. [26] B. Guo, J. Li, H. Zmuda, and M. Sheplak, "Multifrequency microwave-induced thermal acoustic imaging for breast cancer detection," IEEE Transactions on Biomedical Engineering, vol. 54, pp. 2000-10, Nov 2007. [27] G. Ku and L. V. Wang, "Scanning microwave-induced thermoacoustic tomography: signal, resolution, and contrast," Medical Physics, vol. 28, pp. 4-10, Jan 2001. [28] E. G. Moros, J. Penagaricano, P. Novak, W. L. Straube, and R. J. Myerson, "Present and future technology for simultaneous superficial thermoradiotherapy of breast cancer," International Journal of Hyperthermia, vol. 26, pp. 699-709, 2010. 248 [29] P. M. Meaney, M. W. Fanning, D. Li, S. P. Poplack, and K. D. Paulsen, "A clinical prototype for active microwave imaging of the breast," IEEE Transactions on Microwave Theory and Techniques, vol. 48, pp. 1841-1853, Nov 2000. [30] E. C. Fear and M. A. Stuchly, "Microwave system for breast tumor detection," IEEE Microwave and Guided Wave Letters, vol. 9, pp. 470-472, Nov 1999. [31] P. M. Meaney, K. D. Paulsen, A. Hartov, and R. K. Crane, "An Active Microwave Imaging-System for Reconstruction of 2-D Electrical Property Distributions," Ieee Transactions on Biomedical Engineering, vol. 42, pp. 1017-1026, Oct 1995. [32] S. Caorsi, G. L. Gragnani, M. Pastorino, and M. Rebagliati, "An Approach to Focused 2d Tm near-Field Active Microwave Imaging of Dielectric Objects," Microwave and Optical Technology Letters, vol. 8, pp. 302-306, Apr 20 1995. [33] V. Zhurbenko, T. Rubaek, V. Krozer, and P. Meincke, "Design and realisation of a microwave three-dimensional imaging system with application to breast-cancer detection," Iet Microwaves Antennas & Propagation, vol. 4, pp. 2200-2211, Dec 2010. [34] P. M. Meaney, K. D. Paulsen, A. Hartov, and R. K. Crane, "An active microwave imaging system for reconstruction of 2-D electrical property distributions," IEEE Transactions on Biomedical Engineering, vol. 42, pp. 1017-26, Oct 1995. [35] E. C. Fear, X. Li, S. C. Hagness, and M. A. Stuchly, "Confocal microwave imaging for breast cancer detection: localization of tumors in three dimensions," IEEE Trans Biomedical Engineering, vol. 49, pp. 812-22, Aug 2002. [36] L. V. Wang, "Prospects of photoacoustic tomography," Medical physics, vol. 35, pp. 575867, 2008. 249 [37] S. C. Hagness, A. Taflove, and J. E. Bridges, "Two-dimensional FDTD analysis of a pulsed microwave confocal system for breast cancer detection: Fixed-focus and antennaarray sensors," IEEE Transactions on Biomedical Engineering, vol. 45, pp. 1470-1479, Dec 1998. [38] S. C. Hagness, A. Taflove, and J. E. Bridges, "Three-dimensional FDTD analysis of a pulsed microwave confocal system for breast cancer detection: Design of an antenna-array element," IEEE Transactions on Antennas and Propagation, vol. 47, pp. 783-791, May 1999. [39] E. C. Fear, X. Li, S. C. Hagness, and M. A. Stuchly, "Confocal microwave imaging for breast cancer detection: localization of tumors in three dimensions," IEEE Transactions on Biomedical Engineering, vol. 49, pp. 812-22, 2002. [40] J. M. Sill and E. C. Fear, "Tissue sensing adaptive radar for breast cancer detection Experimental investigation of simple tumor models," IEEE Transactions on Microwave Theory and Techniques, vol. 53, pp. 3312-3319, Nov 2005. [41] E. J. Bond, X. Li, S. C. Hagness, and B. D. Van Veen, "Microwave imaging via space-time beamforming for early detection of breast cancer," IEEE Transactions on Antennas and Propagation, vol. 51, pp. 1690-1705, Aug 2003. [42] X. Li, S. K. Davis, S. C. Hagness, D. W. van der Weide, and B. D. Van Veen, "Microwave imaging via space-time beamforming: Experimental investigation of tumor detection in multilayer breast phantoms," IEEE Transactions on Microwave Theory and Techniques, vol. 52, pp. 1856-1865, Aug 2004. [43] J. M. Sill and E. C. Fear, "Tissue sensing adaptive radar for breast cancer detection: study of immersion liquids," Electronics Letters, vol. 41, pp. 113-115, Feb 3 2005. 250 [44] E. G. Fear and J. M. Sill, "Preliminary investigations of tissue sensing adaptive radar for breast tumor detection," in Engineering in Medicine and Biology Society, 2003. Proceedings of the 25th Annual International Conference of the IEEE, 2003, pp. 37873790 Vol.4. [45] M. Klemm, I. Craddock, J. Leendertz, A. Preece, and R. Benjamin, "Experimental and clinical results of breast cancer detection using UWB microwave radar," in Antennas and Propagation Society International Symposium, 2008. AP-S 2008. IEEE, 2008, pp. 1-4. [46] The American Heritage dictionary of the English language, 4th ed. Boston: Houghton Mifflin, 2000. [47] C. Yu, M. Q. Yuan, J. Stang, E. Bresslour, R. T. George, G. A. Ybarra, W. T. Joines, and Q. H. Liu, "Active microwave imaging II: 3-D system prototype and image reconstruction from experimental data," IEEE Transactions on Microwave Theory and Techniques, vol. 56, pp. 991-1000, Apr 2008. [48] X. Millard and Q. H. Liu, "A fast volume integral equation solver for electromagnetic scattering from large inhomogeneous objects in planarly layered media," IEEE Transactions on Antennas and Propagation, vol. 51, pp. 2393-2401, Sep 2003. [49] X. Xu, Q. H. Liu, and Z. W. Zhang, "The stabilized biconjugate gradient fast Fourier transform method for electromagnetic scattering," in IEEE Antennas and Propagation Society International Symposium, 2002, pp. 614-617 vol.2. [50] X. M. Xu and Q. H. Liu, "The BCGS-FFT Method for Electromagnetic Scattering From Inhomogeneous Objects in a Planarly Layered Medium," IEEE Antennas and Wireless Propagation Letters, vol. 1, pp. 77-80, 2002. 251 [51] C. Yu, M. Q. Yuan, Y. J. Zhang, J. Stang, R. T. George, G. A. Ybarra, W. T. Joines, and Q. H. Liu, "Microwave Imaging in Layered Media: 3-D Image Reconstruction From Experimental Data," IEEE Transactions on Antennas and Propagation, vol. 58, pp. 440448, Feb 2010. [52] S. Y. Semenov, A. E. Bulyshev, A. E. Souvorov, A. G. Nazarov, Y. E. Sizov, R. H. Svenson, V. G. Posukh, A. Pavlovsky, P. N. Repin, and G. P. Tatsis, "Three-dimensional microwave tomography: Experimental imaging of phantoms and biological objects," IEEE Transactions on Microwave Theory and Techniques, vol. 48, pp. 1071-1074, Jun 2000. [53] S. Y. Semenov, R. H. Svenson, A. E. Boulyshev, A. E. Souvorov, V. Y. Borisov, Y. Sizov, A. N. Starostin, K. R. Dezern, G. P. Tatsis, and V. Y. Baranov, "Microwave tomography: Two-dimensional system for biological imaging," IEEETransactions on Biomedical Engineering, vol. 43, pp. 869-877, Sep 1996. [54] S. Y. Semenov, R. H. Svenson, A. E. Bulyshev, A. E. Souvorov, A. G. Nazarov, Y. E. Sizov, V. G. Posukh, A. Pavlovsky, P. N. Repin, A. N. Starostin, B. A. Voinov, M. Taran, G. P. Tatsis, and V. Y. Baranov, "Three-dimensional microwave tomography: Initial experimental imaging of animals," IEEETransactions on Biomedical Engineering, vol. 49, pp. 55-63, Jan 2002. [55] A. E. Bulyshev, A. E. Souvorov, S. Y. Semenov, R. H. Svenson, A. G. Nazarov, Y. E. Sizov, and G. P. Tatsis, "Three-dimensional microwave tomography. Theory and computer experiments in scalar approximation," Inverse Problems, vol. 16, pp. 863-875, Jun 2000. [56] S. Y. Semenov, R. H. Svenson, A. E. Bulyshev, A. E. Souvorov, A. G. Nazarov, Y. E. Sizov, A. V. Pavlovsky, V. Y. Borisov, B. A. Voinov, G. I. Simonova, A. N. Starostin, V. G. Posukh, G. P. Tatsis, and V. Y. Baranov, "Three-dimensional microwave tomography: 252 Experimental prototype of the system and vector born reconstruction method," IEEETransactions on Biomedical Engineering, vol. 46, pp. 937-946, Aug 1999. [57] A. E. Souvorov, A. E. Bulyshev, S. Y. Semenov, R. H. Svenson, A. G. Nazarov, Y. E. Sizov, and G. P. Tatsis, "Microwave tomography: A two-dimensional Newton iterative scheme," IEEE Transactions on Microwave Theory and Techniques, vol. 46, pp. 16541659, Nov 1998. [58] S. Y. Semenov, R. H. Svenson, V. G. Posukh, A. G. Nazarov, Y. E. Sizov, A. E. Bulyshev, A. E. Souvorov, W. Chen, J. Kasell, and G. P. Tatsis, "Dielectrical spectroscopy of canine myocardium during acute ischemia and hypoxia at frequency spectrum from 100 kHz to 6 GHz," IEEE Trans Med Imaging, vol. 21, pp. 703-7, Jun 2002. [59] S. Semenov, J. Kellam, P. Althausen, T. Williams, A. Abubakar, A. Bulyshev, and Y. Sizov, "Microwave tomography for functional imaging of extremity soft tissues: feasibility assessment," Physics in medicine and biology, vol. 52, pp. 5705-19, 2007. [60] Q. H. Liu, Z. Q. Zhang, T. H. T. Wang, J. A. Bryan, G. A. Ybarra, L. W. Nolte, and W. T. Joines, "Active microwave imaging I-2-D forward and inverse scattering-methods," IEEE Transactions on Microwave Theory and Techniques, vol. 50, pp. 123-133, Jan 2002. [61] S. Padhi, J. Howard, A. Fhager, and S. Bengtsson, "A PC-controlled microwave tomographic scanner for breast imaging," Review of Scientific Instruments, vol. 82, pp. -, Jan 2011. [62] P. M. Meaney, S. A. Pendergrass, M. W. Fanning, D. Li, and K. D. Paulsen, "Importance of using a reduced contrast coupling medium in 2D microwave breast imaging," Journal of Electromagnetic Waves and Applications, vol. 17, pp. 333-355, 2003. 253 [63] P. M. Meaney, K. D. Paulsen, and T. P. Ryan, "2-Dimensional Hybrid Element ImageReconstruction for Tm Illumination," Ieee Transactions on Antennas and Propagation, vol. 43, pp. 239-247, Mar 1995. [64] V. Semenov, M. Lisak, and D. Anderson, "Electric-Field Enhancement and Power Absorption in Microwave Tr-Switches," IEEE Transactions on Microwave Theory and Techniques, vol. 43, pp. 286-292, Feb 1995. [65] P. M. Meaney, Q. Q. Fang, T. Rubaek, E. Demidenko, and K. D. Paulsen, "Log transformation benefits parameter estimation in microwave tomographic imaging," Medical Physics, vol. 34, pp. 2014-2023, Jun 2007. [66] Q. Fang, "Computational methods for microwave medical imaging " PhD, Thayer School of Engineering, Dartmouth College, Hanover, 2004. [67] Q. Fang, P. M. Meaney, and K. D. Paulsen, "Viable Three-Dimensional Medical Microwave Tomography: Theory and Numerical Experiments," IEEE Trans Antennas and Propagation, vol. 58, pp. 449-458, Feb 1 2010. [68] B. Kaltenbacher, "Some Newton-type methods for the regularization of nonlinear ill-posed problems," Inverse Problems, vol. 13, pp. 729-753, Jun 1997. [69] V. Y. Arsenin and A. N. Tikhonov, Solutions of ill-posed problems. Washington Winston & Sons, 1977. [70] P. C. Hansen, Rank-deficient and discrete ill-posed problems : numerical aspects of linear inversion. Philadelphia: SIAM, 1997. [71] A. Tarantola, Inverse problem theory: methods for data fitting and model parameter estimation. Amsterdam ; New York: Elsevier, 1987. 254 [72] A. H. Golnabi, P. M. Meaney, S. D. Geimer, and K. D. Paulsen, "Comparison of no-prior and soft-prior regularization in biomedical microwave imaging," Journal Medical Physics, vol. 36, pp. 159-70, Jul 2011. [73] D. K. Cheng, Field and wave electromagnetics, 2nd ed. Reading, Mass.: Addison-Wesley Pub. Co., 1989. [74] J. A. Kong, Electromagnetic wave theory. Cambridge, Mass.: EMW Publishing, 2005. [75] S. N. Ghosh, Electromagnetic theory and wave propagation, 2nd ed. Boca Raton, Fla. New Delhi: CRC Press Narosa Pub. House, 2002. [76] MathWorks. (2010, MEX-files Guide. Available: http://www.mathworks.com/support/technotes/1600/1605.html [77] A. H. Golnabi, P. M. Meaney, N. R. Epstein, and K. D. Paulsen, "Microwave imaging for breast cancer detection: advances in three--dimensional image reconstruction," presented at the Conf Proc IEEE Eng Med Biol Soc, 2011. [78] P. C. Hansen, Discrete inverse problems : insight and algorithms. Philadelphia: Society for Industrial and Applied Mathematics, 2010. [79] A. Taflove and S. C. Hagness, Computational electrodynamics : the finite-difference timedomain method, 3rd ed. Boston: Artech House, 2005. [80] T. Rubaek, P. M. Meaney, P. Meincke, and K. D. Paulsen, "Nonlinear microwave imaging for breast-cancer screening using Gauss-Newton's method and the CGLS inversion algorithm," IEEE Transactions on Antennas and Propagation, vol. 55, pp. 2320-2331, Aug 2007. 255 [81] A. H. Golnabi, P. M. Meaney, S. Geimer, Z. Tian, and K. D. Paulsen, "Microwave tomography for bone imaging," in Biomedical Imaging: From Nano to Macro, 2011 IEEE International Symposium on, 2011, pp. 956-959. [82] P. M. Meaney, A. Golnabi, M. W. Fanning, S. D. Geimer, and K. D. Paulsen, "Dielectric volume measurements for biomedical applications," in Antenna Technology and Applied Electromagnetics and the Canadian Radio Science Meeting, 2009. ANTEM/URSI 2009. 13th International Symposium on, 2009, pp. 1-4. [83] P. M. Meaney, M. W. Fanning, R. M. di Florio-Alexander, P. A. Kaufman, S. D. Geimer, T. Zhou, and K. D. Paulsen, "Microwave tomography in the context of complex breast cancer imaging," Conference proceedings : Annual International Conference of the IEEE Engineering in Medicine and Biology Society, vol. 2010, pp. 3398-401, 2010. [84] A. H. Golnabi, P. M. Meaney, S. Geimer, and K. D. Paulsen, "Microwave imaging for breast cancer detection and therapy monitoring," presented at the Biomedical Wireless Technologies, Networks, and Sensing Systems (BioWireleSS), 2011 IEEE Topical Conference on, 2011. [85] B. W. Pogue, T. O. McBride, J. Prewitt, U. L. Osterberg, and K. D. Paulsen, "Spatially variant regularization improves diffuse optical tomography," Applied Optics, vol. 38, pp. 2950-61, May 1 1999. [86] A. Li, E. L. Miller, M. E. Kilmer, T. J. Brukilacchio, T. Chaves, J. Stott, Q. Zhang, T. Wu, M. Chorlton, R. H. Moore, D. B. Kopans, and D. A. Boas, "Tomographic optical breast imaging guided by three-dimensional mammography," Applied Optics, vol. 42, pp. 518190, 2003. 256 [87] H. Dehghani, B. W. Pogue, S. D. Jiang, B. Brooksby, and K. D. Paulsen, "Threedimensional optical tomography: resolution in small-object imaging," Applied Optics, vol. 42, pp. 3117-3128, Jun 1 2003. [88] B. W. Pogue and K. D. Paulsen, "High-resolution near-infrared tomographic imaging simulations of the rat cranium by use of apriori magnetic resonance imaging structural information," Optics Letters, vol. 23, pp. 1716-1718, Nov 1 1998. [89] B. A. Brooksby, H. Dehghani, B. W. Pogue, and K. D. Paulsen, "Near-infrared (NIR) tomography breast image reconstruction with a priori structural information from MRI: Algorithm development for reconstructing heterogeneities," IEEE Journal of Selected Topics in Quantum Electronics, vol. 9, pp. 199-209, Mar-Apr 2003. [90] L. Solymar and D. Walsh, Electrical properties of materials, 8th ed. Oxford ; New York: Oxford University Press, 2010. [91] C. T. Kelley, Iterative methods for linear and nonlinear equations. Philadelphia: Society for Industrial and Applied Mathematics, 1995. [92] R. Pierri and A. Tamburrino, "On the local minima problem in conductivity imaging via a quadratic approach," Inverse Problems, vol. 13, pp. 1547-1568, Dec 1997. [93] J. A. Scales, M. L. Smith, and T. L. Fischer, "Global Optimization Methods for Multimodal Inverse Problems," Journal of Computational Physics, vol. 103, pp. 258-268, Dec 1992. [94] P. M. Meaney, K. D. Paulsen, B. W. Pogue, and M. I. Miga, "Microwave image reconstruction utilizing log-magnitude and unwrapped phase to improve high-contrast object recovery," IEEE Transactions on Medical Imaging, vol. 20, pp. 104-116, Feb 2001. 257 [95] T. M. Grzegorczyk, P. M. Meaney, S. I. Jeon, S. D. Geimer, and K. D. Paulsen, "Importance of phase unwrapping for the reconstruction of microwave tomographic images," Biomedical Optics Express, vol. 2, pp. 315-330, 2011. [96] S. Caorsi and M. Pastorino, "Two-dimensional microwave imaging approach based on a genetic algorithm," IEEE Transactions on Antennas and Propagation, vol. 48, pp. 370373, Mar 2000. [97] M. Donelli and A. Massa, "Computational approach based on a particle swarm optimizer for microwave imaging of two-dimensional dielectric scatterers," IEEE Transactions on Microwave Theory and Techniques, vol. 53, pp. 1761-1776, May 2005. [98] P. Rocca, M. Benedetti, M. Donelli, and A. Massa, "Evolutionary techniques for inverse scattering - Current trends and envisaged developments," in Antennas and Propagation Society International Symposium. APSURSI '09. IEEE, 2009, pp. 1-4. [99] M. Benedetti, M. Donelli, G. Franceschini, M. Pastorino, and A. Massa, "Effective exploitation of the a priori information through a microwave imaging procedure based on the SMW for NDE/NDT applications," IEEE Transactions on Geoscience and Remote Sensing, vol. 43, pp. 2584-2592, Nov 2005. [100] S. Caorsi, G. L. Gragnani, M. Pastorino, and M. Sartore, "Electromagnetic Imaging of Infinite Dielectric Cylinders Using a Modified Born Approximation and Including a-Priori Information on the Unknown Cross-Sections," IEEE Proceedings-Microwaves Antennas and Propagation, vol. 141, pp. 445-450, Dec 1994. [101] W. C. Chew and J. H. Lin, "A Frequency-Hopping Approach for Microwave Imaging of Large Inhomogeneous Bodies," IEEE Microwave and Guided Wave Letters, vol. 5, pp. 439-441, Dec 1995. 258 [102] O. Dorn, E. L. Miller, and C. M. Rappaport, "A shape reconstruction method for electromagnetic tomography using adjoint fields and level sets," Inverse Problems, vol. 16, pp. 1119-1156, Oct 2000. [103] O. Feron, B. Duchene, and A. Mohammad-Djafari, "Microwave imaging of inhomogeneous objects made of a finite number of dielectric and conductive materials from experimental data," Inverse Problems, vol. 21, pp. S95-S115, Dec 2005. [104] R. Ferraye, J. Y. Dauvignac, and C. Pichot, "An inverse scattering method based on contour deformations by means of a level set method using frequency hopping technique," IEEE Transactions on Antennas and Propagation, vol. 51, pp. 1100-1113, May 2003. [105] A. Fhager and M. Persson, "Using a priori data to improve the reconstruction of small objects in microwave tomography," IEEE Transactions on Microwave Theory and Techniques, vol. 55, pp. 2454-2462, Nov 2007. [106] N. Joachimowicz, C. Pichot, and J. P. Hugonin, "Inverse Scattering - an Iterative Numerical-Method for Electromagnetic Imaging," IEEE Transactions on Antennas and Propagation, vol. 39, pp. 1742-1752, Dec 1991. [107] A. Litman, "Reconstruction by level sets of n-ary scattering obstacles," Inverse Problems, vol. 21, pp. S131-S152, Dec 2005. [108] S. J. Osher and F. Santosa, "Level set methods for optimization problems involving geometry and constraints I. Frequencies of a two-density inhomogeneous drum," Journal of Computational Physics, vol. 171, pp. 272-288, Jul 20 2001. [109] L. Crocco and T. Isernia, "Inverse scattering with real data: detecting and imaging homogeneous dielectric objects," Inverse Problems, vol. 17, pp. 1573-1583, Dec 2001. 259 [110] M. El-Shenawee, O. Dorn, and M. Moscoso, "An Adjoint-Field Technique for Shape Reconstruction of 3-D Penetrable Object Immersed in Lossy Medium," IEEE Transactions on Antennas and Propagation, vol. 57, pp. 520-534, Feb 2009. [111] B. Brooksby, S. D. Jiang, H. Dehghani, B. W. Pogue, K. D. Paulsen, J. Weaver, C. Kogel, and S. P. Poplack, "Combining near-infrared tomography resonance imaging to study in vivo and magnetic breast tissue: implementation of a Laplacian-type regularization to incorporate magnetic resonance structure," Journal of Biomedical Optics, vol. 10, pp. -, Sep-Oct 2005. [112] B. W. Pogue, H. Q. Zhu, C. Nwaigwe, T. O. McBride, U. L. Osterberg, K. D. Paulsen, and J. F. Dunn, "Hemoglobin imaging with hybrid magnetic resonance and near-infrared diffuse tomography," Oxygen Transport to Tissue Xxiv, vol. 530, pp. 215-224, 2003. [113] D. R. Lynch, Numerical partial differential equations for environmental scientists and engineers : a first practical course. New York: Springer, 2005. [114] P. K. Yalavarthy, B. W. Pogue, H. Dehghani, and K. D. Paulsen, "Weight-matrix structured regularization provides optimal generalized least-squares estimate in diffuse optical tomography," Medical Physics, vol. 34, pp. 2085-2098, Jun 2007. [115] H. Ojeda-Fournier and C. E. Comstock, "MRI for breast cancer: Current indications," The Indian Journal of Radiology and Imaging, vol. 19, pp. 161-9, 2009. [116] E. Morris, L. Liberman, and SpringerLink (Online service). (2005, Breast MRI diagnosis and intervention [SpringerLink]. Available: http://dx.doi.org/10.1007/0-387-27595-9 MIT Access Only [117] H. J. Shin, H. H. Kim, J. H. Ahn, S. B. Kim, K. H. Jung, G. Gong, B. H. Son, and S. H. Ahn, "Comparison of mammography, sonography, MRI and clinical examination in 260 patients with locally advanced or inflammatory breast cancer who underwent neoadjuvant chemotherapy," The British Journal of Radiology, vol. 84, pp. 612-20, 2011. [118] K. P. McGuire, J. Toro-Burguete, H. Dang, J. Young, A. Soran, M. Zuley, R. Bhargava, M. Bonaventura, R. Johnson, and G. Ahrendt, "MRI staging after neoadjuvant chemotherapy for breast cancer: does tumor biology affect accuracy?," Annals of Surgical Oncology, vol. 18, pp. 3149-54, 2011. [119] M. Lorenzon, C. Zuiani, V. Londero, A. Linda, A. Furlan, and M. Bazzocchi, "Assessment of breast cancer response to neoadjuvant chemotherapy: is volumetric MRI a reliable tool?," European Journal of Radiology, vol. 71, pp. 82-8, 2009. [120] C. K. Kuhl, "Current status of breast MR imaging - Part 2. Clinical applications," Radiology, vol. 244, pp. 672-691, Sep 2007. [121] C. K. Kuhl, "The current status of breast MR imaging - Part I. Choice of technique, image interpretation, diagnostic accuracy, and transfer to clinical practice," Radiology, vol. 244, pp. 356-378, Aug 2007. [122] D. Galanaud, S. Haik, M. G. Linguraru, J. P. Ranjeva, B. Faucheux, E. Kaphan, N. Ayache, J. Chiras, P. Cozzone, D. Dormont, and J. P. Brandel, "Combined diffusion imaging and MR spectroscopy in the diagnosis of human prion diseases," American Journal of Neuroradiology, vol. 31, pp. 1311-8, Aug 2010. [123] A. Fatemi-Ardekani, N. Samavati, J. Tang, and M. V. Kamath, "Advances in multimodality imaging through a hybrid PET/MRI system," Critical Reviews in Biomedical Engineering, vol. 37, pp. 495-515, 2009. 261 [124] B. J. Pichler, A. Kolb, T. Nagele, and H. P. Schlemmer, "PET/MRI: paving the way for the next generation of clinical multimodality imaging applications," Journal of Nuclear Medicine, vol. 51, pp. 333-6, Mar 2010. [125] Y. Yuan, M. L. Giger, H. Li, N. Bhooshan, and C. A. Sennett, "Multimodality computeraided breast cancer diagnosis with FFDM and DCE-MRI," Academic Radiology, vol. 17, pp. 1158-67, Sep 2010. [126] C. M. Carpenter, B. W. Pogue, S. Jiang, H. Dehghani, X. Wang, K. D. Paulsen, W. A. Wells, J. Forero, C. Kogel, J. B. Weaver, S. P. Poplack, and P. A. Kaufman, "Image-guided optical spectroscopy provides molecular-specific information in vivo: MRI-guided spectroscopy of breast cancer hemoglobin, water, and scatterer size," Optics Letters, vol. 32, pp. 933-5, Apr 15 2007. [127] J. M. Golio and J. Golio, The RF and microwave handbook, 2nd ed. Boca Raton: CRC Press/Taylor & Francis Group, 2008. [128] M. Miyakawa and J. C. Bolomey, Non-invasive thermometry of the human body. Boca Raton: CRC Press, 1996. [129] D. Li, P. M. Meaney, and K. D. Paulsen, "Conformal microwave imaging for breast cancer detection," IEEE Transactions on Microwave Theory and Techniques, vol. 51, pp. 11791186, Apr 2003. [130] P. M. Meaney, F. Shubitidze, M. W. Fanning, M. Kmiec, N. R. Epstein, and K. D. Paulsen, "Surface Wave Multipath Signals in Near-Field Microwave Imaging," International Journal of Biomedical Imaging, vol. 2012, p. 11, 2012. 262 [131] P. M. Meaney, N. K. Yagnamurthy, and K. D. Paulsen, "Pre-scaled two-parameter GaussNewton image reconstruction to reduce property recovery imbalance," Physics in Medicine and Biology, vol. 47, pp. 1101-19, 2002. [132] J. Sierpowska, "Electrical and dielectric characterization of trabecular bone quality," Ph.D., University of Kuopio, Kuopio, Finland, 2007. [133] S. R. Smith and K. R. Foster, "Dielectric-Properties of Low-Water-Content Tissues," Physics in Medicine and Biology, vol. 30, pp. 965-973, 1985. [134] A. Peyman, C. Gabriel, E. H. Grant, G. Vermeeren, and L. Martens, "Variation of the dielectric properties of tissues with age: the effect on the values of SAR in children when exposed to walkie-talkie devices," Physics in Medicine and Biology, vol. 54, pp. 227-241, Jan 21 2009. [135] P. Sprent and N. C. Smeeton, Applied nonparametric statistical methods, 4th ed. Boca Raton: Chapman & Hall/CRC, 2007. [136] B. M. Ingle, S. M. Hay, H. M. Bottjer, and R. Eastell, "Changes in bone mass and bone turnover following ankle fracture," Osteoporos Int, vol. 10, pp. 408-15, 1999. [137] P. P. Sarangi, A. J. Ward, E. J. Smith, G. E. Staddon, and R. M. Atkins, "Algodystrophy and osteoporosis after tibial fractures," Journal of Bone and Joint Surgery-British Volume, vol. 75, pp. 450-2, May 1993. [138] T. Ahl, H. E. Sjoberg, and N. Dalen, "Bone mineral content in the calcaneus after ankle fracture," Acta Orthopaedica Scandinavica, vol. 59, pp. 173-5, Apr 1988. [139] M. Lazebnik, L. McCartney, D. Popovic, C. B. Watkins, M. J. Lindstrom, J. Harter, S. Sewall, A. Magliocco, J. H. Booske, M. Okoniewski, and S. C. Hagness, "A large-scale study of the ultrawideband microwave dielectric properties of normal breast tissue obtained 263 from reduction surgeries," Physics in Medicine and Biology, vol. 52, pp. 2637-2656, May 21 2007. 264 265

1/--страниц