Abstract
Microwave tomography (MT) measurement setups for different configurations based on breast compression are compared to classical circular measurement setups. Configurations based on compression allow measuring the evanescent component of the scattered field and lead to a compact measurement setup that allows direct image comparison with a standard mammography system. The different configurations are compared based on the singular value decomposition (SVD) of the radiation operator for a 2D TM case. This analysis allows determining under which conditions the image quality obtained from the reconstructions can be enhanced. These findings are confirmed by a series of reconstructions of breast phantoms based on synthetic data obtained at a single frequency of operation.
1. Introduction
Active microwave imaging for breast cancer detection has been an active research field in the past years [1]. Two different approaches have been considered, namely, radar-type approaches [2] and microwave tomography (MT) [3]. MT is used to estimate the complex permittivity distribution inside the object under test (OUT), by solving an inverse scattering problem when the object is illuminated under several different conditions and the scattered field is measured at different locations. In MT, detailed information of the OUT is obtained at the expense of solving a nonlinear and ill-posed inverse scattering problem [4].
A measurement setup in MT typically includes a series of receiving (Rx) and transmitting (Tx) antennas for measuring the scattered fields needed for the reconstruction algorithms. The antennas can be arranged in planar, cylindrical, or hemispherical surfaces. In the available clinical prototypes, the patient lies prone on a table with the breast pendant in a liquid [3, 5]. This configuration has been proposed for maximizing the comfort of the patient. In such a case, a cylindrical or hemispherical configuration is preferred since the antennas can be conformal with the shape of the breast and can accurately sample the scattered fields. However, compression of the breast may present several advantages over this classical configuration [6, 7]. The Tx and Rx antennas can be placed very close to the OUT, allowing the measurement of the evanescent scattered fields. Also, this procedure can lead to simplifications of the reconstruction algorithm since the thickness of the compressed breast is precisely known. Finally, because both MT and X-ray mammography can use the same compressed configuration a combined approach would take advantage of the benefits of both techniques, namely, the high resolution of X-rays and the quantitative information provided by MT. Such a combined approach has already been used with other techniques, that is, mammography with optical imaging or positron emission tomography with computed tomography [8]. In the case of a configuration using breast compression, a planar configuration for the Tx/Rx array will be better adapted.
Evanescent scattered fields have been exploited for a long time in the area of microwave near-field microscopy [9, 10]. This approach is commonly used in nondestructive evaluation of materials allowing spatial resolutions in the order of the nanometer. The exploitation of evanescent scattered fields in MT was proposed in [11, 12]. The authors showed that an enhancement in the reconstructed images obtained in MT was obtained by using close-proximity setups; however, this enhancement was limited to the region near the edges of the imaging domain.
In this paper, we will analyze how the quality of the reconstructed images can be improved through breast compression. First of all, we will discuss different possibilities for the implementation of a planar MT measurement setup using breast compression. Next, we will shortly introduce the analysis techniques used for the comparison of the different configurations. This will include a review of the mechanical properties of the breast and a simple analysis of the deformation of breast tissues under compression. Finally the results of the analysis will show under which conditions we can improve the image quality by exploiting the measurement of evanescent scattered fields.
2. Different Configurations for the MT Setup Based on Breast Compression
Planar measurement setups for microwave imaging have already been studied in the literature, including the planar microwave camera [13, 14] or the raster scan proposed in [15] that used breast compression. However, no analysis of the effect of breast compression on the quality of the obtained images was discussed.
The configurations considered in this paper are shown in Figure 1. Figure 1(a) shows a classical circular MT setup where both the Tx and Rx probes are placed in a circle surrounding the OUT. Usually, the Tx probes, Rx probes and OUT are all immersed in a matching medium that is properly chosen to maximize field coupling between the probes and the region of interest in the OUT. Different possibilities are available for the implementation of a planar MT setup. The first and more straightforward idea could be to implement a planar array of Tx and Rx antennas as it can be seen in Figure 1(b). The array of antennas could be immersed in the matching media, or it could be placed in the air region if a sealed container is used. This configuration has the drawback that the antennas, being in close proximity of the scatterers, will disturb the field to be measured. A possible solution to this problem is to include the antennas in the formulation of the problem taking into account their presence in the discretization of the MT system of [16, 17]. Another possibility is to use an array of small near-field probes to minimize perturbations in the measured field [18]. Two configurations based on this type of near-field probes are shown in Figures 1(c) and 1(d). The first one, which will be referred to as the camera configuration, is a setup that consists in a series of Tx antennas illuminating the OUT from one side and an array of near-field measurement probes on the other side. This setup arrangement creates a very compact system that could be possibly integrated with a standard mammography apparatus. A drawback of this configuration is the fact that we can only have transmission information through the object, which translates into a limited view of the OUT. The second configuration shown in Figure 1(d) will be called the waveguide configuration. The compound structure formed by the immersion medium and the enclosure formed by the dielectric plates and the surrounding air region can be seen as a multilayer dielectric waveguide. The structure is capable of guiding different TE () or TM ( and ) modes depending on the frequency of operation, the dielectric properties of the immersion media and dielectric plates, and the thickness of the different layers. Combinations of the different modes can be launched by using an array of transmitting antennas and varying the phase shift between them [19]. This configuration has the advantage that the scattered fields can be measured both on top and bottom of the structure, thus doubling the information compared to the camera configuration. This configuration also allows the exploitation of spatial and modal diversity to create the desired number of independent illuminations. In the 2D case this fact can be a limitation since the number of independent illuminations is limited to twice the number of possible modes in the waveguide. In the 3D case this is no longer a limiting fact since the additional dimension gives the freedom in the spatial diversity to create the desired number of independent illuminations.
(a)
(b)
(c)
(d)
3. Numerical Models
In this section the different numerical phantoms used in the reconstructions will be presented. Since we are comparing measurement setups with and without compression, it is important, in order to make a fair comparison between the configurations, to have a model that takes into account, at least approximately, the deformation of the various breast regions resulting from the compression. For this purpose, we performed a simple mechanical analysis of the compressed breast tissues and the results will be briefly discussed here. Then, the numerical phantoms used in the different reconstructions will be presented. These numerical breast models have been inspired from the anatomically realistic numerical breast phantoms presented in [20].
3.1. Mechanical Deformation of Breast Tissues
The mechanical properties of human tissues have been an active research subject in different medical applications. Breast tissue deformation is also studied for image comparison between mammography and other noncompressing techniques such as CT or MRI of the breast [21] or image guide surgery [22]. In all of these applications, breast tissues are generally approximatively modelled as linearly elastic materials. Such materials are characterized by two parameters, the Young’s modulus and Poisson’s ratio. In order to make a realistic comparison between the different configurations, we have analyzed the deformation experienced by a breast-like object composed of fat and glandular tissue using two different mechanical solvers, CST multiphysics and ANSYS mechanical solver. The simulated objects were assigned the mechanical properties given in [21] and reported in Table 1. These properties were obtained from a compilation of different published results both from in vivo and in vitro measurements of human tissues. However, breast tissues are quite heterogeneous and these values could vary greatly among different patients or parts of the breast.
The results of the analysis showed that for the mechanical properties shown in Table 1, the fat will suffer most of the deformation while the glandular tissue remains more or less unchanged. On the other hand, other studies suggest that both types of tissue present a similar value of the Young’s modulus [23]. These later studies have been performed in the framework of the small deformations used in elastography. This is not the case in the levels of compression used in mammography, and therefore a mechanical modelling able to handle large deformations and the expected nonlinearity in the mechanical behavior should be used in future studies.
3.2. Numerical Breast Phantoms
The numerical phantoms that will be used in the different reconstructions performed in Section 5 are shown in Figure 2. They differ in the content of fibroglandular tissue, and they will be identified in the sequel as mostly fatty, scattered fibroglandular, and very dense. All these phantoms consist of a mixture of fat and glandular tissue, whose electrical properties at the frequency of 2.45 GHz are given in Table 1 [24]. In the uncompressed phantoms, the fat portion is a circle of 29 mm in radius and the glandular portion varies from one phantom to the other to reflect the percentage of glandular tissue corresponding to each category. In the previous section, it was discussed how under compression the fat tissue suffers more deformation than the glandular tissue due to its lower stiffness. Taking this into account, it was assumed that the fat portion is deformed into a flattened ellipse keeping the volume of the breast constant, whereas the glandular portion is kept unchanged. The same volume is also maintained for the solution domain .
(a)
(b)
(c)
(d)
(e)
(f)
4. Analysis Techniques
In this section, the inversion algorithm used for the reconstructions will be shortly explained as well as an efficient way to address the inhomogeneous inverse scattering problem of the breast compression based configurations. Finally, the singular value analysis of the radiation operator will also be introduced.
4.1. Inverse Problem
In a general MT problem the OUT is placed inside a solution domain and it is sequentially illuminated using a series of different transmitting antennas. For each illumination the scattered field is measured at points. For all the configurations shown in Figure 1, a 2D TM case is considered; that is, the objects are considered infinite along the direction and it is assumed that the electric field has only a component. A system of two integral equations defines the scattering problem. Following the development in [25], we can discretize the system of equations using the method of moments (MOM) with pulse basis functions and point matching and obtain the following system of matrix equations:where we have . The domain is discretized in pixels (). Vector of size contains the measured scattered fields due to the th illumination. is the incident field vector of length in the domain when no OUT is present. and are the observation and coupling Green matrices of size and , respectively. These matrices are calculated from the discretization of the 2D Green’s function (GF) of the domain .
and are the unknowns of the problem. is a diagonal matrix (), where is a vector of length containing the discretized version of the contrast function which is defined as . is the unknown complex permittivity in the discretized region and is the complex permittivity of the background medium. Finally, is a vector which contains the induced polarization currents in the OUT (), where is the total electric field in the domain . The unknowns and are found through a modified version of the well-known CSI technique presented in [26]. The quality of the reconstructed images will be compared in terms of the mean square error (MSE) for the quantity of interest after a fixed number of iterations of the algorithm. MSE is defined as , in which is the known imposed contrast.
The configurations for the MT proposed in Figure 1 and based on breast compression need some modifications for the solution of the MT inverse problem. In these configurations we have to consider the dielectric medium of the compression plates. Following the general formulation, we should use the homogeneous GF for the air region and define the contrast function for the whole region containing the plates and the immersion medium, which is leading to a considerable computational overhead. However, since we are only interested in imaging the region included between the compression plates it would be quite advantageous to reduce the size of the discretization region. This is quite possible because the geometry and electrical properties of the plates are well known. For this purpose a modified GF and incident field need to be computed to take into account the presence of the inhomogeneous background.
The computation of inhomogeneous GF is widely treated in the literature, especially for the case of infinite multilayered media [27–29]. However, the geometries presented in Figure 1 are laterally closed to contain the matching fluid and allow a compact measurement setup; in such a case our inhomogeneous media cannot be considered as infinite multilayered media. Another possibility is to compute this GF numerically. According to [27], the inhomogeneous discrete version of the GF for a 2D TM scalar case, , can be expressed as where is the GF corresponding to the homogeneous case, and is the contribution of the inhomogeneous part of the structure. can be obtained combining equations (1a) and (1b). As an example, the discrete GF can be expressed as: where and are defined in the smaller region corresponding to the imaging domain defined in Figures 1(c) and 1(d). On the other hand, and are the free space discrete GF and inhomogeneous contrast, and they are defined in a region that contains the whole structure including the dielectric plates. For this reason on the term we only have to take into account the pixels contained in the smaller region . The inhomogeneous contrast has been defined with respect to the air region ; for this reason a constant term is used to obtain the GF relative to the new background permittivity, . The same procedure can be obtained to calculate and .
An example of this procedure is shown here and the results are compared to a simulation using CST Microwave Studio. The geometry shown in Figure 3(a) consists of a slab of glycerine (, ), with a contour of plexyglass (, ). The size of the glycerine slab is 240 by 38 mm and the thickness of the plexyglass contour is 6 mm. The scattering object is a rectangle of 12 by 10 mm, with an electric contrast of 2 referred to the glycerine permittivity, located at ; and . The object is illuminated by 4 electric line sources located at and and the scattered field is measured at at 80 positions along the direction (40 located at and 40 at ) with a spacing of 2 mm between the measurement probes. The results obtained with (1a), (1b), (2), (3a), and (3b) implemented in MATLAB and the finite volume time domain formulation of CST Microwave Studio are shown in Figure 3(b) both for the amplitude and phase of the scattered field. Here the four sets of 80 measurements are shown in sequence to ease visualization. CST Microwave Studio is a 3D solver so appropriate boundary conditions are used in the simulation to represent a 2D geometry. The mesh size is for the MOM simulation () and assigned to in CST at a frequency of 2.45 GHz. As it can be seen from the figure, the results for these two very different techniques agree very well both for the amplitude and the phase. This fact validates the approach used for the calculation of the scattering from an inhomogeneous object using an inhomogeneous Green’s function. This formulation reduces considerably the size of Green’s matrices and and the incident electric field vector . For instance, in the given example, the size of is reduced from 3150 by 3150 elements to 1026 by 1026 elements, considering an imaging region of 54 by 19 pixels. This has a considerable impact on the reconstructions that will be performed in Section 5, since the CSI type algorithms involve a lot of matrix vector products of the quantities previously mentioned.
(a)
(b)
4.2. SVD Analysis
Several tools are available in the literature to a priori analyze the capability of an MT system to reconstruct a given object from its scattered field measurements. Among these techniques we can cite the SVD analysis of the radiation operator [30, 31] or more recently the Fisher information analysis [32], singular value analysis of the Jacobian matrix [33], and the scattering matrix decomposition [34]. In this work we will be restricted to the SVD analysis of , which is a useful tool to qualitatively compare the behavior of different MT configuration of Tx and Rx. This technique has been previously used to study the behavior of different measurement configurations in MT for geophysical exploration [12, 35] or general purpose imaging systems [36].
This analysis is based on the behavior of the singular values of the radiation operator . It has been shown in [12] that for a given value of signal-to-noise ratio (SNR) or measurement accuracy, a higher number of significant singular values will increase the “information content” contained in a set of measured scattered fields. In some cases it is possible to obtain an analytic expression for the singular values; however, in this paper a numerical procedure will be used since the radiation operator is computed numerically in case of the waveguide and camera configurations. The radiation operator of size admits an SVD decomposition of the form , where and are square matrices of size and , respectively, and is a diagonal matrix whose diagonal contains a vector with the singular values of the matrix .
5. Results
5.1. SVD Analysis of the Configurations
Here the results of the SVD analysis of the configurations are presented. The three configurations considered in this section, circular, camera, and waveguide, share the same immersion medium previously defined. The dielectric plates have the same material and thickness as in the example previously presented in Section 4. The results for the different cases are shown in Figure 4. The singular values are plotted for different radii of the circumference where the Rx antennas are placed in the circular configuration and for the waveguide and camera configuration. In both cases, the breast is compressed up to a value of 38 mm which represents a level of compression similar to a standard mammography test [37]. In the circular configuration, a homogeneous GF is used, whereas in the waveguide and camera configuration an inhomogeneous GF is considered.
Some interesting remarks can be made based on Figure 4. The spectrum of the singular values for the circular configuration presents a rate of decay that becomes less abrupt as the probes are placed closer to the imaging region [11]. In the case of the waveguide and camera configurations, the rate of decay of the singular values is much smoother and follows a quasilinear tendency. This effect is due to the discontinuity between the different media and the close distance between the Rx probes and the imaging region [35]. Another important effect is that the camera configuration has much less significant singular values than the waveguide configuration, since the Rx probes are placed on both sides of the imaging region. The largest number of significant singular values from all of the cases is obtained with the waveguide configuration, and thus an enhanced quality in the reconstructed images is expected from this configuration. Another interesting fact is that the change of the rate of decay of the singular value spectrum makes the configurations much more sensitive to the SNR in the measurements. This means that when a high SNR is available in the measurement system, there is an advantage in probing the scattered field much closer to the imaging region. On the other hand, when the SNR decreases, the additional “information” that can be extracted from the scattered field will be hidden by noise or measurement errors.
5.2. Effect of the Compression Plates on the Reconstruction Capabilities
In this section, the effect of the compression plates will be analyzed. Apart from the breast tissues, there are three media present both in the camera and waveguide configuration: the immersion medium, the compression plates, and the surrounding air where the measurement takes place. Usually, there is a large difference in permittivity between the air and the immersion medium, and this medium is specifically chosen to couple the maximum energy inside the breast. In principle it could be possible to enhance the quantity of scattered field transmitted to the air region by an appropriate choice of the thickness and electrical properties of the compression plates. For this purpose the SVD decomposition will be performed for different thickness and permittivity values of the compression plates.
The SVD analysis for the waveguide configuration is shown in Figure 5. Figure 5(a) shows the spectra of singular values for different thicknesses of the compression plates made of dielectric material with and corresponding to the properties of acrylic plastic at 2.45 GHz. Figure 5(b) shows the variation of the singular values for different values of for a constant thickness of 6 mm. The position of the receivers is always kept at a distance of 5 mm of the top and bottom plates. As it can be seen from the figures, the thickness of the plates plays a major role in the reconstruction capabilities of the configuration, whereas the permittivity of the plates produces essentially no variation in the configuration of the singular values. From this result we can conclude that the plates should be kept as thin as possible in order to place the receivers as close as possible to the OUT.
(a)
(b)
5.3. Reconstruction of the Numerical Phantoms
In order to confirm the results obtained from the SVD analysis of the operator, this section is showing the reconstructions of the three proposed phantoms with and without compression and for different levels of SNR. The synthetic measured data are generated using a MOM code and Gaussian noise is added to represent experimental conditions. The CSI algorithm of [26] is used and the results of the reconstructed images are compared after 2000 iterations. This number of iterations assures that the algorithm has reached convergence. In all the cases, the starting point of the algorithm is a zero contrast, and no a priori information about the shape of the object is used. All the details about the implementation of the inversion algorithm can be found in the cited reference.
As it was previously mentioned in Section 2, the number of independent illuminations that can be generated with the waveguide configuration is limited to the number of independent modes that can propagate in the structure. The number of modes can be calculated using the transmission matrix method described in [38]. Wave propagation is considered along the direction with the electric field oriented in the direction (i.e., TE modes in the formulation of the cited reference). For a multilayer structure composed of a slab of the immersion medium described in Table 1 of thickness 38 mm, compression plates of plexyglass ( and ) of thickness 6 mm, and the outer air region, there is a total of three independent modes propagating (i.e., above cut-off) at the frequency of operation, 2.45 GHz. From this result the number of transmitters is fixed to (i.e., three modes can be launched on each side of the structure) and the number of receivers is fixed to , with half of them located at and the other half at for the waveguide structure. These receivers are located 5 mm above and below the compression plates. For the circular structure, the numbers of transmitters and receivers are the same; however, in this case they are equally spaced around the different circumference radii.
The reconstructions were performed for the three different phantoms and for different levels of SNR. The MSE of the different reconstructions as a function of the SNR is shown in Figure 6 for the three considered phantoms. These results agree with the previous conclusions obtained about the spectrum of singular values. As it can be seen from the figures, the quality of the reconstructions is greatly enhanced when the probes are placed in close proximity of the imaging region. In the mostly fatty case, the waveguide is clearly better than all the circular setups considered. In the case of scattered fibroglandular and very dense breast model, the waveguide setup has a better performance for large SNR values while the performance degrades rapidly as the noise level increases. In general, it can be concluded that an SNR value of around 40 dB is needed to see a significant enhancement of the quality of the reconstructed images. This value of SNR has been previously achieved in state-of-the-art measurement systems [14, 39] in which the field measurement is done in close proximity of the OUT. When the SNR is reduced, the evanescent component of the scattered field is masked by the noise, and the quality of the reconstructed images starts to degrade. As an example, the reconstructions of the phantoms shown in Figure 2, obtained with and without compression, are shown in Figures 7 and 8. The enhancement can be clearly seen in Figures 7(c) and 7(d), where in the latter case the 3 glandular regions can be clearly distinguished whereas in the former case it is not possible to localize the edges of the glandular regions.
(a)
(b)
(c)
(a)
(b)
(c)
(d)
(e)
(f)
(a)
(b)
(c)
(d)
(e)
(f)
6. Conclusion
In this paper we have analyzed different configurations for the MT problem based on breast compression. These configurations were designed in the context of creating a compact measurement setup that uses the same compression geometry as standard mammography apparatus. The different configurations were compared based on the SVD of the radiation operator and through the reconstruction of several breast-like objects with and without compression using synthetically generated measurements for different levels of SNR.(i)It has been shown that quantitative and qualitative enhancements of the reconstructed images can be obtained using breast compression compared to a classical circular MT measurement system for the same measurement conditions (i.e., equal number of Tx and Rx antennas and same SNR level). This enhancement is dependent on the SNR level. In our case, it was found that an SNR of 40 dB or better was needed for an enhancement in the quality of the reconstructions. However due to the nonlinearity of the inverse scattering problem, this value may be dependent on the geometry of the configuration, the frequency of operation, and the particular algorithm used for the solution of the inverse scattering problem. (ii)The waveguide configuration has shown a higher potential of reconstruction than the camera configuration since this latter can only provide information about transmission through the OUT, resulting in a limited view of the scattered field. (iii)An efficient way to treat the inhomogeneous background of the breast compression based configurations based on the calculation of the inhomogeneous Green’s function has been presented. The computation of this Green’s function is performed numerically and its reduced size allows a much smaller cost both in terms of memory and speed of the solution of the inverse scattering problem. (iv)The reconstruction of simple breast-like models was performed with and without compression where only the fat portion of the breast was deformed, taking into account the results from a simplified mechanical simulation of the breast tissue deformation. An interesting study could be performed by using realistic numerical phantoms available in the literature and applying a nonlinear formulation for the mechanical deformation that allows large deformations.
The absence of painful breast compression is often cited as an advantage of radar and MT based techniques. Deciding whether the image enhancements obtained with compression and the associated potential for better tumor detection overcome the benefits in terms of the patient’s comfort is of course beyond the scope of this paper, as it would require detailed clinical studies. The results nevertheless suggest that compression should not be neglected in future MT developments, as a potential for better diagnostics.
Acknowledgments
The authors would like to thank Mohamed Aissiou, Professor Delphine Pernie-Curier, and Professor Martin Levesque from the Mechanical Engineering Department at Ecole Polytechnique for the useful discussions on the mechanical deformation of breast tissues. This work was performed with the financial support of the FQRNT and CRSNG.