† Corresponding author. E-mail:
Hyperthermia has proven to be beneficial to treating superficial malignancies, particularly chest wall recurrences of breast cancer. During hyperthermia, monitoring the time–temperature profiles in the target and surrounding areas is of great significance for the effect of therapy. An ultrasound-based temperature imaging method has advantages over other approaches. When the temperature around the tumor is calculated by using the propagation speed of ultrasound, there always exist overshoot artifacts along the boundary between different tissues. In this paper, we present a new method combined with empirical mode decomposition (EDM), similarity constraint, and continuity constraint to optimize the temperature images. Simulation and phantom experiment results compared with those from our previously proposed method prove that the EMD-based method can build a better temperature field image, which can adaptively yield better temperature images with less computation for assistant medical treatment control.
Hyperthermia is proved to be beneficial to treating chest wall recurrences of breast cancer.[1–4] Researches have proved that the combination of hyperthermia and other medical treatments, especially radiotherapy in local control[5,6] and chemotherapy,[7,8] can yield better medical effects. In order to achieve the best treatment effect in hyperthermia, temperature ranges in the target area are typically tuned close to 42 °C or 43 °C ± 1 °C. Therefore, precisely monitoring temperature is of great importance to assure safe medical treatment.
Invasive temperature measurement methods have been used primarily by embedding heat-sensitive sensors in the tissues. These methods give accurate local results, low cost, and easy operation. However, spatial sampling is too limited to make corrections for most of the tissue heterogeneities in a complex tissue like breast and it brings pain and risk of infection to the patient.[9] Common medical measurement methods with non-invasive attributes have previously been studied including: magnetic resonance (MR),[10,11] ultrasound,[12,13] and microwave[14] technologies. The three-dimensional (3D) temperature imaging with MR methods based on the diffusion coefficient,[15] the relaxation time, or proton resonance frequency (PRF)[16] of tissue has also been well reported. A hybrid PRF approach to simultaneously measuring time variations of PRF shift temperatures in water-based tissue and in fat-based tissues is proposed,[17] which makes temperature changes accurately monitored in both water-based and fat-based tissues. However, it costs quite a lot in the cycle of multiple treatments when the MR method is reused. While a temperature measurement method based on microwaves has proved that the microwave can penetrate into tissue in only a limited depth and can be easily influenced by the environment.[18] Ultrasound is a convenient, non-ionizing, and inexpensive method with less calculation requirements for temperature monitoring. The methods of using ultrasound as a temperature monitor involve measuring the sound attenuation coefficient,[19,20] measuring echo-shift deviation in tissue thermal expansion,[21,22] and measuring the backscattered energy difference of tissue in homogeneities.[23] Several research centers in the world are working on constructing a prototype of ultrasonic transmission tomography for women's breast examination. In our experiment, a clinical ultrasound ring array scanner for breast cancer diagnosis was built to record both reflected ultrasound energy and transmitted ultrasound energy.[24,25] However, the prior information about propagation speed of ultrasound and coefficient of thermal expansion are necessary for the echo-shift method.[26]
Maintaining uniform temperature over one or more large target regions requires the compensation for a number of factors. The Pennes bioheat equation[27,28] allows fairly good calculation provided that each of these variables is reasonably a function of time, but it is not so easy to determine the values of these factors for hyperthermia treatments and measurements as required. Some of these factors, particularly perfusion and absorption coefficients are substantially different in different tissues. Therefore, temperature often changes rapidly near tissue boundaries. This phenomenon will inevitably be wrongly reflected as overshoot artifacts around tissue boundaries, such as the greatly different thermal coefficients of speed of ultrasound between fat and muscle or connective tissues, when building temperature images with ultrasonic transmission tomography.
Researches on reducing overshoot artifacts normally consider the root cause of the artifacts and then perform an inverse operation.[29,30] Our previous contribution[31] proposed a method based on a physical mode. When applying our previous method, we solve a differential equation based on heat conduction rules with boundary conditions. Although the result is satisfactory for future biomedical treatment, it is characterized with a heavy computation burden and lack of the adaptive property.
In this paper, an alternative method of combining de-convolution, continuity constraint and empirical mode decomposition (EMD)[32,33] is proposed to build a temperature image based on ultrasonic transmission speed tomography, which can yield better border information adaptively with less calculation. The proposed method can provide accurate and fast measurements for assistant medical treatment control by effectively alleviating the negative effect from overshoot artifacts in temperature images. The rest of this paper is organized as follows. In Section
EMD, introduced by Huang et al.,[33] is an excellent data analysis method for non-stationary signals. The kernel of EMD is to decompose a signal into a group of intrinsic mode functions (IMF), which are the fundamental representing of the signal. As the decomposition process is adaptive, the EMD is especially suitable for analyzing signals with nonlinear or non-stationary properties.[34] An IMF must satisfy two requirements. The first one is that the number of zero crossings and the number of extreme points should differ by no more than one. The second one is the local average, which is normally equal to the mean of the maximum and minimum envelopes around the signal, obtained by using cubic spline interpolation through the corresponding local extrema, and must be zero.
With an ultrasound transmission speed image and the profiles between the speed of sound and temperature for different homogeneous tissues, we can build a basic temperature image. As for each homogeneous tissue, which can be segmented from the speed of sound image, we can calculate temperature from these profiles directly. While for the mixture area between two tissues, on the assumption that multi-scattering is neglected, the relationship between the transmission ultrasound speed and temperature can be expressed as a linear transition between these two tissues.
To refine the speed of sound image segmentation in the mixture area, homogeneous areas, and the corresponding mixture area are separated by the following rule:
![]() |
![]() |
![](cpb_26_6_064301/cpb_26_6_064301_eqn003.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn004.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn005.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn006.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn007.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn008.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn009.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn010.gif)
With the profiles between the speed of sound versus temperature in
![]() |
![](cpb_26_6_064301/cpb_26_6_064301_eqn014.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn015.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn016.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn017.gif)
The propagation speed of sound in the mixture area of two homogeneous tissues is a multi-variable function caused by factors like the applied intensity absorption, the distribution of metabolic heat generation, heat conduction in the tissues, heat convection by flowing blood, heat loss by radiation, convection, and evaporation at the surface, and heat exchange between veins and arteries. Although the mapping from speed of sound to temperature in the mixture area is a multi-variable function, we can still simplify it by a linear interpolation PSF since the thickness of the mixture area is relatively small compared with that of the homogeneous area. Then, we can obtain
![]() |
![](cpb_26_6_064301/cpb_26_6_064301_eqn019.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn020.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn021.gif)
![]() |
Substituting Eq. (
![]() |
![](cpb_26_6_064301/cpb_26_6_064301_eqn024.gif)
Since the EMD is characterized with an adaptive property for a non-stationary signal, the IMFs from the temperature EMD reflect the versatile intrinsic relationship between the ground truth values and artifacts. Thus, attaching a proper weight coefficient on each IMF and reconstructing it under different constraint can yield a signal similar to that applied to a certain filter.[35] In this contribution, after building the weight coefficient through ultrasonic transmission tomography, we adopt EMD on the temperature image without its DC component to generate IMF and then adaptively select the weight coefficient under continuity and similarity constraints for image optimization. The method is described in detail as follows.
The overshoot artifacts reflect the overshoot of temperature around the boundary, which is not correct in physics. Using the same method as in our previous study,[31] we take the mean variation of adjacent values from optimized temperature image as an evaluation criterion, which is defined as
![]() |
![](cpb_26_6_064301/cpb_26_6_064301_eqn046.gif)
![](cpb_26_6_064301/cpb_26_6_064301_eqn047.gif)
![]() |
Application of our proposed method to the abnormal data (marked by the red curve in Fig.
Reflection mode ultrasound image (Fig.
From the original un-optimized three-dimensional (3D) temperature (Fig.
![]() | Table 1.
Mean temperature variation in different angles. . |
When our previous method is applied, the temperature distribution of the mixture area is calculated by the finite element method with using Comsol® Multi-Physics 3.5a with 25 repeated iterations and temperature of the other area is also calculated. The total calculation time is over 5 minutes when using an 8-core Intel ® X5570 HP XW8600 workstation with 64-GB memory. In addition, for comparison, the temperature field of the whole area optimized with the current method is also computed using an nVidia ® Quadro K5200 GPU card since it contains no iteration. The total calculation time is less than 400 mil-seconds, which proves the advantages in reducing the computation burden with our proposed method in this paper. Furthermore, the whole optimization and reconstruction procedure is characterized with an adaptive property, which does not need any prior information or assumption.
Empirical mode decomposition is an excellent data analysis method for non-stationary signals with adaptive property. Signals reconstructed from IMFs with different weights can prevail favorably over those from a complicated filter. Here in this paper, our adopted EMD approach combined with similarity constraint and continuity constraint, can yield better optimization temperature images adaptively with less computation time than our previously proposed method. The phantom in this study mimics the cancerous and benign masses embedded in glandular tissue of human breast. In future human clinical study, temperature imaging through ultrasound transmission tomography for localized feedback of hyperthermia can present the more precise details to enhance the feasibility with acceptable computation time when stable ultrasound transmission tomography is applied.
[1] | |
[2] | |
[3] | |
[4] | |
[5] | |
[6] | |
[7] | |
[8] | |
[9] | |
[10] | |
[11] | |
[12] | |
[13] | |
[14] | |
[15] | |
[16] | |
[17] | |
[18] | |
[19] | |
[20] | |
[21] | |
[22] | |
[23] | |
[24] | |
[25] | |
[26] | |
[27] | |
[28] | |
[29] | |
[30] | |
[31] | |
[32] | |
[33] | |
[34] | |
[35] | |
[36] | |
[37] | |
[38] | |
[39] | |
[40] |