Thermometry and ablation monitoring with ultrasound

Abstract

In this review we present the current status of ultrasound thermometry and ablation monitoring, with emphasis on the diverse approaches published in the literature and with an eye on which methods are closest to clinical reality. It is hoped that this review will serve as a guide to the expansion of sonographic methods for treatment monitoring and thermometry since the last brief review in 2007.

• Ablation monitoring
• thermal therapies
• thermometry
• ultrasound

Introduction

Image-guided thermal therapies are now routinely applied in a variety of clinical settings. For laser, RF, microwave, and high-intensity focused ultrasound (HIFU) therapies, several medical imaging modalities have been proposed, tested, and employed in the guidance of these therapies in research, culminating now with MRI-guided HIFU (MR-HIFU) for patients in the clinic. Despite the positive clinical experience with MR-HIFU and the ongoing clinical trials for new systems and indications, image-guidance with ultrasound remains particularly attractive in simplicity, portability, accessibility, and cost. While many groups have explored the potential for both quantitative ultrasound thermometry and qualitative ablation monitoring with ultrasound, ultrasound-guided focused ultrasound surgery (USgFUS) has many outstanding questions and technical challenges before clinical use of the technology to monitor and control thermal therapies is routine. Here, we review the current status of ultrasound thermometry and ablation monitoring, with emphasis on the diverse approaches published in the literature and with an eye on which methods are closest to clinical reality. Although ultrasound thermometry would have applications to all thermal therapies, this paper refers to HIFU as the specific energy-based thermal therapy modality of interest.

A review of the English-language literature indexed in PubMed from 1967 to the present was completed for the following search expressions: (i) ultrasound thermometry, and (ii) (‘ultrasonic’ OR ‘ultrasound’) AND (‘temperature’ OR ‘thermal’) AND (‘estimation’ OR ‘measurement’ OR ‘imaging’). In some cases, these citations led to additional citations in Scopus and Web of Science, primarily in physics, acoustics, and signal processing journals not indexed by PubMed. In most cases, the related term ‘thermography’ refers to methods for measuring the temperature at the surface of the breast, scrotum, intravascular vulnerable plaque, or other sites using thermal or optical methods. While these methods may provide information about lesions at depth, the thermography literature (which is older than thermometry) is not reviewed here, except for the exceptional papers where the term ‘ultrasound thermography’ was noted. It is hoped that this review will serve as a guide to the expansion of sonographic methods for treatment monitoring and thermometry since the last brief review included in Rivens et al. 8 years ago [1].

Why measure temperature?

The measurement of temperature during thermal therapy provides an accurate estimate of the region of thermal damage that will occur in tissue. The relationship between time, temperature, and cell survival has been confirmed in numerous studies ranging from cells to excised tissue, animals, and patients. The influence of temperature on cells has been observed to follow an Arrhenius relationship from mild heating to coagulative temperatures [2]. Histological evaluation of tissues following ablative treatments supports the notion that temperature shows a strong correlation with regions of coagulative necrosis for short-duration exposures as typically encountered during HIFU [3]. As a result, the ability to image the spatial temperature distribution during treatment provides a form of ‘in vivo dosimetry’ that is unavailable with other energy-based therapies such as ionising radiation. In addition to predicting regions of thermal damage, in vivo temperature imaging can also be used to monitor sensitive areas to protect them from unwanted damage, thereby enhancing the safety of a treatment.

Over 35 years ago, Christensen [4] noted that temporal resolution requirements for thermometry were reduced by the large thermal mass of treated tissues. The temporal resolution required to effectively monitor a thermal therapy depends on the rate of heating and the time response of tissue governed by the bioheat equation. For interstitial treatments where the heating duration can last up to 15 min, a temporal sampling of 10–30 s may be sufficient. For HIFU treatments where a single exposure reaches coagulative temperatures in 10–30 s, a temporal resolution of a few seconds is necessary to adequately track and control treatment. Spatial resolution requirements are also tightly related to the type of heating. Interstitial treatments or hyperthermic exposures that gradually heat a volume many cubic centimeters in size typically have gentler thermal gradients and can be monitored effectively with a spatial resolution of 2–3 mm. HIFU treatments have very sharp thermal gradients and heat volumes much less than 1 cm³ in a single exposure; hence a spatial resolution of 1 mm is often desirable (but usually not attainable). In general external HIFU treatments place the strictest requirements on spatial and temporal resolution for thermometry. The final consideration is the temperature precision required for a thermometry method. In the case of ablative treatments, elevating temperature beyond a critical temperature is often sufficient to achieve irreversible tissue destruction. Beyond this threshold, the effects are largely the same until the temperature exceeds the boiling point in tissue. Therefore, a precision of 2–3 °C for thermometry is often acceptable for these types of treatments. On the other hand, for hyperthermic treatments, a difference in 1 °C results in a requirement for double the exposure time to achieve the same biological end point in tissue. Therefore, for these treatments a much higher precision, ideally less than 1 °C, is necessary for adequate monitoring and control. The reader will appreciate that spatial, temporal and temperature resolution are intrinsically linked and are traded off for one another to customise for a particular heating application. HIFU treatments typically try to achieve the highest spatial and temporal resolution with a trade-off on the temperature precision, which is acceptable for coagulative treatments. In hyperthermia treatments a lower spatial and temporal resolution can be tolerated, with a concomitant gain in temperature precision. A few studies have looked at these trade-offs through simulations and have evaluated their impact on treatment accuracy and precision [5–7].

Gold standards for in vivo thermometry

Invasive temperature sensors

The gold standard for in vivo temperature measurement is the utilisation of invasive temperature sensing probes [8]. Temperature sensors generally fall into two classes: thermocouples and fibre-optics. The ideal invasive sensor would not modulate the thermal therapy field and would be insensitive to strains due to thermal expansion and/or patient motion. Significant distortion of the HIFU field was reported for metallic thermocouple probes larger than , but probes smaller than produced only local perturbations that are thought to have negligible effect on the overall thermal distribution [9] but may create a discrepancy between the temperature at the probe and the surrounding tissue [10]. Fibre-optic probes are generally insensitive to environmental electromagnetic interference, but can be influenced by light in the case of laser thermal therapy [11].

The use of invasive temperature probes reduces the advantage of non-invasive thermal therapies such as HIFU. Further, accurate temperature measurements are only obtained at the spatial location of the sensor, precluding estimation of the spatial distribution of temperature in biological tissue due to heterogeneous blood flow and energy absorption. As a result, multiple temperature sensors are typically required for adequate monitoring of spatial heating, which is often impractical clinically because of the inconvenience or anatomy (i.e. brain, bone).

MR thermometry

The current clinical standard for non-invasive temperature measurement in the body is with magnetic resonance imaging (MRI). Generally speaking, the advantage of MRI for thermometry is that it provides quantitative temperature measurements in the body without the need for a priori calibration in the target tissue. Furthermore, it is non-invasive, non-ionising, and can be acquired in multiple planes or volumetrically. The spatial, temporal and temperature resolutions are acceptable for clinical thermal therapies such as hyperthermia and thermal ablation. The primary disadvantage is the cost and lack of portability associated with the method, as well as the need to design custom therapy systems to work within the strong magnetic field of the scanner. Excellent reviews can be found on this topic and its application in medicine [12–15]. Almost all the tissue parameters associated with MRI exhibit some form of temperature dependence including relaxation parameters (T1, T2), water diffusion, magnetisation, and proton resonant frequency. Exploiting these physical phenomena, MRI can be used to measure relative temperature changes or absolute temperature in vivo. The most mature and utilised technique for MR thermometry is the proton resonant frequency (PRF) shift method [16,17], which provides relative temperature measurements in soft tissues (excluding fat).

The PRF shift method is based on the physical phenomenon that in a magnetic field, water protons experience an approximately −0.01 ppm/°C shift in the Larmor frequency as the temperature changes. The origin of this frequency shift is attributed to the temperature dependence of hydrogen bonding in water, which affects the local microscopic magnetic field experienced by protons on water molecules. The change in magnetic field is reflected as a change in the resonant frequency of these protons when excited with a radiofrequency pulse in a static magnetic field. With the development of MRI, efforts were successful in implementing imaging strategies to exploit this phenomenon based on the phase information in MR images [16,18]. Since the change in phase due to temperature is much smaller than other changes in phase due to factors including susceptibility variations in tissue and static field inhomogeneities within the scanner, a subtraction of a reference image is usually performed to remove these static variations (leaving behind only changes in the magnetic field due to temperature). While this subtraction is generally robust, it renders the method a relative temperature technique, and the temperature of the reference image must be known or assumed. Early studies of this technique demonstrated that the PRF shift was approximately the same across different tissue types, and in native and coagulated tissues [19]. Further, it was sensitive to other factors including magnetic field drift, motion, and susceptibility variations from sources including lung volume changes, and moving devices in the vicinity. In addition it was observed that fat did not exhibit the same temperature-dependent frequency shift, and the technique could not measure temperature in this tissue. As a result, applying PRF shift methods in fat-rich organs such as breast has been challenging and typically requires alternative methods based on T1 or T2 relaxation [20,21]. Furthermore, PRF shift thermometry is challenging to apply in the abdomen due to the proximity of the lungs and subsequent motion of organs due to respiration. Techniques to overcome the challenges in performing simple subtraction-based thermometry include the referenceless method of thermometry and multi-baseline MR thermometry [22,23]. Nonetheless, this technique enables measurement of relative temperature changes in the body with a resolution of a few mm, temporal resolution of seconds, and a temperature precision of approximately 1 °C in human subjects on a clinical scanner. The temperature dependence of phase changes is independent of patient, organ, and whether tissue is treated or untreated, when the optimal conditions are met. Finally, numerous studies have been performed to validate the quantitative accuracy and precision of this method. It is for these reasons that the technique is used in conjunction with approved medical devices for the delivery of thermal therapy in the body using MRI-guided high-intensity focused ultrasound therapy [24,25].

Other non-ultrasonic approaches to thermometry

Despite the desire to minimise ionising radiation, investigations on CT-based thermometry date back to shortly after clinical introduction of the modality in the 1970s. From a physics perspective, the linear attenuation coefficient is a product of the mass attenuation coefficient of the tissue and the density. Since density changes slightly with temperature, the linear attenuation coefficient μ, and consequently the Hounsfield units for a heated pixel, will change compared to baseline. Unfortunately, the thermal sensitivity coefficient for CT thermometry is in the order of −0.50 HU per °C which limits the thermal resolution. For mild hyperthermia, this sensitivity is most likely limiting. However, the spatial resolution of underlying modality is quite high and offers the prospect for high-resolution thermometry with the cost of increasing the dose. An excellent historical and current status review of CT-based thermometry was covered in a recent review article on thermometry methods for laser-induced thermotherapy (LITT) [8].

Using temperature measurements on the surface (thermography), numerical inversion of the bioheat equation has also been attempted [26]. In this ambitious study the investigators also attempt the simultaneous reconstruction of not only the temperature increase map in the subsurface volume, but also the acoustic absorption coefficient and the thermal diffusivity. As inverse problems for models related to the modified Helmholtz equation are in general ill-posed, parameter estimation with the bioheat equation is extremely challenging: the expected resolution at depth cannot compete with image-based guidance. Since ill-posed inverse problems are commonly regularised using prior information from an ancillary imaging modality, it is doubtful that this approach will have advantages over ultrasound thermometry, even in the developmental stage.

Ultrasound considerations for HIFU monitoring

In this section we review the types of ultrasound acquisition geometries employed in ultrasound thermometry and ablation monitoring, the various types of raw data that are incorporated in the subsequent temperature estimation problem, and the potential forms of ultrasound contrast that can be exploited.

Ultrasound acquisition geometries

Passive

Passive acoustic mapping is a recent area of interest for monitoring HIFU. At present, it appears to be only applicable to detecting non-thermal signals from processes such as cavitation. However, since the passive acoustic signals from cavitation and thermal strain are also affected by temperature-dependent changes in speed of sound (SOS) and attenuation, it may be possible to estimate relative temperature changes using these signals.

In passive imaging the instrumentation for data acquisition is relaxed since transmit control is not required. In some passive applications a single, focused ultrasound transducer is used as a sensor to receive acoustic emissions generated in the HIFU treatment zone. Confocal HIFU and diagnostic transducer probes work in this regard. With larger array-based apertures and receiver-side beamforming, larger volumes may be monitored, and in practice diagnostic transducers may be used in this application, but they may not be optimally matched in terms of bandwidth to the acoustic emissions. In some case the same transducer is used for both passive and active cavitation detection [27]. While circular apertures associated with tomographic imaging could be utilised for passive mapping, they have not been fully explored due to the limited number of organs outside the breast and brain that can be investigated with ultrasound tomography.

Pulse–echo

Pulse–echo is the most common ultrasound geometry due to its routine use in diagnostic ultrasound instrumentation. Exploiting the principle of reciprocity, an ultrasound transducer is used to both transmit an acoustic wave into tissue and receive waves that arise through interaction processes inside the tissue. If the interaction of the ultrasound wave with tissue is temperature-dependent, then the signals received as echoes will exhibit some modulation from which temperature could ideally be estimated. As with the passive geometry, single element and array-based transducers are both commonplace. Depending upon the nature of the temperature estimation method, entire portions of the normal diagnostic ultrasound image formation chain can be utilised as is.

Pitch–catch

If the transmitter and receiver transducers are distinctly off-set from each other, then a pitch–catch geometry is invoked. Unlike pulse–echo, pitch–catch geometries facilitate transmission of acoustic waves through tissue sections if the transducers are facing each other. The information content for transmission signals is in many ways complementary to reflected signals in pulse–echo. For example, the transmitted acoustic wave is more directly influenced by SOS (as seen in propagation time) than the subtle effects of SOS on pulse–echo data. While pulse–echo is sensitive to changes in attenuation due to scattering, the transmitted wave is directly modulated by attenuation due to both scattering and absorption. The most complete pitch–catch geometry is of course the full aperture tomographic geometry, where transmitted waves are captured along with waves scattered at all angles.

Hybrid

Hybrid geometries combine ultrasonics with another imaging modality to produce an imaging system with increased performance. The hybrid imaging geometry of most interest in thermometry is photoacoustics: high spatial resolution but low contrast ultrasonics is combined with techniques from optical imaging where the resolution is poor but the contrast due to endogenous chromophores is good. In the photoacoustic geometry, a brief but intense pulse of light is produced by a nanosecond laser capable of outputting more than 10 mJ per pulse in the order of 10 times per second. Due to differential optical absorption, small temperature increases in tissue lead to thermal expansion and the subsequent generation of an acoustic wave. Similar to passive imaging above, ultrasound receivers can be synchronised to the illuminating light, and the photoacoustic waves can be detected. Reconstruction of the sources of the photoacoustic waves can be accomplished either by raster-scanning a focused transducer or through receiver-side beamforming with arrays.

Data types

The numerous ultrasound temperature estimation algorithms found in the literature span the entire image formation chain. At the lowest level, a temperature estimation algorithm can work with the raw RF data from the ultrasound transducer or array (the so-called channel or pre-beamformed RF data). Data at the level of the individual transducer is only available on select ultrasound research platforms (e.g. Verasonics). More commonly, the post-beamformed RF data is available, and in some cases this has been demodulated at the centre frequency and filtered into in-phase/quadrature-phase (I/Q) data. I/Q data are also available on some systems for pre-beamformed channel data. In clinical ultrasound systems, data up to this point in the image formation chain are available only in a few select systems from a minority of vendors. At the highest level, temperature estimation algorithms can take as input the processed images produced by the scanner, including B-mode, Doppler, and in the latest generation of clinical devices elastography images. The speed at which these images can be derived from a clinical device is highly variable and an important consideration for research purposes. It should also be noted that recent ultrasound research platforms also allow significant transmitter-side control, including arbitrary waveform generation and the integration of HIFU interleaved with imaging in custom probes.

