A Feasibility Study of 2-D Microwave Thorax Imaging Based on the Supervised Descent Method

In this paper, the application of the supervised descent method (SDM) for 2-D microwave thorax imaging is studied. The forward modeling problem is solved by the finite element-boundary integral (FE-BI) method. According to the prior information of human thorax, a 3-ellipse training set is generated offline. Then, the average descent direction between an initial background model and the training models is calculated. Finally, the reconstruction of the testing thorax model is achieved based on the average descent directions online. The feasibility using One-Step SDM for thorax imaging is studied. Numerical results indicate that the structural information of thorax can be reconstructed. It has potential for real-time imaging in future clinical diagnosis.


Introduction
Microwave plays an important role in non-invasive detections, such as geophysics explorations and industrial monitoring. Its capability in biomedical imaging is developing gradually. Compared with present clinical imaging modalities including computed tomography (CT), magnetic resonance imaging (MRI), and ultrasound, microwave takes advantage of non-ionizing radiation, portable, and cost-efficient devices. As an emerging biomecial imaging approach, microwave has already been proven to be effective in breast cancer detection [1,2]. Its application for brain imaging has also been studied [3][4][5][6][7]. Microwave can also be used for thorax imaging and monitoring the human respiratory system. The respiratory system is of great importance to human health. Human thorax is a complex environment as the structure is complicated, consisting of lungs, heart, ribs, and thoracic vertebrae. Various physiological functions are completed in the thorax. Some previous work focused on vital sign capturing [8,9], detections of heart failure [10] and pulmonary edema [11], and the measurement of lung water [12][13][14] and cancer [15][16][17]. Various systems [18,19] and sensors [16,[20][21][22][23][24][25][26] are proposed for experimental and clinical use. The feasibility of microwave respiration monitoring has already been studied [27] and sensor antennas for corresponding applications have been proposed [28,29]. Due to the complexity of the thorax, a complete and accurate microwave thorax imaging system is a challenging task, of which researchers are addressing.
Microwave biomedical imaging can be formulated as an inverse problem. It is usually ill-posed and nonlinear. The ill-posedness and nonlinearity are exacerbated by the complex environment in the human thorax. Researchers applied the synthetic radar imaging method in detection of human torso fluid [19] and brain stroke [30]. The Born iterative method and Gauss Newton method have also been adopted in brain stroke detection [7,31]. Some other optimization approaches are also applied to biomedical inverse problems, which are usually computationally expensive.
Machine learning techniques are rapidly developing in recent decades, which can largely improve the imaging quality by integrating prior information and exploiting more connections between data and images during the learning stage. Deep learning networks can be qualified for microwave image reconstruction with a better performance than the distorted Born iterative method and phase confocal method [32]. A convolutional neural network (CNN) is used to learned the nonlinearity in breast tumour detection [33]. In [34], a CNN is added after the non-linear inversion stage to achieve super resolution imaging for breast phantom imaging. Deep neural networks can also improve the performance of tumour classification in microwave biomedical imaging [35]. Learning-by-examples strategies have proven to be successful in the prediction of lung dimensions from electrical impedance tomography [36]. To the authors' limited knowledge, there is not much research being done on the application of machine learning techniques to microwave thorax imaging.
The supervised descent method (SDM) is a learning-based technique. Average descent directions of cost functions are learned from a series of training samples and then used to update the model and minimize the data misfit between measured and simulated data. SDM is widely applied to face alignment [37]. Present work shows its feasibility in microwave imaging [38][39][40]. SDM is able to reconstruct structures accurately in data sparse case by recovering coefficients of compactly supported radial basis functions, which is shown in [41]. SDM has also been applied to 2-D magnetotelluric inversion [42,43]. The performance of SDM in pixel-based inversion for transient electromagnetic data is studied as well [44]. Moreover, SDM has been used in electrical impedance tomography thorax imaging [45].
In this paper, the application of SDM for 2-D microwave thorax imaging is discussed, which is the first time for SDM applied in this field. The formulations of forward problem and SDM are given. Then, numerical experiments based on thorax models are conducted and discussions are made according to numerical results of iterative SDM. In order to achieve a high temporal resolution, the capability of One-Step SDM is tested. Its potential for real-time thorax imaging is analyzed. The influence of different training sets is given as well. In the experiments, 433 MHz and 915 MHz are adopted as they both belong to the industrial, scientific, and medical band (ISM). This paper is organized as follows. Formulation of SDM is given in Section 2. Numerical experiments including the training set and testing thorax models are described in Section 3. Conclusions are given in Section 4.

Formulations
A 2-D inverse scattering problem can be depicted as Figure 1. Transmitters and receivers are located around the domain of interest (DoI) D, in which a 2-D object is embedded. The forward problem is solved by the finite element-boundary integral (FE-BI), which inherits advantages of the finite element method in dealing with complex medium and integral equations in dealing with propagation of electromagnetic wave in free space. With transverse magnetic (TM) sources, this problem can be formulated as: where E z and J z represent the z-component of electric field and current density respectively. G 0 ( ρ, ρ ) denotes the Green's Function in free space. k( ρ), j = √ −1, ω and µ denote the wavenumber, imaginary unit, angular frequency, and permeability in D respectively. Γ denotes the boundary of D. To solve the forward problem, D together with its boundary is firstly discretized into triangular meshes. E z in D and its normal derivative ∂E z /∂n on the boundary of D can be expressed as: where V and V b denote the number of nodes in and on the boundary of D. Furthermore, E I n denotes the value of the z-component of electric field at node n in D and N I n ( ρ) denotes the corresponding interpolation function; E b n and E b n denote the value of the z-component of electric field and its normal derivative at node n on the boundary of D; and N b n ( ρ) denotes the corresponding interpolation function. The E I n , E b n , and E b n are unknowns to be solved. Then multiplied with a testing function w m ( ρ) and integrated over D, a weak-form representation of Equation (1) can be obtained as: Using the interpolation functions as testing functions, namely Galerkin's method, Equation (5) can be transformed into (V + V b ) equations. Similar steps can be applied to Equation (2) as: Another V b equations can be obtained with point matching testing functions. After substituting Equations (3) and (4) into Equations (5) and (6), these equations can be combined into a compact matrix equation as: where A 1 and A 2 denote the coefficients calculated from Equations (5) and (6) respectively. E and E are column vectors consisting of the unknowns in Equations (3) and (4), namely E I n , E b n , and E b n . Finally, E z can be calculated by solving the matrix equation. The above forward modeling process can be compactly written as: where d denotes the field measured by receivers, m denotes the distributions of dielectric properties in DoI, and F(·) denotes the operator of forward modeling. Solving an inverse problem can be considered as an optimization procedure to find out a reconstruction model m that can minimize the data misfit between the measured data d obs and the simulated data d. The cost function can be defined as: Based on the Taylor expansion, Equation (9) at m + ∆m can be written as: where J C and J F are the Jacobian matrices of C(m) and F(m) respectively. ∆m is an update to the model. The minimum of C is reached when the gradient of C with respect to ∆m is equal to zero, that is: so the update of the model based on Equation (11) can be obtained as: which defines the Gauss Newton method. K sets the model update direction together with the data misfit. As the inverse problem of thorax imaging is always nonlinear due to multiple scattering, iteration process is adopted to realize the inversion. In the Gauss Newton method, the update direction needs to be calculated during each iteration, which is usually time consuming. In order to avoid excessive time cost, SDM provides an alternative approach to calculate the update direction K offline in advance. Firstly, a number of training samples m j is generated with corresponding forward modeling data d j . Then, an average update direction K t is calculated between a certain initial model and the training models.
In each iteration, a corresponding average update direction is learned. The average update direction K t in the t-th iteration can be calculated as: where N is the number of training samples and: where m t j is the updated model of the j-th training sample and d t j is the corresponding forward modeling data in the t-th iteration. It can be transformed into a matrix form as: ∆M t together with ∆D t is an N-Row matrix in which each row corresponds to a training sample. In order to solve Equation (16), a least square method is adopted and a regularization procedure is introduced to stabilize the solution in case of the singularity of ∆D t : where α can be set proportional to the largest eigenvalue of ∆D t . After the update direction K t is obtained, the model can be updated by: The model misfit and data misfit can defined as: After the training stage, a set of update directions K 1 , K 2 , K 3 ... can be used for online imaging and the reconstruction can follow a similar process as Equation (18). During the inversion process, the data misfit should be checked in each iteration and the reconstruction could be stopped if the data misfit is less than a predefined criteria. But it should also be noted that this iterative approach can not guarantee an adequate temporal resolution. To solve this problem, a compromised solution is to use One-Step SDM, which means only the average descent direction of the first step is learned. During reconstruction, a one-step inversion is adopted. This One-Step SDM is a qualitative imaging method, yet it can still incorporate sufficient prior information through training stage and achieve a reasonable reconstruction, with a largely saved reconstruction time.

Description of Thorax Model
The goal of thorax imaging is the inversion of the distributions of dielectric properties in the human thorax. The main task is to reconstruct the lungs and heart. For simplicity, a homogeneous background is used as an initial model in the following discussions, which is shown in Figure 2. 2-D delta sources are placed on the boundary between free space and muscle. As it is infinitely long in the z-direction, the electric field in is the TM mode. The dielectric properties of human tissues including muscle, lung, and heart are given in [46]. The relative permittivity and conductivity of muscle, inflated lung, and heart at 433 MHz and 915 MHz are shown in Table 1. The dielectric properties of muscle are used for background. The numerical experiments including forward modeling and inversion are all carried out with Matrix Laboratory (MATLAB R2017b). When solving the forward problem, the model is partitioned into 4862 triangular elements with 150 elements on the boundary, with a maximum side length less than 1 cm, which is 0.013 of the wavelength of 433 MHz in vacuum.

Training Set
Training samples are usually set according to prior information. As lungs and heart are the main reconstruction targets in thorax imaging, a 3-ellipse training model is adopted in this application. Two random ellipses represent two lungs and one random circle in the middle represents the heart. The position, size, dielectric properties and rotation angle of different, training samples vary in a predefined range. Figure 3 shows an example of training samples at 433 MHz.  Table 2. It should be mentioned that the major axis of the two ellipses is in y-direction before rotation. So the rotation angle describes the angle between the major axis and y-direction. It is positive if the rotation is clockwise. Moreover, in clinical diagnosis, abnormalities may occur in only one of the two lungs, hence they may have different parameters. So the two ellipses share same variation range but the values of their parameters are different. SDM takes the advantage of prior information, so 3000 training samples are generated according to the rules above. Figure 4 shows the training misfit at two frequencies. Both converge after 10 iterations using SDM.

Iterative SDM
After the training stage, a thorax model is used to test the algorithm, which is shown in Figure 5a,d. The two ellipses represent two lungs and the circle in the middle represents the heart. They are in different colors as they have different dielectric properties. During reconstruction, 5% white noise is added to the testing data. The reconstructed relative permittivity distributions of the thorax model at 433 MHz are shown in Figure 5b,e. The red and yellow parts in Figure 5b represent reconstructed size and the real part of relative permittivity of lungs and hearts respectively. The orange and black parts in Figure 5e show the reconstructed imaginary part of relative permittivity. The data misfit of reconstruction process is shown in Figure 6. With 10 iterations, a thorax model is reconstructed and the reconstruction process converges. The position, size, and rotation of the lungs and heart can be reconstructed accurately. The reconstruction at 915 MHz is also tested, which is shown in Figure 7.

One-Step SDM
As mentioned before, the iteration process is adopted to overcome the nonlinearity in thorax imaging, which can be computationally expensive and difficult to implement for real-time imaging. In clinical diagnosis, respiration is continuous, which means the dielectric properties of lungs are changing all the time. As a result, thorax imaging with a high temporal resolution is necessary. The convergence of data misfit during training and testing shows that the data misfit decreases significantly at the first iteration, which means keeping only the first step of SDM may also have a reasonable performance. One-step SDM retains the advantage of iterative SDM in the incorporation of prior information and saves reconstruction time. Figure 5c,f show the reconstruction results using one-step SDM. Compared with iterative SDM, a similar performance of a reconstruction size, position, and rotation using one-step SDM is achieved. Lungs and heart can be distinguished from the reconstruction result. Although the value of reconstructed dielectric properties using one-step SDM are not so satisfying, the structural information of lungs and heart it provides can help doctors to monitor human respiration in clinical diagnosis, such as detecting level of atelectasis. Similar phenomenon can be seen for 915 MHz. The above results indicate that microwave real-time throax imaging may be feasible with the help of one-step SDM in the future.

Discussions
According to the results in numerical experiments, the feasibility of SDM on 2-D microwave thorax imaging is basically verified. With a training set containing models of three ellipses to mimic lungs and hearts, average descent directions can be learned for the inversion of a testing thorax model. SDM shows great capability in solving nonlinear inverse problem, which also shows prior information can play a significant role in inversion. Training sets need to be specially designed to achieve a better performance. In common applications, training samples can be generated simply according to possible shapes, sizes, positions, and dielectric properties of the real model. However, in thorax imaging, the relative permittivity and conductivity of heart are larger than those of background, while the opposite is true for lungs. The average descent direction for heart and lungs is different. As a result, the dielectric properties of the circle and ellipses are set to be larger and smaller than background respectively. Moreover, it should be noted that, the heart always lies between the two lungs in the front of the thorax. So the position of circles is set accordingly. More importantly, the overlap of position variation range between circle and ellipses should be relatively small to avoid learning and reconstruction mistakes in dielectric properties. If the ellipse is so large that it exceeds the thorax boundary, the exceeding parts should be removed.
The reconstruction results using 433 MHz and 915 MHz are not the same. The outline of heart and lungs is described more completely using low frequency. On the other hand, the separation between two lungs is clearer using 915 MHz, which means high frequency reconstruction focuses more on details. Reconstruction using two frequencies simultaneously may recover both outline and details at the same time. The results of dual-frequency one-step SDM are shown in Figure 8. Compared with single-frequency reconstruction, the structural information is more precise including both outline and details. Microwave at two frequencies shares the same structure during training and testing. As a result, inversion with dual-frequency can better reconstruct structural information than single-frequency reconstruction. In future, data with more frequencies may further improve the inversion.
This work mainly focuses on the reconstruction of dielectric properties changes in thorax caused by respiration. However the lung volume changes caused by thorax expansion and contraction during the respiration process may also affect the received signals. These factors should be further studied quantitatively in future studies. (c) (d)

Conclusions
In this paper, the application of one-step SDM in 2-D microwave thorax imaging was studied. A training set consisting of three ellipses was introduced to train the average update direction. After training stage, a numerical thorax model was used to test the algorithm. The results of numerical experiments verified its feasibility in thorax imaging. The results based on one-step SDM showed its potential to realize thorax imaging with a high temporal resolution, retaining adequate structural information for clinical diagnosis. In future, human respiration experiments will be conducted and one-step SDM will be adopted in experiments. The reconstruction of dielectric property distribution in human thorax will be further tested based on measured data using SDM.

Abbreviations
The following abbreviations are used in this manuscript: