Microwave Breast Imaging Using Compressed Sensing Approach of Iteratively Corrected Delay Multiply and Sum Beamforming

Microwave imaging (MI) is a consistent health monitoring technique that can play a vital role in diagnosing anomalies in the breast. The reliability of biomedical imaging diagnosis is substantially dependent on the imaging algorithm. Widely used delay and sum (DAS)-based diagnosis algorithms suffer from some significant drawbacks. The delay multiply and sum (DMAS) is an improved method and has benefits over DAS in terms of greater contrast and better resolution. However, the main drawback of DMAS is its excessive computational complexity. This paper presents a compressed sensing (CS) approach of iteratively corrected DMAS (CS-ICDMAS) beamforming that reduces the channel calculation and computation time while maintaining image quality. The array setup for acquiring data comprised 16 Vivaldi antennas with a bandwidth of 2.70–11.20 GHz. The power of all the channels was calculated and low power channels were eliminated based on the compression factor. The algorithm involves data-independent techniques that eliminate multiple reflections. This can generate results similar to the uncompressed variants in a significantly lower time which is essential for real-time applications. This paper also investigates the experimental data that prove the enhanced performance of the algorithm.


Introduction
Microwave imaging (MI) has been verified to be a reliable health monitoring technique that can play a prime role in diagnosing anomalies in breast tissue. A microwave transceiver can recognize small signal variations with the changes in the electrical properties of human tissues. Recent advances in microwave imaging show that it can be a viable solution as preclinical detection tools in human body imaging [1][2][3][4]. The on-phantom experiments shown in [4] describe the functionality of a microwave sensor for the successful localization of the target. The iterative reconstruction is referred to as post-processing of the acquired data. The clinical realization of MI is discussed with recent advancements in [3], where the MI technology is selected as a complementary modality of the existing screening systems for breast cancer. The potential role of MI and the clinical advances of this imaging is discussed in [1,2] where MI is described as a rapid and inexpensive method of diagnosis as compared to others.
The reliability of biomedical imaging-based diagnosis is substantially dependent on the imaging algorithm. There are three categories of MI, active, passive and hybrid [5]. Beamforming is the primary part of image reconstruction of active radar-based MI [6,7]. Among the two types of beamformers, data-independent beamformers take less time in image reconstruction while the data-adaptive beamformers have excessive computational complexity but produce better images [8].
The differential microwave diagnosis for breast [9], DMAS image reconstruction algorithm [10] and DAS algorithm have been suggested to improve imaging performance like the robust, interference-canceling, balanced deblurring. Several data-independent algorithms were utilized for microwave-based image processing, like Tissue Sensing Adaptive Radar (TSAR) [11] and confocal microwave diagnosis algorithm [12,13]. However, DAS and its variations are much more popular owing to their robust performance and simplicity. Most of these algorithms are applied to the received signal for the reconstruction of breast interior, feature extraction and post-processing to locate abnormalities. Another data-independent beamformer is DMAS which provides improved contrast and resolution in the reconstructed images by trading the complexity of computation. Several variants of DAS were introduced besides DMAS, including coherence factor DAS [14], improved DAS [15], and the iterative variant of DAS [16]. The improvements of DMAS over DAS has been proven in [13,17]. The iterative variant of DAS and DMAS is proposed in [18] where the authors removed the background clutter in the reconstructed images. The reconstructed images of the DMAS algorithm were used to detect the tumor after normalizing the maximum scatterer [19]. The compressed sensing approach has been reported in several image reconstruction techniques in MI [20,21]. A Bayesian CS approach to MI of an inhomogeneous object based on mesh discretization and electric field integral equation in the imaging region was proposed [22]. The CS was applied to separable surrogate functionals (SSF) optimization to reduce the sampling data volume in [23]. None of them is applied to the most popular DAS or variants of DAS algorithms. An iterative three-dimensional non-linear inverse scattering technique was proposed for the reduction of computational time using a nondispersive material model [24].
In this paper, a compressed sensing iterative variant of the DMAS algorithm is introduced for increasing the performance and reducing the runtime for microwave imaging applications using iteratively corrected DMAS beamforming. The proposed CS-ICDMAS works on the basis of computing all the channel power and eliminating comparatively lower power channels. This reduces the number of observations and saves the calculation time without decreasing the image quality. The practical imaging was tested, and the numerical analysis is also presented.