Contrast

The predominant source of ultrasound contrast is variation in acoustic impedance, the product of density and ultrasound velocity. In the vast majority of ultrasound applications, variations in density are essentially ignored. Even in ultrasound tomography, modelling and image formation using the constant density wave equation is much more common than the variable density wave equation. While direct detection of temperature-dependent changes in density are not possible with ultrasound, the underlying thermal expansion or contraction can be detected using strain-based techniques, as described below.

The SOS in mammalian tissues is also temperature dependent [28]. For fat, the ultrasonic speed decreases monotonically with increasing temperature. For most tissues, however, the SOS increases, levels off near 60 °C, and then decreases at higher temperatures. It is unfortunate that the regime where the SOS sensitivity to temperature is lowest is also in the range where tissue destruction occurs in ablative treatments. So while methods that estimate SOS may be applicable to hyperthermia temperatures below 50 °C, temperature estimation in the ablative zone is challenging.

While the attenuation of water is known to decrease monotonically with temperature [28,29], attenuation measurements in various tissues produce mixed results. Complicating the attenuation measurement is the fact that attenuation is also frequency dependent. In bovine tissues attenuation decreases with temperature until a point above 50 °C where the attenuation begins to rise. For canine muscle, liver, and kidney, ultrasound attenuation was reported to double between 50 °C and 65 °C [30]. Bamber and Hill [28] did not observe this attenuation rise in human liver, but Techavipoo reported similar sigmoidal behaviour in canine liver [31,32]. Qualitative 3D maps of relative ultrasound attenuation changes have also been reported as indicative of tissue coagulation for other thermal therapies [33–35].

Since ultrasound is a coherent imaging modality, it is sensitive to multiple scattering centres in a single resolution cell. The subsequent speckle has proven to be a useful tool for estimating thermal strain and for motion estimation that must be dealt with in many in vivo applications. For that reason, image formation chains that utilise techniques to reduce speckle (e.g. spatial and frequency compounding tomography) may not be optimal for temperature estimation methods that rely on speckle tracking.

Bubbles in the thermal treatment zone have one of three sources: cavitation (both inertial and non-inertial), gross tissue boiling, or an exogenous source in the form of a contrast agent. Since boiling is a clear indication that certain temperature thresholds have been reached, efforts have succeeded in monitoring the onset of boiling using passive acoustic detectors. Cavitation bubbles occur before boiling and are a function of ultrasound power, not temperature. In many cases, cavitation is a nuisance to be minimised, since it distorts most ultrasonic ablation monitoring techniques. Some groups have reported both active and passive detection of cavitation bubbles for therapeutic monitoring. Doppler imaging modes, which are sensitive to both flow and exogenous contrast agents such as microbubbles, have recently been used in HIFU studies. Li et al. [36] more recently described the combination of colour, pulse inversion, and decorrelation Doppler to detect pulsed HIFU-induced cavitation bubbles with short transient times. Colour Doppler has also be used to detect hollow silica nanoshells delivered to a tumour to enhance the HIFU therapy [37].

At temperatures below boiling, other phase transitions can occur in tissue due to chemical and structural modification of proteins and other cellular components. These changes are generally classified as coagulation, but surprisingly these changes are not reliably monitored using B-mode ultrasound. Coagulation also leads to destruction of tissue blood flow and perfusion. Using B-flow imaging, a non-Doppler mode for imaging blood flow, Andrioli and Valcavi [38] reported monitoring of RF ablation in the thyroid. This method has not been tested with HIFU ablation though.

The ultrasound non-linearity parameter B/A has also been reported to exhibit temperature dependence. Using a technique that appears similar to but antecedes vibroacoustography, Ueno et al. [39] reported a hyperthermia thermometry system targeting this parameter with validation studies in pig thigh. Van Dongen and Verwij [40] also reported using the non-linearity parameter to estimate temperature in a glycerol-based phantom using the Burgers equation as the model.

To test ultrasound thermometry and ablation monitoring methods, appropriate phantoms are required [41]. In many cases phantoms and standards adopted from the HIFU literature are utilised [42]. Recently, a dual modality (ultrasound and CT) phantom specifically for non-invasive ultrasound thermometry calibration was reported using agarose mixed with variable amounts of fat-mimicking glyceryl trioleate [43]. Gelatin-based phantoms that match both the acoustic and thermal properties of tissue are also possible [44].

Specific strategies for ultrasound thermometry and ablation monitoring

The ideal ultrasound ablation monitoring system would provide accurate, precise and quick measurements of tissue temperature in the range from 37 °C to boiling. Unlike MR thermometry, however, nature is not kind and conspires against the use of ultrasound thermometry above 50 °C through a variety of physical and acoustic complexities, including non-linear temperature dependencies of contrast parameters, tissue phase transitions, and the stochastic generation of cavitation bubbles. In general, ultrasound thermometry is only applicable for tissue hyperthermia. For HIFU tissue ablation, non-thermometric ultrasound methods are required for therapy monitoring.

B-mode imaging

In general, the diagnostic ultrasound imaging cannot take place while the HIFU treatment is active, but rather must be interleaved. Even with a 10–15 ms delay after cessation of HIFU, a 31.25 ms interference band artefact was reported on B-mode imaging [45]. A recent report, however, described a potential method for simultaneous HIFU therapy at 4 MHz with ultrasound imaging in the 3–9 MHz bandwidth using pulse compression and the application of second order infinite impulse response notch filters on the received data to eliminate cross-talk from the HIFU [46]. Pulse inversion has also been described as a method to facilitate simultaneous US therapy and imaging [47].

The correlation between increased echogenicity and temperature in tissue is slight. In one study, a correlation coefficient of r = 0.455 was reported [48]. It is assumed by many groups that the relationship between B-mode echogenicity and temperature is complicated by the fact that the structures that result in increased echogenicity (and in some cases post-lesion shadowing) also severely modulate the HIFU beam and the subsequent temperature distribution, potentially leading to over-treatment. While B-mode ultrasound is available on clinical HIFU systems popular in China, the feedback from the imaging modality is often qualitative [49,50].

Alvarenga, Teixeira and collaborators [51–53] have studied the possibility of using the Gray-level co-occurrence matrix (GLCM), the GLCM-derived entropy, and the average grey-level content (AVGL) to estimate temperature in vivo. Since there is excellent correlation between temperature and these B-mode image-derived indices up to 44 °C, these easy to implement methods could be utilised for quantitative monitoring of hyperthermia. Unfortunately proof-of-principle studies have been limited to only phantoms and tissue samples.

As described below in thermal strain imaging, a subset of systems and groups use the B-mode image for speckle displacement estimate, but the trend is clearly toward using lower level RF data from the scanner.

Nakagami imaging

Ultrasound Nakagami imaging is a recently developed mode of imaging that spatially maps the shape parameter for an assumed probability density function for the backscatter envelope [54–56]. The Nagakami distribution is defined by two parameters: m and a scaling parameter Ω: where Γ is the gamma function, and U is the Heaviside step function. The Nakagami parameter m can be directly estimated from statistical averages of the envelope R of the backscattered ultrasound signal (available either indirectly from RF data or directly in B-mode images)

The scaling factor is simply the time-average intensity:

