Microwave breast cancer detection using time-frequency representations

Abstract

Microwave-based breast cancer detection has been proposed as a complementary approach to compensate for some drawbacks of existing breast cancer detection techniques. Among the existing microwave breast cancer detection methods, machine learning-type algorithms have recently become more popular. These focus on detecting the existence of breast tumours rather than performing imaging to identify the exact tumour position. A key component of the machine learning approaches is feature extraction. One of the most widely used feature extraction method is principle component analysis (PCA). However, it can be sensitive to signal misalignment. This paper proposes feature extraction methods based on time-frequency representations of microwave data, including the wavelet transform and the empirical mode decomposition. Time-invariant statistics can be generated to provide features more robust to data misalignment. We validate results using clinical data sets combined with numerically simulated tumour responses. Experimental results show that features extracted from decomposition results of the wavelet transform and EMD improve the detection performance when combined with an ensemble selection-based classifier.

Appendix:

1.1 A.1 Tumour response simulation

We adopted the strategy outlined in [24] to simulate tumour responses for each volunteer, based on the transmitted pulses from the antennas and the dielectric properties of breast tissue. For one antenna pair and a tumour position p ₀, the frequency domain representation of the tumour response R ^t(p ₀,ω) was modelled as:

$$ R^{t}(p_{0},\omega)={\Gamma} R(\omega)e^{-j(k_{im}(d^{t}_{im}-d^{d}_{im})+k_{br}(d^{t}_{br}-d^{d}_{br}))}, $$

(9)

where R(ω) is the frequency domain representation of received signal, $d^d_{im}$ and $d^d_{br}$ are the lengths of the direct path for this antenna pair through the immersion medium (ultrasound gel) and breast tissue, respectively. $d^t_{im}$ and $d^t_{br}$ are the lengths of the shortest path between the antenna pair via the tumour position p ₀, in the immersion medium and the breast tissue, respectively. Γ is a constant that can be used to introduce additional attenuation in the tumour response. In this paper, we concentrate on the case Γ = 1. k _{i m} and k _{b r} are the wavenumbers for the immersion medium and breast tissue, respectively, and these have the following expressions:

$$\begin{array}{@{}rcl@{}} k_{im}(\omega) &=& \frac{2\pi}{\lambda_{im}(\omega)} =\sqrt{\epsilon_{im}(\omega)}\frac{\omega}{c}\,, \end{array} $$

(10)

$$\begin{array}{@{}rcl@{}} k_{br}(\omega) &=& \frac{2\pi}{\lambda_{br}(\omega)}= \sqrt{\epsilon_{br}(\omega)}\frac{\omega}{c}\,. \end{array} $$

(11)

Here, 𝜖 _{i m} is the relative permittivity of the immersion medium, c is the speed of light, and 𝜖 _{b r} is the average breast tissue relative permittivity. The latter is specified by the Debye model [15]:

$$ \epsilon_{br}(\omega)=\epsilon_{\infty}+\frac{\Delta\epsilon}{1+j\omega\tau}+\frac{\sigma_{s}}{j\omega\epsilon_{0}}\,, $$

(12)

where 𝜖 ₀ = 8.854 × 10⁻¹² F/m is the permittivity of free space, and 𝜖 _∞ is the dielectric constant of the material at infinite frequency. Δ𝜖 = 𝜖 _s − 𝜖 _∞ , and 𝜖 _s is the static dielectric constant. The pole relaxation constant is τ and the static conductivity is σ _s . The model parameters are chosen to approximate the dielectric properties of breast tissue.

We refer readers to [24] for a more complete description and discussion of the tumour response simulation procedure. We notice that when Γ = 1 and the tumour is close to the direct path between the transmitting antenna and the receiving antenna, the amplitude of simulated tumour response may be too large. We introduce a hard threshold to constrain the ratio of the amplitude between the tumour response and the received signal.

We do stress that this is not a model that assumes homogeneous breast tissue. By including the received signal R(w), rather than the transmitted signal, the model incorporates distortions and delays caused by inhomogeneous tissue. In Fig. 11, we present one example of a tumour response-injected signal.

Fig. 11

An example of received pulses before and after tumour response injection. The pulse is collected by the antenna pair A1A2 from the left breast of volunteer one during her first visit

1.2 A.2 Parameter values

As discussed in [24], signals collected by the antenna pairs located on the opposite sides of the breast can be highly distorted and vary significantly among different volunteer visits. Thus, we discarded signals from any antenna pair whose median peak amplitude of all the training data is less than a threshold of 20 mV. The length of original signal is 1024 and the sample rate is 40GHz [31]. We first windowed data from 60 samples before the peak amplitude of measurement from antenna pair A1A2 measurement. Since antenna A1 and A2 are one of the closest antenna pairs, we do not expect any tumour responses to occur before the peak amplitude of the measurement from A1A2. Scans from the same antenna pair were then aligned based on maximal correlations with a reference scan and were windowed between the 61^st sample to 604^th sample, corresponding to an actual time period of 13.6 ns. A window of 544 samples allows the discrete wavelet transform of 5 levels and covers all regions of significant tumour responses as shown in [24]. The maximum ratio of the amplitudes between the simulated tumour response and the received signal was set to 0.3.

The candidate values of the 2ν-SVM hyper-parameters used for cross validation are listed in Table 4. We tested the detection performance with the γ value chosen from the candidate set γ = {2⁻¹⁵, 2⁻¹³,..., 2⁵} using a small subset of the data, and observed that the ensemble classifier almost always chose the γ values of {2⁻¹, 2¹}. We further tested a range of fixed γ values between 2⁻² and 2² and observed that they lead to similar performance. To reduce the computational cost during training, we set γ = 1 for the experiments detailed here.

Table 4 Candidate values of hyper-parameters of 2ν-SVM in cross validation

Thus, there are 1 × 18 × 18 = 324 different 2ν-SVM hyper-parameter combinations. These were used to produce a model library consisting of M × 324 base models, where M is the number of retained antenna pairs in each data set. The ensemble classifier selected 100 base models, choosing those with the smallest Neyman-Pearson measure when applied to the training data, to perform classification on the test data.

The mother wavelet used in the DWT was the Daubechies 5 (db5) [12]. The filter used in the first step of the DTCWT was the Farras filter [1], and Kingsbury Q-shift filters[19] with six taps were used for subsequent stages. 5-level decomposition was performed for DWT and DTCWT as we observed that the magnitudes of the decomposition results at higher levels were too small to contain useful information. This leads to 24 features for each antenna pair as we compute statistics for outputs from both high-pass filter and low-pass filter in the last decomposition level. We performed 5-level EMD decomposition for the same reason, and this leads to 20 EMD-based features for one antenna pair.

The magnitudes of the features derived from decomposition results from higher levels were often significantly smaller in magnitude. For the measurements which were collected by the same antenna pair, we rescaled each feature of the training data into the range [0, 1]. The testing data was scaled by the same ratio. The process was repeated for each antenna pair.