Imaging System Setup
The detailed flowchart of the research methodology is presented in Figure 1. The imaging system consisted of a total of 16 antenna arrays with a mechanical rotating platform. A set of heterogeneous phantoms was tested. The system setup and phantom fabrication procedures can be found in the literature [15]. The system comprised eight transmitters and eight receiver antennas. A rotating turntable held the antenna array and data were taken using a MATLAB-based image processing platform using a vector network analyzer (VNA). The turntable rotated 360 degrees at a step of 7.2 degrees that resulted in 50 equal datapoints. As a result, the total number of channels was 8 × 8 × 50. Figure 2 shows the imaging system setup while collecting data from the system in the presence of a breast phantom.

Antenna Performance Analysis
The antenna used in the proposed system is a modified Vivaldi antenna. The design evaluation of the antenna can be found in [17]. The antenna has dimensions of 77.72 mm × 60 mm and a bandwidth from 2.70-11.20 GHz. The reflection coefficient is considered in free space and the heavy load of the phantom in this work. The setup for taking measurements in the presence of a phantom is shown in Figure 3. Figure 4 illustrates the measured reflection coefficient in free space and a phantom load. It was observed that the antenna could firmly maintain its bandwidth in both scenarios with minimal distortions. The mutual coupling between the adjacent antennas was also investigated. Figure 5a shows the setup of Tx and Rx antennas considered for measuring the isolation between the adjacent antennas. The results are presented in Figure 5b It was observed that the antennas had sound isolation and were maintained between −26 dB and −60 dB.

Antenna Performance Analysis
The antenna used in the proposed system is a modified Vivaldi antenna. The design evaluation of the antenna can be found in [17]. The antenna has dimensions of 77.72 mm × 60 mm and a bandwidth from 2.70-11.20 GHz. The reflection coefficient is considered in free space and the heavy load of the phantom in this work. The setup for taking measurements in the presence of a phantom is shown in Figure 3. Figure 4 illustrates the measured reflection coefficient in free space and a phantom load. It was observed that the antenna could firmly maintain its bandwidth in both scenarios with minimal distortions. The mutual coupling between the adjacent antennas was also investigated. Figure 5a shows the setup of Tx and Rx antennas considered for measuring the isolation between the adjacent antennas. The results are presented in Figure 5b It was observed that the antennas had sound isolation and were maintained between −26 dB and −60 dB.

Beamforming Method
The proposed algorithm was analyzed based on the microwave signal contrast. The scattered microwave signals from healthy and unhealthy breast phantoms were compared. First, microwave signal contrast comparison took place between the reference microwave signal using numerical simulation of the full-wave time domain model and the scattered signal from the computational phantom. The successive approximation method was obtained to get the exact position and size of the tumor. There was no coupling medium between the phantom and the antenna; as a result, we needed to eliminate skin reflections. Thus, for artifact removal, the rotation subtraction was applied, which relied on the comparison between the original scattering and at least a single-rotated scattering [15]. The array of antennas is placed around the Region of Interest (ROI) in such systems. The array was rotated once, and the offset data were collected and recorded as a reference of initial radiance.