During HIFU, both cavitation and boiling may produce bubbles that produce acoustic shadowing. Zhang et al. [57] at Xi’an Jiaotong University (China) explored the potential of using Nakagami imaging to overcome the artefacts in traditional B-mode images. By inserting a digitiser at the appropriate point in a modified diagnostic ultrasound scanner (DP9990, Mindray, Shenzhen, China), they captured all post-beamformed RF lines and estimated the Nakagami parameter along each A-line using a 1D sliding window. 2D Nakagami images were then derived from the collection of A-lines. They found that the Nakagami parameter was statistically higher in the zone of thermal ablation compared to any regions where bubbles formed. However, the lesion size predicted by the Nakagami image underestimated the true lesion size as confirmed in polyacrylamide+BSA phantoms. More recently, the same group has proposed a lesion-to-bubble ratio (LBR) derived from Nakagami (and other imaging modes) to monitor thermal ablation [58]. They reported that a ratio of 2.0 on Nakagami images was indicative of protein coagulative necrosis despite shadowing of the therapy lesion by bubbles. However, the authors acknowledge that this approach may only be valid for small thermal lesions.

The Tsui group [59] at Chang Gung University (Taiwan) have also explored Nakagami imaging for monitoring both RF and ultrasound ablation. Compared to B-mode imaging, the Nakagami image exhibited better contrast-to-noise ratio in the lesion vs. background. Although using microwave ablation for their investigation, the Nakagami images were produced from a Terason commercial ultrasound system (Burlington, MA, USA) using a Nakagami window three times the width of the transmitted pulse. In a related study of the temperature-dependence of the Nakagami parameter and derived parameters such as absolute regional ratio of change in Nakagami parameter, the same groups reports that ‘the magnitude of change in the statistical distribution of the backscattered-signal envelope is a monotonic function of temperature over a temperature range from 37 °C to approximately 45 °C [60]. Since these methods involve a polynomial fit to the derived data, this latter work includes a discussion of the potential differences in appropriate polynomial for HIFU vs. RF and microwave ablation.

In a recent report, investigators in Toronto, Canada and Tehran, Iran reported similar experiences with Nakagami imaging using a confocal arrangement of a HIFU transducer with an endocavity convex array probe applied to ex vivo porcine muscle [61].

Speckle tracking and thermal strain imaging

In an imaging scenario with no motion, the apparent displacements of the ultrasound speckle due to localised heating is correlated with both changes in the SOS and thermal expansion of the tissue. Local displacements can be estimated from the time shifts Δt in the RF data, and given a known SOS c₀ at baseline temperature T₀, the temperature change can be related to the rate of change of time shifts with depth z: where α(°C⁻¹) is the thermal expansion coefficient and the SOS sensitivity to temperature is known (). This linear relationship is most often assumed to be valid for hyperthermia below 50 °C. In the range above 50 °C up to phase transitions such as coagulation, the relationship is non-linear and may be difficult to calibrate. A more complete derivation of the basic thermal strain theory can be found in a recent review [62].

Since both the thermal expansion and the change in ultrasound velocity lead to shifts, it is important to understand the relative contribution of these two competing processes. In general, the change in SOS dominates thermal expansion by an order of magnitude [63–65]. Souchon and colleagues modelled the contributions of thermal expansion and SOS change in a one-dimensional model [66]. As shown in Figure 1, the largest echo strain is shifted to lower temperatures compared to the temperature of maximum speed of sound. Thus, there is the potential for quantitative temperature measurements above 60 °C provided that the speckle pattern is not destroyed by coagulation or other destructive phenomena.

Figure 1. 1D simulation of relative contributions of speed-of-sound change and thermal expansion to observed strain. Reprinted, with permission from Elsevier [66].

In pioneering work, Seip and Ebbini [67] modelled the RF A-line signal as an autoregressive (AR) model of discrete tissue scatterers that are shifted by thermal expansion and pseudo-shifted by changes in SOS. Through careful analysis, a relationship between the temperature change and changes in the power spectral density (PSD) for this model were developed and tested in phantoms. Of particular note, the authors commented that standard data analysis techniques using the fast Fourier transform would not be able to detect temperature changes in the order of 1 °C. By utilising the autoregressive PSD approach, however, arbitrary sub-degree temperature resolution is enabled, subject to noise considerations.

In practice though, cross-correlation techniques quickly supplanted the original PSD approach [64,68]. To estimate displacements in either RF data or B-mode images, the most common algorithm today is the maximum normalised cross-correlation: where S₀ and S₁ are the baseline and subsequent ultrasound data, with the bar quantities denoting spatial averages. A straight two-dimensional cross-correlation is known to be sensitive to changes in intensity. The kernel or window area for the cross-correlation must also be considered. We note in passing that these methods can also be applied to regions-of-interest (ROIs) outside the HIFU intent-to-treat zone. In those cases, displacements are related to gross motion of the animal or patient, and this information can be used for motion correction or compensation with other temperature estimation methods (see the backscatter section below). For real time temperature imaging, speckle tracking algorithms have been implemented on graphical processing units (GPUs) [69–72]. The 500 Hz temporal resolution reported in Liu and Ebbini [69,71] is probably the best reported in the literature. This extreme temporal resolution is not useful from a heat transfer perspective, but it may be advantageous for minimising artefacts due to respiratory and cardiac cycles in limited fields of view and for rapidly detecting bubble formation.

The other important aspect of this model is that the temperature change is related to the gradient of the time shifts as a function of depth. This gradient can be numerically calculated by a variety of methods. In Lai et al. [73] the Sobel operator is recommended since it has a smoothing effect which offsets the increase in noise associated with numerical differentiation. An alternative approach also adopted from elastography is to estimate the displacement using a local parametric model, avoiding differentiation. Ye et al. [74] reported the implementation of an iterative Kalman filter to estimate a model-based temperature map from B-mode images. A typical processing pipeline for thermal strain imaging as described in [73] is reproduced in Figure 2.

Figure 2. Prototypical thermal strain image formation chain using displacement (disp.) estimation. © 2010 IEEE. Reprinted, with permission, from Lai et al. [73].

Adolhassani and colleagues reported that isotropic speckle-tracking kernels from 2 mm to 4 mm had similar performance when maximum normalisation cross-correlation was applied to B-mode images [75]. In later work, Mehrabani et al. [76] utilised the Horn-Schunck method for visualising optical flow to estimate the speckle displacement. Their B-mode ultrasound thermometry technique was ultimately tested in a sheep model where occlusion of the renal segmental artery led to a decreased renal lobular temperature [77].

When strain imaging was used for qualitative ablation monitoring during HIFU, the contrast for strain imaging was superior to B-mode imaging. The Chapelon group [78] reported a mean contrast ratio between treated and untreated ROIs in porcine liver to be 22.0 ± 10.3 dB for strain images compared to 9.3 ± 5.5 dB. The B-mode image analysis was complicated by the fact that the treated zone exhibited both hyper- and hypoechoic subregions. When analysed separately, neither subzone contrast was better than strain imaging. Ultrasound in this work was, however, applied at rather high frequency (12 MHz). Miller et al. [79] analysed the minimum signal-to-noise ratio (SNR) in strain images that would be required to visualise a thermal lesion of a known size: in this case, an ellipsoid 10 mm long and 1 mm wide. The authors estimated that strain images with a SNR greater than 20 dB exceed the threshold for detection. In the case that the thermal lesion size is unknown, they reported a lower bound on the SNR threshold for size estimation of 26 dB. Miller and colleagues also analysed echo strain imaging for coaxial and perpendicular geometries relative to the HIFU transducer [80]. In an atypical coaxial geometry, the imaging transducer was opposite the therapy transducer. Although the coaxial geometry is more appropriate for clinical utility, the perpendicular geometry actually exhibited a higher contrast-to-noise ratio (2.00 ± 0.72 vs. 0.37 ± 0.24). In more recent work, this group attempted to calibrate the relationship between strain and temperature deviations up to 8 °C using a second order polynomial model [81].

Thermal strain imaging has been one of the most prolifically studied methods for estimating temperature with ultrasound [82–91]. Some continued work on temperature estimation using frequency shifts in the spectral domain due to thermal strain has been reported [92,93]. Frequency domain estimation of pulse–echo time shifts has also been reported using a zero-crossing detector [94–96]. Recently, axial-shear strain elastograms have been proposed as an improvement to traditional axial-strain elastograms in thermal strain imaging, especially when cavitation bubbles are present [97]. Using strong reflectors in tissue (e.g. implanted thermocouples), Teixeira and colleagues in Portugal utilised neural networks to estimate temperature from echo-shift RF lines [98–100]. In their initial work, the neural network input consisted of five data channels for five scatterers. It is not clear if this approach would scale to systems for speckle tracking.

In a separate approach, ultrasound-estimated strain has also been used to estimate the thermal diffusivity K (m²/s) in tissue. Anand and colleagues showed that the rate of change of a secondary parameter derived from the estimated strain and displacement profile widths is directly proportional to K [101]. They used this non-invasive procedure to estimate K for the bioheat equation inverse approach described in the Passive section below. The authors proposed that a spatial map of K over a larger region of interest containing the target zone would be feasible.

The achievable spatial resolution of thermal strain imaging is limited by the same physical considerations as in ultrasound elastography [102]. The size of the spatial kernel utilised in speckle tracking directly influences both the SNR and the spatial resolution of the strain image, but with opposite effects. In practice, Shah et al. reported using an anisotropic kernel (0.6 mm axial vs. 2.1 mm lateral) as a compromise between these competing outcomes [103], but the optimal kernel may be tissue and depth dependent. In their phantoms, the lateral spatial temperature resolution was reported to be in the order of several millimetres, while in the axial direction the resolution was sub-millimetre. For further investigation of the effects of SNR on echo strain estimation for thermometry, the reader is also directed to early work in Miller et al. [79].

Tissue stiffness monitoring and elastography

The mechanical properties of soft tissue are temperature dependent, and both bulk and shear elastic properties exhibit changes with temperature [104–106]. The temperature dependence is non-linear, with a large transition occurring when tissues undergo thermal coagulation [107]. A number of groups have investigated the possibility of exploiting this physical property of tissue to monitor or control a thermal ablation therapy with both ultrasound and MRI. Most efforts have focused on utilising this effect for mapping or visualising the region of thermal coagulation in tissue after treatment is completed due to the large change in mechanical properties that occur after reaching this end point in tissue. Indeed, before the development of thermal strain imaging, classical ultrasound elastography using strain estimation was reported to be sensitive to changes in strain contrast due to both acoustic and optical thermal delivery [108–110]. However, some methods have attempted to track the changes in mechanical properties dynamically during heating to track temperature of stiffness changes as they occur.

There are obvious links between tissue stiffness monitoring and the previously described thermal strain imaging methods. In some cases, similar data analysis methods are used for both. For example, Sumi et al. used estimate displacement and strain in paired RF data frames to solve the relative change in shear modulus after coagulation in fresh calf liver [111]. Quasi-static elastography with a compression plate can also be used to monitor lesion stiffness after ablation [112], but both the pre–post contrast-to-noise ratio (CNRe) and tissue modulus are non-linear with temperature.

Shear wave thermometry

Imaging of shear wave propagation in tissues is possible with both MRI and ultrasound [113,114]. In the case of ultrasound this can be accomplished by generating tissue displacements remotely in tissue through the radiation force principle or low frequency mechanical stimulation. The displacement of tissues results in shear waves propagating from the location of the displacement, at a speed of approximately 2–8 m/s. Ultra-fast ultrasound imaging can be used to track the propagation of these shear waves, which is related to the underlying shear elastic properties of the tissue. One method developed for this type of imaging is referred to as supersonic shear wave imaging [115] in which a series of localised focal regions are generated in a linear fashion with depth in tissue, generating a cylindrical shear wave around these points. This serves to cover an entire imaging field of view with shear waves, enabling mapping of the elastic properties throughout a region of interest. One of the first applications investigated with this technique was the visualisation of thermally induced lesions in tissue due to the inherent change in elastic properties upon undergoing thermal coagulation [116]. A factor of three-fold increase in shear modulus was measured in coagulated chicken breast using this technique in this preliminary study. Subsequent to this work, the same group demonstrated the ability to visualise thermal lesions after RF ablation in vivo in a pig model using shear wave imaging [117]. In this study, the mean elasticity of lesions was 38.1 ± 2.5 kPa vs. 6.4 ± 0.3 kPa in normal liver. The authors concluded that by using a threshold of 20 kPa, the presence of coagulation could be detected with a sensitivity of 0.8 and a positive predictive value of 0.83 in vivo, which represents some of the best data to date reported for ultrasound imaging of coagulated tissue. In addition to visualisation of thermal coagulation in tissues, Arnal et al. [118] investigated the temperature dependence of the shear modulus using shear wave imaging, to evaluate the potential to use this imaging method for thermometry. The changes in shear modulus with temperature were compared to the temperature also measured with ultrasound imaging (using the method of Seip [67]), and a strong correlation was observed for heating below 40 °C. These measurements in ex vivo tissue samples indicated a change of approximately 1 kPa/°C in muscle and the authors suggested this may have application for monitoring temperature changes below 45 °C. They also claimed that this method was very robust to motion, which overcomes one of the limitations of traditional echo-shift approaches for ultrasound thermometry. Similar trends in the shear modulus with temperature and coagulation have been observed in other studies using 1D transient elastography [119,120]. Sapin-de Brosses et al. performed a series of studies investigating the relationship between shear modulus as measured with shear wave imaging and thermal dose in excised tissues and in vivo [104,105]. As highlighted in Figure 3, these studies identified a strong dependence of shear modulus with temperature and thermal dose, with different trends observed across tissue types and in the case of muscle depending on the orientation of the measurement with respect to fibre orientation. The in vivo studies, performed in rat muscle, demonstrated a consistent increase in shear modulus of at least five times upon reaching a thermal dose of approximately 150 cumulative equivalent minutes at 43 °C (CEM43). When the time required to cause this stiffness increase was plotted against the temperature of the muscle, an exponential relationship that overlapped with previous measurements of the same relationships for cells was observed. Through this careful analysis, the authors were able to conclude that the changes in stiffness observed with shear wave imaging were indicative of cellular necrosis. This work represents some of the more thorough evaluations of the use of ultrasound imaging for monitoring temperature and thermal damage.

Figure 3. Variation of shear wave velocity as a function of temperature, as assessed using shear wave elastography. © 2011 Institute of Physics and Engineering in Medicine. Reproduced by permission of IOP Publishing from Sapin-de Brosses et al. [105]. All rights reserved.

In shear elasticity estimation using speckle tracking, the axial spatial resolution is dependent upon the correlation time window widths and delays [121]. Using a total data acquisition time of 23 ms, a pulse repetition frequency of 3500 Hz, and 2.4 μs correlation windows with 60% overlaps, Benech and Negreira reported a 0.7 mm axial spatial resolution in a 1D transient elastography geometry [119].

Local harmonic motion (LHM) imaging