Channel Compression on Skin Reflection Removal
Different methods of artifact removal including rotation subtraction [25] and singular value decomposition (SVD) were proposed [26]. In the rotation subtraction method, it was assumed that the tumor was out of the axis of the rotation and the response received from the tumor would be well-preserved even after the subtraction of the scattered rotated position signals. In the SVD method, the breast shape is restricted as uniform and submerged in a coupling medium and the antenna array is maintained at a certain distance from the skin. In this work, the artifact removal was performed using the rotation subtraction method. The antenna array completed full rotation surrounding the phantom at a 7.2 • step with a total of 8 × 8 × 50 channels. The equally sliced datapoint was denoted Nϕ = 50. The system collected complex S-parameters in the frequency domain, S (f, t x , r x , ϕ), where f is the frequency, t x = 1, 3, 5, 7, 9, 11, 13, 15 and r x = 2, 4, 6, 8, 10, 12, 14, 16 are the transmitting and receiving antennas, respectively, and the alignment of the rotating platform angle is ϕ. A total of 201 datapoints was recorded for each rotation within the band from 2.7 to 8.0 GHz. Once the data were recorded for the original illumination, the array s rotated around the phantom to get offset data. The difference of the data from the different illuminations was computed to remove the effect of the skin.
In this work, the sample S (f, t x , r x , ϕ) was parted into two matrices based on ϕ becoming odd and even, or S odd (f,t x ,r x ,ϕ odd ), and S even (f,t x ,r x ,ϕ even ), individually, where ϕ odd = odd number of datapoints (1, 3, 5, . . . N ϕ -1) and ϕ even = even number of data points (2, 4,6, . . . N ϕ ). Thus, S odd can be considered initial illumination and S even -"offset" illumi- nation. Finally, rotation subtraction was implemented by simply calculating the difference between the two matrices according to Equation (1): The power of all the channels was calculated. The low-power channels were subtracted based on the compression factor Finally, rotation subtraction was implemented by simply calculatin tween the two matrices according to Equation (1): The power of all the channels was calculated. The low-power tracted based on the compression factor 03D7. The measured 8 × 8 sca presence of a breast phantom is illustrated in Figure 6. The differenc and low-power channels was visible. Afterwards, the total skin reflect were calculated based on Equation (2): ϗ × n (Rskin_removed) channels with the highest power were conside from Sskin_removed to produce compressed data as Scompressed(f,tx,rx, φodd), whe number of channels and ϗ is the compression factor of 0 ≤ ϗ ≤ 1. The lo were eliminated to reduce the execution time and make the process fa the skin reflections, the signals were converted to the time domain us rier transform to create Γ(t, tx, rx, φodd). Then, the data in the Γ(t, tx, rx, φ via the delay multiply and sum (DMAS) algorithm for the reconstr [27,28]. Normally, DAS-based methods struggle in noisy settings wh tions from various scattering sources need to be considered. In that c forms DAS by generating sharper contrast images by utilizing the corr However, an earlier study reported that DMAS struggles in high-nois in multiple target situations [30]. Therefore, in this research work, a pressed sensing (CS), we revised the typical DMAS algorithm by imp correction of the delay calculation. The correction was done by severa until convergence to the desired level of accuracy.

Compressed Sensing in DMAS
Three-dimensional Cartesian coordinates of all points inside th were represented in the i by 3 matrices, C, where i = 1,2,…, and I is points. Then, the I was generated from C by the I matrix, PC-C, conta . The measured 8 × 8 scattering matrix in the presence of a breast phantom is illustrated in Figure 6. The difference between the high-and lowpower channels was visible. Afterwards, the total skin reflection data (R skin_removed ) were calculated based on Equation (2): Finally, rotation subtraction was implemented by simply calculating the difference between the two matrices according to Equation (1): The power of all the channels was calculated. The low-power channels were subtracted based on the compression factor 03D7. The measured 8 × 8 scattering matrix in the presence of a breast phantom is illustrated in Figure 6. The difference between the highand low-power channels was visible. Afterwards, the total skin reflection data (Rskin_removed) were calculated based on Equation (2): ϗ × n (Rskin_removed) channels with the highest power were considered for the imaging from Sskin_removed to produce compressed data as Scompressed(f,tx,rx, φodd), where n(Rskin_removed) is the number of channels and ϗ is the compression factor of 0 ≤ ϗ ≤ 1. The lower power channels were eliminated to reduce the execution time and make the process faster. After adjusting the skin reflections, the signals were converted to the time domain using the inverse Fourier transform to create Γ(t, tx, rx, φodd). Then, the data in the Γ(t, tx, rx, φodd) were processed via the delay multiply and sum (DMAS) algorithm for the reconstruction of the image [27,28]. Normally, DAS-based methods struggle in noisy settings where multiple reflections from various scattering sources need to be considered. In that case, DMAS outperforms DAS by generating sharper contrast images by utilizing the correlation process [29]. However, an earlier study reported that DMAS struggles in high-noise settings as well as in multiple target situations [30]. Therefore, in this research work, after applying compressed sensing (CS), we revised the typical DMAS algorithm by implementing iterative correction of the delay calculation. The correction was done by several iteration processes until convergence to the desired level of accuracy.
× n (R skin_removed ) channels with the highest power were considered for the imaging from S skin_removed to produce compressed data as S compressed (f,t x ,r x , ϕ odd ), where n(R skin_removed ) is the number of channels and The power of all the channels was calculated. The low-power tracted based on the compression factor 03D7. The measured 8 × 8 sca presence of a breast phantom is illustrated in Figure 6. The differenc and low-power channels was visible. Afterwards, the total skin reflect were calculated based on Equation (2) ϗ × n (Rskin_removed) channels with the highest power were conside from Sskin_removed to produce compressed data as Scompressed(f,tx,rx, φodd), whe number of channels and ϗ is the compression factor of 0 ≤ ϗ ≤ 1. The lo were eliminated to reduce the execution time and make the process fa the skin reflections, the signals were converted to the time domain us rier transform to create Γ(t, tx, rx, φodd). Then, the data in the Γ(t, tx, rx, φ via the delay multiply and sum (DMAS) algorithm for the reconstr [27,28]. Normally, DAS-based methods struggle in noisy settings wh tions from various scattering sources need to be considered. In that c forms DAS by generating sharper contrast images by utilizing the corr However, an earlier study reported that DMAS struggles in high-nois in multiple target situations [30]. Therefore, in this research work, pressed sensing (CS), we revised the typical DMAS algorithm by imp correction of the delay calculation. The correction was done by severa until convergence to the desired level of accuracy. The power of all the chann tracted based on the compressio presence of a breast phantom is and low-power channels was vis were calculated based on Equati _ skin removed R ϗ × n (Rskin_removed) channels w from Sskin_removed to produce compr number of channels and ϗ is the were eliminated to reduce the ex the skin reflections, the signals w rier transform to create Γ(t, tx, rx, via the delay multiply and sum [27,28]. Normally, DAS-based m tions from various scattering so forms DAS by generating sharpe However, an earlier study repor in multiple target situations [30 pressed sensing (CS), we revised correction of the delay calculatio until convergence to the desired ≤ 1. The lower power channels were eliminated to reduce the execution time and make the process faster. After adjusting the skin reflections, the signals were converted to the time domain using the inverse Fourier transform to create Γ(t, t x , r x , ϕ odd ). Then, the data in the Γ(t, t x , r x , ϕ odd ) were processed via the delay multiply and sum (DMAS) algorithm for the reconstruction of the image [27,28]. Normally, DAS-based methods struggle in noisy settings where multiple reflections from various scattering sources need to be considered. In that case, DMAS outperforms DAS by generating sharper contrast images by utilizing the correlation process [29]. However, an earlier study reported that DMAS struggles in high-noise settings as well as in multiple target situations [30]. Therefore, in this research work, after applying compressed sensing (CS), we revised the typical DMAS algorithm by implementing iterative correction of the delay calculation. The correction was done by several iteration processes until convergence to the desired level of accuracy.

Compressed Sensing in DMAS
Three-dimensional Cartesian coordinates of all points inside the region of interest were represented in the i by 3 matrices, C, where i = 1,2, . . . , and I is the total number of points. Then, the I was generated from C by the I matrix, P C-C , containing the Euclidian distances between each possible pair of points in the imaging domain. Matrices A Tx and A Rx contained the optimized three-dimensional Cartesian coordinates of the transmitting and receiving antennas, respectively. Since the imaging domain was stationary, the antennas changed their distance from the points to be reconstructed as they rotated. Thus, A Txϕodd and A Rxϕodd were generated by determining all the antenna positions considered in the original orientation. Then, P Txϕodd-C and P C-Rxϕodd , containing the distances from each point to the transmitting and receiving antennas, were calculated from C, A Txϕodd , and A Rxϕodd . Finally, the delays needed for focusing each channel to focus on point i in C were evaluated by summing the distances of the transmitting and receiving antenna positions concerned to the point being focused. The wave propagation velocity in the background medium, air in this study (theoretically, speed of light c), was divided by the total distance to yield the appropriate delay, τ (i, rx, ϕ odd ), as indicated in Equation (3): where ε b is the dielectric constant of the coupling medium. The delay τ was calculated from the estimated shortest distance and the reflected signal from C(i), which are the straight-line Euclidean distances calculated previously. The delays were further added to the signals for delivering a proper delayed signal. After multiplying by the paired delayed signal, they were summed for determining the scatter intensity at the allotted point in the ROI as shown in Equation (4). The paired signal multiplication compensated coherent signals together with higher values and hence improved the performance of DMAS over DAS [15,27].

Iterative Correction of CS-DMAS
The dielectric materials diminish the propagation speed; this results in increased time delay. Thus, the higher value of Γ in the C region can be implied as to the higher dielectric region. The additional time adjustment can be made by properly increasing the distance while calculating τ. Here, we introduce an iterative approach to estimate the best-fitted delay as well as the scattered intensity map (SIM) evaluation. At the same time, a greater estimation of SIM is dependent on the proper adjustment of τ. However, the iterative process can lead to vulnerability and noise-sensitivity if we use Γ directly. An inverse distance weighted integral averaging was introduced to produce a smoothed SIM, Γ'(i). This method can be considered analogous to applying the 3D Green's function in the distorted Born iterative method.
An extra time, P C-C(i,j) , was adjusted by appropriately increasing the distances. Then, the following equation was used to calculate modified delay: Here, c is the speed of light.
The new set of delays was used to reconstruct the SIM. Finally, the completion benchmark was inspected for convergence. The iteration of Equations (7) and (8)  The iteration was terminated when E Γ decreased the expected level of accuracy after achieving convergence. In this paper, E Γ was < 10 −5 .

Imaging Results
The simulation of DMAS was performed on an HP Workstation Mid 2018 with a 3.7 GHz Intel Xeon E3-1245 processor and 32 GB 2400 MHz DDR4 RAM.

Contour Images
Two sets of the lab-made heterogeneous phantom were used to evaluate the performance of the proposed beamformer. The first one had a single tumor inside (phantom 1) and another one was designed with two tumors (phantom 2). The heterogeneous phantoms 1 and 2 had four layers of skin, fat, gland and tumor to be as realistic as possible [31]. The power of all the channels tracted based on the compression fa presence of a breast phantom is illu and low-power channels was visible were calculated based on Equation ϗ × n (Rskin_removed) channels with from Sskin_removed to produce compress number of channels and ϗ is the com were eliminated to reduce the execu the skin reflections, the signals wer rier transform to create Γ(t, tx, rx, φod via the delay multiply and sum (D [27,28]. Normally, DAS-based meth tions from various scattering sourc forms DAS by generating sharper co However, an earlier study reported in multiple target situations [30]. T pressed sensing (CS), we revised th correction of the delay calculation. T until convergence to the desired lev The power of all the channels was calculated. The low-power channels were subtracted based on the compression factor 03D7. The measured 8 × 8 scattering matrix in the presence of a breast phantom is illustrated in Figure 6. The difference between the highand low-power channels was visible. Afterwards, the total skin reflection data (Rskin_removed) were calculated based on Equation (2): ϗ × n (Rskin_removed) channels with the highest power were considered for the imaging from Sskin_removed to produce compressed data as Scompressed(f,tx,rx, φodd), where n(Rskin_removed) is the number of channels and ϗ is the compression factor of 0 ≤ ϗ ≤ 1. The lower power channels were eliminated to reduce the execution time and make the process faster. After adjusting the skin reflections, the signals were converted to the time domain using the inverse Fourier transform to create Γ(t, tx, rx, φodd). Then, the data in the Γ(t, tx, rx, φodd) were processed via the delay multiply and sum (DMAS) algorithm for the reconstruction of the image [27,28]. Normally, DAS-based methods struggle in noisy settings where multiple reflections from various scattering sources need to be considered. In that case, DMAS outperforms DAS by generating sharper contrast images by utilizing the correlation process [29]. However, an earlier study reported that DMAS struggles in high-noise settings as well as in multiple target situations [30]. Therefore, in this research work, after applying compressed sensing (CS), we revised the typical DMAS algorithm by implementing iterative correction of the delay calculation. The correction was done by several iteration processes until convergence to the desired level of accuracy.

Numerical Imaging Data
The numerical assessments in terms of signal-to-mean ratio (SMR) and execution time is described in this section. The SMR states the proportion of the maximum tumor energy to the mean energy from backscattered signals at a similar sample [32].
Numerical imaging data in terms of SMR and execution time are shown in Table 1. The graphical representation of the signal to mean ratio and execution time against the compression factor is presented in Figure 9a,b. The SMR of CS-ICDMAS remains decreased with the compression factor. The SMR shows a proportional relation with the execution time. But the drop of SMR concerning execution time is not huge. The proposed algorithm shows exponential improvement in reducing execution time by maintaining image quality and detection. These results validate the capability of CS-ICDMAS to decrease computation time and reduce computational complexity by eliminating the weak channels by preserving detection performance. A comparison of the execution time for the DAS, CF (Coherence Factor)-DAS, DMAS and ICDMAS with the proposed method is provided in Table 2. It is noted that the proposed algorithm can reduce the time required for computations and simultaneously gives similar imaging output.
The power of all the channels was calculated. The low-power channels were subtracted based on the compression factor 03D7. The measured 8 × 8 scattering matrix in the presence of a breast phantom is illustrated in Figure 6. The difference between the highand low-power channels was visible. Afterwards, the total skin reflection data (Rskin_removed) were calculated based on Equation (2): ϗ × n (Rskin_removed) channels with the highest power were considered for the imaging from Sskin_removed to produce compressed data as Scompressed(f,tx,rx, φodd), where n(Rskin_removed) is the number of channels and ϗ is the compression factor of 0 ≤ ϗ ≤ 1. The lower power channels were eliminated to reduce the execution time and make the process faster. After adjusting the skin reflections, the signals were converted to the time domain using the inverse Fourier transform to create Γ(t, tx, rx, φodd). Then, the data in the Γ(t, tx, rx, φodd) were processed via the delay multiply and sum (DMAS) algorithm for the reconstruction of the image [27,28]. Normally, DAS-based methods struggle in noisy settings where multiple reflections from various scattering sources need to be considered. In that case, DMAS outperforms DAS by generating sharper contrast images by utilizing the correlation process [29]. However, an earlier study reported that DMAS struggles in high-noise settings as well as in multiple target situations [30]. Therefore, in this research work, after applying compressed sensing (CS), we revised the typical DMAS algorithm by implementing iterative correction of the delay calculation. The correction was done by several iteration processes until convergence to the desired level of accuracy.

Compressed Sensing in DMAS
Three-dimensional Cartesian coordinates of all points inside the region of interest were represented in the i by 3 matrices, C, where i = 1,2,…, and I is the total number of

Conclusions
This work introduces a new compressed sensing approach of iteratively corrected delay multiply and sum beamforming over the well-known DMAS in terms of microwave imaging of the breast. The main disadvantage of conventional DMAS is the high computation load and runtime. An array of 16 antennas was used to receive the backscattered microwave signal from a breast phantom. The compressed sensing approach was proposed by finding low-power channels across all the observations and eliminating them based on the compression factor. Then, an iterative approach was applied till convergence. The artifact removal was performed using the rotation subtraction method. Imaging data from the experimental system were tested with the proposed approach and the reconstructed image and numerical analysis were presented. The proposed CS-ICDMAS algorithm shows exponential improvement in reducing execution time as compared to DAS, CF-DAS, DMAS and ICDMAS by maintaining image quality and detection which validate the capability of CS-ICDMAS to decrease computation time and reduce computational complexity by eliminating the weak channels and preserving detection performance.
Institutional Review Board Statement: Not applicable.

Conflicts of Interest:
The authors declare no conflict of interest.