This method, also referred to as harmonic motion imaging (HMI), is more of a coagulation monitoring technique using ultrasound. The principle is based on the idea that an amplitude modulated acoustic wave can be used to simultaneously heat and perturb tissue in a harmonic fashion. The amplitude modulation is typically less than 100 Hz in most studies and produces motion in the 10–50 μm range. The harmonic motion generated in tissues is dependent on the local absorption of ultrasound energy as well as the underlying mechanical properties of the tissue. The technique monitors the displacement of the tissue element undergoing harmonic oscillation and relies on the premise that upon coagulation the changes in absorption and mechanical properties produce a measurable change (typically a reduction) in this displacement. This change in displacement can be used to indicate when coagulation occurs and can be used to terminate a sonication in the absence of a direct temperature measurement. The method was first described by Konofagou and Hynynen in 2003 [122] as an alternative to ultrasound-stimulated acoustic emission imaging, and was later evaluated for its diagnostic applications. In parallel, efforts were also made to evaluate the application of this technique for monitoring and controlling focused ultrasound surgery. A series of studies were performed in muscle and implanted tumours in rabbits comparing the signals obtained using LHM during HIFU with temperature measurements acquired using MR thermometry to validate the technique against an accepted standard [123,124]. These studies demonstrated a repeatable change in LHM amplitude that was frequency dependent and typically higher at lower frequencies (<100 Hz). They also indicated that the difference in LHM between native and coagulated tissue was significantly different and occurred at temperatures and thermal doses associated with thermal coagulation. In a study in implanted tumours it was observed that control of HIFU exposures using LHM resulted in successful coagulation in 69% of the exposures delivered [123,125]. Similar results were observed in excised tissues, gel phantoms and in vivo by Maleke and Konofagou [126], with these authors presenting M-mode imaging of the HIFU-exposed region as a function of time to estimate the extent of thermal coagulation depth in tissue [127]. These results were obtained using a confocal imaging and therapy system enabling B-mode monitoring of LHM during HIFU [128]. The results in vivo from this group indicated an LHM amplitude before and after coagulation of a mammary tumour from 27.34 ± 1.34 to 20.98 ± 1.82 μm across 16 lesions. Using GPU processing, real-time harmonic motion imaging has been reported [129]. A review of results and progress with this technique was recently published by Konofagou et al. [130].

Stimulated acoustic emissions

The method of ultrasound-stimulated acoustic emissions (USAE), or vibroacoustography, was introduced in 1998 [131], and is a technique for probing the mechanical properties of a medium through the use of two incident ultrasound waves with a slight frequency difference. The technique is based on the idea that the frequency difference between the beams (typically in the low kHz range) results in an oscillatory mechanical stimulation of tissue. This results in an acoustic emission that can be measured with a hydrophone or acoustic sensor, and is reflective of the underlying mechanical properties of the tissue. Building on earlier simulation work, Konofagou investigated the temperature dependence of this signal in excised muscle and fat tissue samples and observed a correlation between the USAE signal and temperature measured by a thermocouple at the focus of an ultrasound beam in both tissues [132,133]. The change in the USAE signal was measured to be approximately 0.8–1.0%/°C, and was most linear between the range of room temperature to approximately 45 °C. Above the coagulation threshold, non-linear responses were observed, and the response was different between fat and muscle. One challenge with this technique is that the signal depends both on the change in absorption and stiffness of tissue with temperature, and decoupling these two effects from the signal is impossible. Later investigations suggested that a strong frequency response exists in the temperature-dependent USAE response, and that the shift in frequency may be a more robust method for temperature measurement in tissue using this technique. Since this body of work, the technique has not progressed, and no results in vivo have been presented.

Backscattered energy and quantitative ultrasound

Scattering cross-section depends on differences in the compressibility and density of a sub-wavelength scatterer with respect to the surrounding medium. The compressibility of the two materials changes with temperature, based on their respective temperature-dependent changes in SOS, which results in a change in the scattering coefficient. Changes in scattering properties can be quantified by measuring the change in backscattered energy (CBE) received from a heated tissue volume, relative to that received at resting temperature [134].

Simulations from Straube and colleagues suggest that changes in backscattered energy will be monotonic with temperature for individual scatterers [134] and in images of random distributions of scatterers [135]. Experimental results in ex vivo tissue samples confirmed that CBE is monotonically positive (increasing power) with increasing temperature for lipid-based scatterers, and negative (decreasing power) for aqueous scatterers. For collections of heterogeneous individually tracked scatterers in 1D [136], or for image regions in 2D [137] or 3D [138], the overall standard deviation of the CBE was monotonically positive in the range of 37–50 °C. At 7.5 MHz, the calibrated temperature sensitivity of CBE was 0.300 ± 0.016 dB/°C (linear correlation coefficient 0.991–0.999) in ex vivo turkey breast, which allowed temperature measurement with 0.5 °C accuracy and 1 cm³ resolution in other samples of the same tissue type [138].

The accuracy of this technique depends on tissue composition, region size, SNR, and apparent motion. Because the technique aims to compare backscattered energy for individual scattering regions, the RF data must be compensated for changes in apparent scatterer location caused by changes in SOS or motion. In 1D, echoes from individual scattering regions were manually identified in the RF waveform and compensated [136]. For 2D and 3D CBE quantification with a linear phased array diagnostic probe, automatic echo-tracking was accomplished by maximising cross-correlation of the RF signal between consecutive images using either rigid [137] or non-rigid [138] registration. During homogeneous water bath heating and heterogeneous HIFU hyperthermia, 2D non-rigid registration identified and corrected between −0.1 and 0.7 mm of displacement in the axial direction where SOS changes are expected, and −1.1 to 1.1 mm (mean 0.4 mm) in the lateral direction which are primarily due to motion [139]. In three mouse tumours heated from 37.5–43 °C with HIFU, after motion compensation the means of image pixels with positive CBE monotonically increased by a mean of 3.5 dB (3.1 to 3.6 dB); the negative CBE monotonically decreased by −3.4 dB (−2.0 to −4.5 dB). Li et al. reported similar findings for small transducer displacements in both elevation and axial directions [140]. Breaking from the motion- and echo shifted-compensation CBE paradigm though, Tsui and colleagues proposed in 2012 that CBE imaging without echo shifting may have improved sensitivity to increased temperature detection [141].

Seo and colleagues reported using an L7-4 transducer with the Verasonics data acquisition system to monitor HIFU treatments in tumour bearing mice [142]. While they utilised CBE for temperature estimation, speckle tracking was also employed for a motion-mapped reference model. Motion compensation reduced the RMSE for temperature estimation from 4.2 °C to 1.1 °C.

In a related method, techniques from quantitative ultrasound for tissue characterisation have been used to assess the effects of different thermal doses. The Oelze group [143] measured both attenuation and backscatter coefficient (BSC). Using a parametric model of the backscattered signal, the effective scatter diameter (ESD) and effective acoustic concentration (EAC) were also estimated. In the hyperthermia range below 50 °C, attenuation changes up to 10–20%, BSC changes of an order of magnitude, and 30% changes in ESD and EAC were all reported as more sensitive than thermally induced changes in SOS. In their most recent work [144], the relationship between ESD, EAC, and temperature was extended up to 80 °C in experiments with in situ rat tumours. Good linear regression with temperature was observed in three out of four rats. The interested reader should also consider exploring the existing literature on the use of quantitative ultrasound for detecting other forms of cell death in vivo [145].

In another non-thermal strain technique that uses post-beamformed I/Q data, Mast and colleagues at the University of Cincinnati have explored echo decorrelation imaging, which is dependent upon the formation of an integrated backscatter image [146]. Echo decorrelation imaging is reported to be more robust to tissue motion and bubble activity than thermal strain imaging, but at higher temperatures the images are qualitative, not quantitative. Initial application of this technique was for RF ablation, but the technique should be compatible with HIFU [147].

In terms of temporal resolution, proponents of CBE report temperature maps upgraded every 1 to 2 s, which is competitive with MR thermometry [138]. While the temporal resolution is limited by the specific imaging hardware that implements CBE and the volume of interest, the authors developed a five-image 3D image acquisition that ran at 10 Hz. The ultimate temporal resolution was reduced by subsequent 3D motion compensation for the central slice and CBE processing.

Phase-contrast and thermo-acoustic lensing

Changes in the SOS due to temperature modulation can also be viewed as another form of phase aberration. In a transmission pitch–catch geometry, the received signal is sensitive to these slight aberrations, and this idea for thermometry was explored as early as 1977 by Sachs and Janney in fixed human brains [148]. Of particular note in this early work, the authors reported a spatial resolution better than 3 mm. Extending this approach, Clement and Hynynen developed a frequency-domain phase contrast filter to directly image the aberrations created by an orthogonal HIFU transducer [149,150]. To test this method they simulated a continuous wave signal transmitted through a treated tissue and then detected using an array with time-averaged field intensity sensitivity. Since this algorithm is sensitive to attenuation, frequencies below 1 MHz were simulated. At a localised heating event at the interface of fat and soft tissue, they report the successful mapping of a 5 °C temperature anomaly in both tissue compartments. The spatial resolution for ultrasound phase contrast thermal imaging can potentially be superior to MR thermometry. Farny and Clement recently reported 0.3 mm isotropic spatial resolution with a 1.68-MHz centre frequency [150]. Phase contrast in transmission ultrasound has also been reported to have sufficient sensitivity to detect sub-100 mK temperature increases associated with radiation therapy [151].

The phase aberration due to this thermal lens effect was also studied by Le Floch et al., since it also produces self-defocusing in the therapy beam [152] and is known to corrupt other temperature estimation methods. For pulse–echo temperature estimation methods, spatial compounding has been proposed and investigated as a method to reduce the effect of phase aberrations due to thermo-acoustic lensing in ultrasound temperature estimation [153]. With seven steered plane waves utilised in the spatial compounding, the variance in the temperature estimation was reduced to half.

Ultrasound tomography

Given the explicit relationship between SOS and temperature in water and soft tissue, it is perhaps not surprising that thermometry was proposed as an application of ultrasound computed tomography (USCT) shortly after its introduction in the 1970s [4]. In the 1980s, time-of-flight tomography based on linear propagation paths was investigated for internal temperature monitoring in industrial processes [154–157]. At present, a sole vendor markets a US Food and Drug Administration (FDA)-approved breast ultrasound tomography system in the USA [158,159]. The original device approval was based on reflection imaging predicates, but recently a colour elastographic image (qualitative for tissue stiffness) derived from the quantitative SOS and attenuation images was also cleared. While the underlying SOS images are sensitive to temperature changes, the quantitative images are evidently not available to the clinical observer. Fortunately, there is much interest in the field for both FDA approval and clinical indication for quantitative transmission imaging.

Ultrasound tomography employs a 360° aperture to collect either reflection (pulse–echo) and/or transmission (pitch–catch) data. From a k-space perspective, the reflection data fills an annulus in the frequency domain and is therefore high-resolution but non-quantitative. Since the reflection ultrasound tomography image is essentially compounded over 360°, there is significantly reduced speckle and shadowing compared to B-mode imaging. The transmission data, however, contains quantitative information on both the SOS and attenuation. With a suitable image reconstruction algorithm applied to the tomographic data, a map of the complex index of refraction can be generated for the tissue of interest, with the real component representing the spatially varying SOS and as the imaginary, the attenuation coefficient. In advanced ultrasound tomography systems, the SOS map can be used to improve the reflection tomography image [160–163] which may have artefacts due to beamformers that incorrectly assume a constant SOS. In other cases, full wave inversion can be performed using all of the scattering data. It should also be noted that reflection and transmission tomography can be performed with the same transducer if a strong reflecting surface is placed on the opposite side of the object [164–166].

The most significant disadvantage of ultrasound tomography-based thermography is the inferior temporal resolution. Most systems use a gantry to rotate reflection and transmission transducers around the object to produce a 2D slice. Motion of the gantry in elevation then produces a stack of images in three dimensions. To generate new SOS images in intervals of less than 5 s would require the engineering of non-rotating, high channel count tomographic systems. The alternative is to attempt the tomographic reconstruction using less data. Zhu et al. proposed and simulated limited angle USCT for temperature measurements in thermal engineering [167]. In their approach, the missing projection data is extrapolated using the X-ray transform consistency conditions for a time-of-flight (TOF), no-refraction ultrasound propagation model. It may also be possible to exploit the sparsity of temperature changes in the field of view by using an approach similar to keyhole MRI or compressive sensing [168].

Using a prototype USCT system, Azhari recently investigated the feasibility of HIFU guidance using both ultrasound tomography and projection-based transmission images [196]. For tomographic guidance, a baseline SOS image is formed and then segmented into three compartments based on SOS values: fat, non-fat, and tumour (in order of typically increasing SOS). Immediately after therapy, a second USCT image is acquired, and the SOS differential from baseline (ΔSOS) is calculated on a pixel-by-pixel basis. Using values taken from the literature, the temperature change is related to the change in SOS by a compartment-dependent coefficient:

Since this linear relationship is inexact, this method is only applicable to non-ablative hyperthermia less than 45 °C. For ablation application, Azhari proposes that the SOS is tracked over time, and when the SOS begins to decrease in non-fat tissues one may infer that ablation has succeeded. Given the aforementioned poor temporal resolution, this ∂SOS/∂t approach may be difficult in practice. Azhari proposes a faster approach by abandoning rotation for tomography and instead forming a 2D TOF projection in a manner similar to phase contrast techniques above. This acquisition requires only a movement in elevation, and in a first approximation the change in TOF is a function of the integrated temperature change along a line of propagation. The quantitative accuracy of this approach will be a function of volume of the zone of treatment with best performance with the smallest HIFU therapy volumes.

Passive monitoring

Passive monitoring strategies focus on the detection of bubble activity associated with ultrasound ablation. In general, there are thought to be two sources of bubble activity: cavitation and boiling. A less explored source of acoustic emission is broadband emissions associated with other phase transitions such as tissue coagulation.

In a groundbreaking study, Mast and colleagues [48] investigated the temperature dependence of three different types of acoustic emissions received by a 1 MHz unfocused passive cavitation detector (PCD) during HIFU in fresh bovine liver at 3.1 MHz. The PCD was oriented orthogonal to the HIFU beam produced by an array, and acoustic emissions were filter bank processed for a subharmonic band at 1.55 MHz, a 0.3–1.1 MHz broadband, and low frequency band below 30 kHz. The authors conjectured that the subharmonic band would be indicative of non-linear bubble vibrations from stable cavitation activity, while the broadband and low frequency channels would monitor inertial cavitation and tissue boiling respectively. Acoustic emissions captured by the PCD were correlated with increasing temperature measured using an invasive thermocouple. While the low frequency and broadband acoustic emission were generally observed to increase as a function of temperature, the subharmonic emissions were not a reliable indicator of ablation rate or tissue temperature. The following year, Mast and colleagues reported passive acoustic emissions with a 192-element array operating around 520 kHz [169]. Using a novel frequency domain beamformer (Equation 10 in previous reference [169]), they passively imaged sources in simulations of PBS solution, bovine liver, and in an in vitro flow phantom with echogenic liposomes as cavitation nuclei [170]. Although the linear array was 192 elements, experiments were limited to sub-apertures of 64 elements. The authors proposed that changes in temperature during HIFU therapy may alter the thresholds for inertial cavitation and the related processes such as rectified diffusion, and they offer a method for generating a ‘worst case’ estimate for temperature elevation.

Rabkin et al. compared the appearance time of B-mode hyperechogenicity with inertial cavitation measured with broadband PCD (0.5–5.5 MHz). They defined the onset of inertial cavitation as ‘the time of the first HIFU pulse in which the maximal amplitude of the broadband noise was greater than three Standard Deviations (SD) above mean broadband noise measured during recordings of acoustic emissions before each HIFU exposure’ [27, p. 950]. Above HIFU intensities of 800 W/cm², they reported the inertial cavitation detection with PCD to precede the visualisation of a hyperechoic region on B-mode by up to 0.5 s.

The Coussios group at Oxford [171,172] also considered PCD, but using a 128-element, 38-mm aperture, 5–10-Mhz diagnostic probe (L10-5, Zonare Medical Systems, Mountain View, CA) as the passive cavitation detector. Rather than utilising dynamic receive beamforming with a null transmit profile, Coussios and colleagues employed a novel time-exposure acoustics (TEA) beamforming algorithm that is well-suited for passive imaging of broadband and near-field sources. They evaluated the PCD using both non-cavitating (active) sources and single bubble inertial cavitation before tests with larger contiguous and disjoint cavitating regions in agar phantoms. Compared to real human tissue where the energy deposition tends to migrate towards the HIFU transducers, the Coussios gel phantom experiments exhibited bubbles that tended to migrate away from the HIFU transducer, according to the authors due to acoustic radiation force and phantom melting. In further studies in freshly excised, degassed bovine livers, Coussious and colleagues determined an 84% correction detection rate for lesion formation with all HIFU exposures over 5.4 MPa [173]. In comparison, only 53% of B-mode images with hyperechogenicity correctly mapped the lesion formation. While the lateral resolution of the PCD point spread function (PSF) for their custom beamformer was good (1 mm), the axial PSF was predicted and observed to be in the order of 5 mm. The authors expressed a hope to reduce this limitation to below 2 mm through future developments. Indeed, in more recent work using compressive sensing beamformers [174] and the robust Capon beamformer [175–177], improved axial spatial resolution has been achieved. Recent implementation of the renamed passive acoustic mapping (PAM: formerly called PCD) approach on an NVIDIA GTX260 graphics processing unit (GPU) produced real-time cavitation monitoring at 5 Hz [178].

In addition to ultrasonics, efforts have been made to detect acoustic emissions solely below 20 kHz due to boiling. Anand and Kaczkowski [101,179] used a medical stethoscope coupled by a standard PC sound card to a microphone sampled at 44.1 kHz. Using this approach, time-domain output from the stethoscope was acquired for 45 s continuously, even during small time periods when the HIFU system was turned off to acquire B-mode images. After 30 s, the onset of acoustic emission below 2 kHz was indicative of boiling and a time to onset of boiling (t_boil) was estimated. Using this boiling time, the known spatial HIFU beam profile, and the heat transfer equation (no cooling due to perfusion) in alginate and turkey breast phantoms, the local heating rate was then estimated. Since the stethoscope is an unfocused sensor, spatial mapping of the local heating rate was not possible. While Anand and colleagues utilised linear acoustics in their modelling of acoustic emission due to boiling, the Coussios group [180] has used the non-linear Khokhlov–Zabolotskaya–Kuznetsov (KZK) equation to model cavitation-enhanced heating, including the experimentally observed sharp temperature rises.

Lastly, since MRI is generally assumed to be incapable of detecting cavitation, Arvanitis and McDannold have reported attempts to co-register both MR temperature imaging and PAM in order to simultaneously monitor cavitation and temperature during HIFU [181]. To accomplish this, a custom 10 -m extension cable for a 128-element clinical imaging probe previously determined to be MR compatible was fabricated.

As stated above, passive cavitation detection and source localisation with tomographic apertures has received limited attention [182–185].

Hybrid methods

A variety of hybrid imaging modalities have been explored as potential solutions to the thermometry and ablation monitoring problems. In work by the Emelianov group at the University of Texas [103,186], the contrast in photoacoustic (optoacoustic) imaging (PAI) was theoretically taken to be a function of the dimensionless Grüneisen parameter Γ which in turn is a function of both the thermal coefficient of volume expansion β and the squared SOS. Both of these latter parameters are temperature dependent in a nearly linear fashion from 37–55 °C. Emelianov and colleagues then showed that a change in tissue temperature T would be proportional to changes in the normalized photo acoustic signal P:

Furthermore, they argued that the use of photoabsorbers such as gold nanoparticles as contrast agents would not change this relationship [187]. In comparison to the small changes seen on ultrasound strain images with temperature, changes of up to 40% in the photoacoustic signal were observed for 9 °C temperature changes. Wang and colleagues reported similar success applying photoacoustic thermometry in phantoms [188] and chicken tissue [189]. In a dual modality experiment with ultrasound transmission tomography, a group at the University of Twente (Netherlands) combined a photoacoustic tomography (PAT) image with SOS tomography to measure temperature increase of up to 9 °C [190].

For ablation monitoring, a reduction in the photoacoustic signal within the ablated zone has been reported with a preclinical PAI system [191]. The authors attribute the reduced signal to a ‘reduced concentration of blood, tissue dehydration, denaturation of proteins and porphyrins, and reduction of thermoacoustic efficiency in the thermally treated tissue’ [191]. However, other researchers have reported an increased PA signal in images acquired at 720 and 845 nm in HIFU-induced thermal lesions in chicken breast [192].

One potential advantage of PAT thermometry is high spatial resolution compared to other methods [189]. With a 3D tomographic geometry, Chitnis et al. reported an impressive isotropic spatial resolution of 0.5 mm [191]. In terms of time, though, the required signal averaging produces a 2-s temporal resolution which is not the fastest thermometry method, but sufficient for in vivo applications [188]. Due to optical attenuation, imaging at depths greater than 1 cm may require further trade-off in temporal resolution.

In acousto-optics (AO), limited studies have demonstrated a reduction in AO signal during formation of a thermal lesion, but the mechanism of modulation is not well-understood [193].

In an experimental system capable of measuring both electrical impedance and ultrasound velocity in phantoms undergoing mild hyperthermia, Islam et al. describe mixed reports on the accuracy of temperature estimation using combined measurements, especially when the percentage of fat in tissue was in error [194].

Conclusions

As summarised in Table 1, ultrasound methods have made significant progress in the areas of thermometry for hyperthermia below 50 °C and in qualitative ablation monitoring which targets tissue coagulation and bubble formation. There is, however, a dearth of evidence that robust ultrasound thermometry exists for temperatures above 50 °C. This shortcoming is at present adequately filled by the availability of mature MR thermometry techniques. While Nakagami imaging, axial shear strain elastograms, and echo decorrelation imaging have all shown promise for improved performance at higher HIFU power levels and temperatures, it is unknown whether any of these new techniques has the quantitative precision required for real-time control.

Table 1. Summary of various methods for ultrasound thermometry and ablation monitoring. For our purposes, current clinical equipment includes ultrasound research systems such as Verasonics and Ultrasonix that utilise FDA-approved probes. Techniques that require accessory equipment such as mechanical stimulators, lasers, or custom transducers are not included in that category even though they may include a clinical or research ultrasound system in the method.

Methods	Relative thermometry	Absolute thermometry	Quantitative. 37–50 °C	Quantitative. >50 °C	Ablation monitoring	Detects cavitation and boiling	Implementable usingcurrent clinical equlpment
PRF MR thermometry	✓		✓	✓			✓
Spectroscopic MR thermometry	✓	✓	✓	✓			✓
CT thermometry	✓			✓			✓
B-mode imaging					✓	✓	✓
Nakagami imaging					✓		✓
Thermal strain imaging	✓		✓				✓
Shear wave thermometry	✓		✓		✓
Local harmonic motion					✓
Stimulated acoustic emissions	✓		✓
Backscattered energy and QUS	✓		✓				✓
Transmission imaging and ultrasound tomography	✓	✓	✓		✓
Passive monitoring					✓	✓	✓
Photoacoustics	✓		✓

One possible route for improvement is the combination of two or more techniques into a multi-parametric analysis. For example, Chen et al. [195] recently described a hybrid method for temperature measurement that combines changes in B-mode images with thermal strain imaging. For temperature deviations less than 10 °C they reported a linear relationship between their combined index and temperature, but above 60 °C non-linearity was still unavoidable in pork tenderloin. It is unknown at this time, however, whether any combination of the methods described can expand the quantitative thermometry regime to treatment temperatures above 50 °C. With elastography modes now routinely available on clinical scanners, ultrasound thermometry based on thermal strain imaging may be the most appropriate technique for future clinical evaluation for thermal therapies above 50 °C, especially if coupled with a secondary system for cavitation detection.

Declaration of interest

M.A.L. acknowledges support from NIH-NCI 1R01CA138536-01A2. R.M.S. is a paid employee of Philips Research. R.C. is a Cancer Prevention & Research Institute of Texas (CPRIT) Scholar in Cancer Research supported by NIH-NCI Grant R1308 from the Cancer Prevention and Research Initiative of Texas and has also received generous support from the M.R. and Evelyn Hudson Foundation. The authors alone are responsible for the content and writing of the paper.