Hyperthermia is a cancer treatment that uses the biological effect of artificially induced heat to negatively impact cancer growth . The fundamental aim of a microwave hyperthermia treatment consists in selectively elevating the tumor temperature to 42–43 C by means of a proper antenna system, while keeping the heat in the surrounding healthy regions in a tolerable range. The targeted increase of the tumor temperature has been demonstrated to sensitize cancer cells to radiotherapy as well as to enhance the delivery of chemotherapy drugs to the tumor, proving particularly effective in the treatment of recurrent cancers [2, 3, 4].
In the state of the art for internal tumors, the microwave heating of the target is achieved by exploiting constructive wave interference from a certain number of antennas, laid out around the patient [5, 6, 7]. Due to loss mechanism of wave physics, the heating of the skin is inevitable; to avoid this, a plastic or silicon bag filled with circulating demineralized water (the so called water bolus) is introduced between the applicator and the body; this has the dual function of coupling the electromagnetic field into the patient and to cool the body surface, preventing superficial overheating [8, 9].
The effectiveness of a hyperthermia treatment is strongly dependent on the quality of the heating process, especially for challenging anatomical sites such as the head and neck (H&N) region; this has prompted advances in treatment planning to optimally set the amplitudes and phases of the applied signals. To this purpose, patient specific 3D models, based on the acquisition and segmentation of CT or MRI scans, together with the applicator model, are imported in appropriate numerical solvers; the optimization on the antenna feedings is then performed on the basis of the simulate energy deposition in the patient.
Maximizing the temperature inside the tumor is the ultimate goal of an hyperthermia treatment ; however, a temperature-based optimization requires a very high computational burden, involving the (numerical) solution of the bioheat equation at each step. Because of this, treatment planning usually relies on the optimization of the specific absorption rate (SAR), which is the heat source (volume density of deposited power). This significantly eases the required computational cost; because of its practicality, and its being an indirect measure of the temperature reached, a systematic study can be found in the literature where the relationship between SAR coverage indicators and clinical outcome is assessed . For these reasons, different SAR-based optimization strategies have been proposed [12, 13, 14], and actually adopted in the clinic for applicators in the H&N region .
The correlation between SAR and the corresponding temperature peak locations has also been studied theoretically ; that study is based on the assumption of a homogeneous infinite medium, which is reasonable in many cases. However, the thermal boundary conditions imposed on the skin-water bolus interface may significantly impact on the temperature profile in realistic scenarios, thus leading to a shift in the spatial peaks of temperature and SAR distribution.
As reported below, we have analyzed this aspect in both simplified and realistic 3D models and simulations; we have found indeed temperature-SAR shifts, attributable to water bolus boundary conditions. We propose here an easy-to-implement strategy for correcting this, while still employing SAR optimization. This is achieved by finding the SAR focusing location that leads to a temperature focusing in the desired region. The presented approach can be used to improve the performance of the existing SAR-based optimization processes, such as , or the recently proposed focusing via constrained power optimization [15, 16], especially when applied to sub-superficial tumors .
The paper is organized as follows. In Section II, the equation governing the heat transfer in tissues is briefly introduced. In Section III, we describe the implemented numerical testbed, concerning a tumor in the H&N region, the field focusing technique adopted to maximize the power deposition on the target and the corresponding thermal analysis. The proposed corrective procedure on the SAR focusing is presented in Section IV. The achieved results are reported in Section V, both for a simple model of the human neck and for a more realistic numerical phantom. The conclusions follow in Section VI.
2 Bioheat Transfer and Boundary Conditions
The temperature increase in a tissue caused by the exposure to an external heating source can be described as usual using the Pennes’ bioheat equation ; as common in the hyperthermia planning literature , we will be concerned with the steady-state version of the equation:
where is the temperature, is the tissue-specific thermal conductivity, and is the source term, given by:
In (2), is the external heat source, and are the blood density and specific heat, is the arterial blood temperature, and is the volumetric blood perfusion rate. The metabolic heat generation term has been neglected in (2) as usual. In the presence of an electromagnetic source (such as an antenna system), the heat source can be written as:
where is the electric field (peak value), and are the tissue-specific mass density and electrical conductivity, and SAR denotes the specific absorption rate.
The thermal model has to be completed by the boundary conditions (B.C.) for the temperature , specifying how the system interacts with the outside environment.
For the realistic case of a finite region, two basic types of thermal boundary conditions can be identified: specified temperature (isothermal condition) and specified heat flux. The isothermal condition is appropriate when the body is in contact with a good heat conductor maintained at a constant temperature
(such as an ice pad) and it is given by:
where indicates the boundary of the considered volume region . The effects of an exterior fluid cooling the surface of the body are instead better described by a convective heat flux boundary condition:
is the unit vector normal to the boundary,is the local thermal conductivity, is the external temperature and is the heat transfer coefficient.
3 SAR focusing and temperature shift
3.1 Numerical testbed
To address the issue of the temperature focusing shift in an hyperthermia treatment, we implemented a numerical model concerning the application of the microwave heating to a tumor placed in the human neck.
For the sake of illustrating issues and the proposed method, we first consider a simplified computational phantom of the human neck, where the vertebrae, the trachea, the spinal cord, as well as the neck shape are modeled by means of cylinders, while the tumor is represented by a sphere of diameter mm; this is shown in Fig. 1. We will then test our method to the realistic case described in Section 5.2.
A phased array applicator is very effective for the heating of H&N tumors [6, 7], and the reference design is reported in literature . Hence, we considered a circular array made of patch antennas with water substrate, immersed in a water bolus and properly matched to operate at the frequency MHz . The dielectric and thermal properties were assigned according to the material database in , using a tumor perfusion value 1.85 times higher than the rest value for the muscle perfusion  (see Table 1).
The model boundary conditions were set as per the current literature . A convective heat flux boundary condition (5) was applied on the tissue-water bolus interface, with C and , and on the internal boundaries of the laryngotracheal canal, with C  and . Moreover, the isothermal condition C was imposed on the upper and lower surfaces of the neck, being reasonable to assume that the temperature of the head and the body trunk remains approximately constant during treatment in real applications.
3.2 SAR optimization
We considered the SAR-based optimization strategy currently usually employed in the literature . The implemented approach seeks the antenna feedings that maximize the SAR in the gross target volume (GTV), i.e. the tumor, minimizing the risk of hotspots in the surrounding healthy tissues. The corresponding cost function is the target-to-hotspot SAR quotient (THQ) :
3.3 SAR-Temperature shift
The resulting temperature map is reported in Fig. 2.c. As can be observed, a significant temperature shift outside the GTV occurs. This shift can be attributed to the thermal boundary conditions dictated by the external scenario, i.e. the water bolus cooling.
4 Temperature shaping via SAR focusing
In this section, we concentrate on effectively shaping the temperature field, avoiding the shift between SAR focusing and temperature distribution described in the previous section. We would like to keep the SAR focusing approach, though, because it is numerically effective, widespread and assessed; we will thus proceed using the SAR focusing, but carrying out the SAR focusing task in a manner that achieves the desired temperature distribution.
We do this by finding where the SAR should be focused in order to have the desired temperature coverage around a given wanted position. In its general setting, this amounts to an optimization problem; hence, the first issue to be addressed is the definition of a meaningful objective function.
In practical applications, the requirements on the temperature distribution involves also uniformity in addition to positioning the maximum; the corresponding quality estimator is the T90 parameter[28, 26], defined as the temperature exceed by of the points in the GTV. Hence, we will use this quality indicator in our procedure.
It is apparent that a good shaping involves both temperature “focusing” (i.e., maximizing in the desired GTV) and temperature uniformity; this is achieved by optimizing a single parameter: the T90 normalized to in the GTV, i.e via:
Expression (7) will be our objective function.
Another important issue is to decouple this part of the optimization from the SAR optimization (e.g. as in Section 3.2); this can be effectively done recognizing that in the GTV the electric field distribution is well approximated by a Gaussian function. We will set the width and peak value of the Gaussian approximating field by fitting them to the squared amplitude of the field distribution of the initial SAR optimization, and leaving only the center position as the optimization parameter, with the above-defined objective function .
We observe that the Gaussian approximation is done on the squared amplitude of the electric field and not on the SAR itself. This allows a greater flexibility, as the temperature-SAR shift may involve crossing regions with different electrical conductivities.
The Gaussian approximation of the electric field distribution is of course valid only locally, but is a good approximation in the region of interest; this is seen by comparing the value of our objective function, , for actual SAR distributions and their Gaussian approximation, as reported in Fig. 3.
As for the optimization, any local (convex) approximate optimization would do, as we are looking for a correction on an initial global optimization for the SAR. In this work, we have chosen the crudest approach, that has the advantage of simplicity and can be used as an add-on onto any existing treatment planning system. We look for the optimal SAR peak location in a region around the desired temperature peak point. This will be called the refinement region , surrounding the GTV. In this limited region, we discretize the space, compute the objective function for all points in the grid, and select the best point, i.e. the one with the highest value of the objective function.
We now describe the process in more detail.
4.2 Gaussian fitting
As alluded above, we fit the squared amplitude of the electric field, , with a multivariate Gaussian function with height , center , and width :
This will allow to insert a source term in the bioheat equation (1) to perform the optimization of the SAR focusing position, without the need to perform SAR optimization each time.
The suitability of this local approximation is readily assessed by checking the final metric of the overall optimization exercise, i.e. via the objective function of the optimization, , defined above in (7). This can be seen in Fig. 3, where the value of the objective function is evaluated for the case where the bioheat source term is the actual SAR optimized for a given position, and when it is substituted for the Gaussian approximation in (8). As it is seen there, the objective function value is well preserved by the Gaussian approximation of the electric field squared amplitude; the differences between the temperature profiles (corresponding to the exact squared amplitude of the E-field and its Gaussian approximation) only occur moving away from the region of interest.
4.3 Optimum focusing point search
Thanks to the Gaussian approximation of the electric field squared amplitude, using (8), the center of the Gaussian SAR function can be easily moved within a region surrounding the GTV and the corresponding temperature map can be calculated accordingly, solving the bioheat equation (1); for each position of the SAR center, one thus computes the temperature distribution, and from this the objective function (7). From the list of these values one picks the SAR peak location which yields the largest value of the objective function.
In this exhaustive search, it is necessary to establish a) the extension of the search region, called refinement region , and b) the sampling density in this region. The distance between the SAR and temperature peaks appears as a meaningful characteristic length of our problem.
It is also reasonable to consider a refinement region with the shape of a ball centered around the tumor centroid , with a diameter , being the radius of the smallest sphere enclosing the GTV. The volume of can be discretized with a sampling distance among points, where indicates the sampling density.
Using an exhaustive search approach, the SAR map:
with and fixed by the Gaussian fit described in Section 4.2, is generated for each point , i.e., the Gaussian function center is moved on the grid of sampling points inside , the bioheat equation (1) is solved, and the parameter is evaluated. The effect of this search, here limited to a two-dimensional ball surrounding the GTV on the plane, is shown in Fig. 4.a.
5 Results and Discussion
In this section, we report the results achieved applying the temperature-SAR shift corrective procedure described in Section 4; we consider both the above-introduced simple phantom and a realistic model of the human neck region.
5.1 Simple neck model
For the case under study, we test the proposed strategy limiting the refinement region to a two-dimensional ball on the plane of the SAR-temperature shift; with our coordinate system this is on the plane and centered at mm, with a diameter cm, being mm the tumor sphere diameter and mm the SAR-temperature shift magnitude (see Fig. 3.a).
Following the approach described in Section 4.2, a bi-variate Gaussian mask with V/m and cm has been created to fit the squared amplitude of the E-field corresponding to the SAR distribution reported in Fig. 2.b. The normalized Gaussian fitting function and the corresponding temperature profile, obtained by solving the bioheat equation (1), are shown in Fig. 3.a. The good agreement among the curves within the region of interest (i.e., in proximity of the tumor) proves the effectiveness of the implemented fit.
Computation of the objective function was then performed for SAR center points placed at the nodes of a grid of 347 evenly spaced points on the circular refinement region , with a inter-node spacing mm. For each grid point, the temperature map was computed solving the bioheat equation (1), and the corresponding parameter was evaluated.
Fig. 4.a shows the parameter as function of the centers of the Gaussian mask inside the refinement region. The position corresponding to the maximum provides the Gaussian SAR center that optimizes the uniformity of the tumor temperature coverage. In the considered case, it was found: mm.
The antenna feedings optimization, performed on the 3D model according to Section 3.2, maximizing the SAR in a target sphere of diameter centered around , provides the temperature map displayed in Fig. 4.b on the plane. By comparing Fig. 2.c and Fig. 4.b, it appears evident how the implemented shift correction has significantly improved the temperature coverage of the GTV.
Fig. 4.c reports the temperature maps on the other planes ( and ) cutting the tumor along the canonical directions before and after the application of the above-described temperature refinement. For the considered example, although the implementation of the shift correction was limited to the plane, a significant improvement in the temperature coverage of the GTV was also obtained in the plane.
5.2 Realistic neck model
The temperature shift correction procedure was then repeated for the more realistic neck computational model visualized in Fig. 5. This phantom has been realized using the 3D mesh models of the vertebrae (from C3 to C7), the spinal cord and the laryngotracheal canal provided by the Visible Human Project (VHP Female Version 2.2 ). In this example, the GTV has been simulated as a mass with irregular shape, centered at mm.
As per our proposed method, in the first step, the SAR-based optimization aimed at maximizing the SAR coverage of the tumor was performed.
Then, using the same thermal boundary conditions as those reported in Section 2 for the simple model, the corresponding temperature map was obtained (see Fig. 6.a), showing a shift in the temperature focusing outside the GTV. The corrective strategy described in Section 4 was applied again on a two-dimensional refinement region on the plane, centered at , with a diameter cm and a sampling distance of 1.3 mm among points.
The search for the optimal Gaussian SAR peak location, i.e. providing the largest value of the objective function , is reported in Fig. 6.b; the position corresponding to the maximum was found to be mm.
To ensure the therapeutic effectiveness of a hyperthermia treatment it is important to optimally focus the temperature increase on the tumor, with minimum collateral damage to the surrounding tissues. In the clinical practice, this is currently achieved via SAR-based sophisticated optimization routines, employing accurate numerical simulations where both the applicator and the CT-based patient model are included.
In this paper, we highlighted how SAR optimization on the target region could be not enough to ensure an equally satisfying temperature focusing, due to the thermal boundary conditions imposed by the external scenario (e.g. water bolus).
To address this problem, we proposed a simple strategy to compensate the temperature peak shift outside the tumor via a proper displacement of the SAR focusing target. The implemented procedure has been applied to improve the heating of a tumor mass in the head and neck region, providing satisfying results both in the case of a simple neck phantom, and for a more realistic model. For both the considered examples, a two-dimensional search for the optimal shift proved to be enough for the current correction purposes. However, a three-dimensional search could lead to further improvements in the final temperature coverage of the tumor for antenna arrays that allow a major control of the focusing process along the vertical axis (e.g., when the antennas are organized in more than one ring).
It is worth noting that the proposed strategy represents a way to improve the temperature focusing on the target in the context of a SAR-focusing procedure, both in pre-treatment and during in-treatment adjustments, potentially leading to an overall enhancement of the treatment quality with a small additional computational effort.
This work has been supported by the Italian Ministry of Research under PRIN “Field and Temperature Shaping and Monitoring for Microwave Hyperthermia FAT SAMMY”.
-  M. W. Dewhirst, B. L. Viglianti, M. Lora-Michiels, M. Hanson, and P. J. Hoopes, “Basic principles of thermal dosimetry and thermal thresholds for tissue damage from hyperthermia,” Int. J. Hyperthermia, vol. 19, no. 3, pp. 267–94, 2003.
-  H. H. Kampinga, “Cell biological effects of hyperthermia alone or combined with radiation or drugs: a short introduction to newcomers in the field,” Int. J. Hyperthermia, vol. 22, no. 3, pp. 191–6, 2006.
-  N. Cihoric et al., “Hyperthermia-related clinical trials on cancer treatment within the ClinicalTrials.gov registry,” Int. J. Hyperthermia, vol. 31, no. 6, pp. 609–14, 2015.
-  N. R. Datta et al., “Local hyperthermia combined with radiotherapy and-/or chemotherapy: Recent advances and promises for the future,” Cancer Treat. Rev., vol. 41, no. 9, pp. 742–53, 2015.
-  M. M. Paulides et al., “The HYPERcollar: A novel applicator for hyperthermia in the head and neck,” Int. J. Hyperthermia, vol. 23, no. 7, pp. 567–76, 2007.
-  J. Crezee et al., “Improving locoregional hyperthermia delivery using the 3-D controlled AMC-8 phased array hyperthermia system: a preclinical study,” Int. J. Hyperthermia, vol. 25, no. 7, pp. 581–92, 2009.
-  T. Drizdal, M. M. Paulides, N. van Holthe, and G. C. van Rhoon, “Hyperthermia treatment planning guided applicator selection for sub-superficial head and neck tumors heating,” Int. J. Hyperthermia, vol. 34, no. 6, pp. 704–13, 2018.
-  M. L. van der Gaag, M. De Bruijne, T. Samaras, J. van der Zee, and G. C. van Rhoon, “Development of a guideline for the water bolus temperature in superficial hyperthermia,” Int. J. Hyperthermia, vol. 22, no. 8, pp. 637–56, 2006.
-  Z. Rijnen, P. Togni, R. Roskam, S. G. van de Geer, R. H. Goossens, and M. M. Paulides, “Quality and comfort in head and neck hyperthermia: A redesign according to clinical experience and simulation studies,” Int. J. Hyperthermia, vol. 31, no. 8, pp. 823–30, 2015.
-  H. P. Kok, P. M. van Haaren, J. B. van de Kamer, J. Wiersma, J. D. van Dijk, and J. Crezee, “High-resolution temperature-based optimization for hyperthermia treatment planning,” Phys. Med. Biol., vol. 50, no. 13, pp. 3127–41, 2005.
-  H. K. Lee et al., “Superficial hyperthermia and irradiation for recurrent breast carcinoma of the chest wall: prognostic factors in 196 tumors,” Int. J. Radiat. Oncol. Biol. Phys., vol. 40, no. 2, pp. 365–75, 1998.
-  E. Zastrow, S. C. Hagness, B. D. van Veen, and J. E. Medow, “Time-multiplexed beamforming for noninvasive microwave hyperthermia treatment,” IEEE Trans. Biomed. Eng., vol. 58, no. 6, pp. 1574–84, 2011.
-  Z. Rijnen et al., “Clinical integration of software tool VEDO for adaptive and quantitative application of phased array hyperthermia in the head and neck,” Int. J. Hyperthermia, vol. 29, no. 3, pp. 181–93, 2013.
-  G. Cappiello et al., “Differential evolution optimization of the SAR distribution for head and neck hyperthermia,” IEEE Trans. Biomed. Eng., vol. 64, no. 8, pp. 1875–85, 2017.
-  D. A. M. Iero, L. Crocco, and T. Isernia, “Thermal and microwave constrained focusing for patient-specific breast cancer hyperthermia: A robustness assessment,” IEEE Trans. Antennas Propag., vol. 62, no. 2, pp. 814–21, 2014.
-  G. G. Bellizzi, T. Drizdal, G. C. van Rhoon, L. Crocco, T. Isernia, and M. M. Paulides, “The potential of constrained SAR focusing for hyperthermia treatment planning: analysis for the head & neck region,” Phys. Med. Biol., vol. 64, no. 1, 2018.
-  H. H. Pennes, “Analysis of tissue and arterial blood temperatures in the resting human forearm,” J. Appl. Physiol., vol. 1, no. 2, pp. 93–122, 1948.
-  M. M. Paulides et al., “Simulation techniques in hyperthermia treatment planning,” Int. J. Hyperthermia, vol. 29, no. 4, pp. 346–57, 2013.
B. Gao, S. Langer, and P. M. Corry, “Application of the time-dependent Green’s function and Fourier transforms to the solution of the bioheat equation,”Int. J. Hyperthermia, vol. 11, no. 2, pp. 267–85, 1995.
-  P. A. Hasgall et al., IT’IS database for thermal and electromagnetic parameters of biological tissues, 2018. [Online]. Available: http://www.itis.ethz.ch/database
-  D. Andreuccetti, R. Fossi, and C. Petrucci, An Internet resource for the calculation of the dielectric properties of body tissues in the frequency range 10 Hz - 100 GHz. IFAC-CNR, Florence (Italy), 1997. [Online]. Available: http://niremf.ifac.cnr.it/tissprop/
-  M. M. Paulides, J. F. Bakker, N. Chavannes, and G. C. van Rhoon, “A patch antenna design for application in a phased-array head and neck hyperthermia applicator,” IEEE Trans. Biomed. Eng., vol. 54, no. 11, pp. 2057–63, 2007.
-  J. Lang, B. Erdmann, and M. Seebass, “Impact of nonlinear heat transfer on temperature control in regional hyperthermia,” IEEE Trans. Biomed. Eng., vol. 46, no. 9, pp. 1129–38, 1999.
-  R. F. Verhaart et al., “Accurate 3D temperature dosimetry during hyperthermia therapy by combining invasive measurements and patient-specific simulations,” Int. J. Hyperthermia, vol. 31, no. 6, pp. 686–92, 2015.
-  E. R. McFadden Jr. et al., “Thermal mapping of the airways in humans,” J. Appl. Physiol., vol. 58, no. 2, pp. 564–70, 1985.
-  R. A. Canters, P. Wust, J. F. Bakker, and G. C. van Rhoon, “A literature survey on indicators for characterisation and optimisation of SAR distributions in deep hyperthermia, a plea for standardisation,” Int. J. Hyperthermia, vol. 25, no. 7, pp. 593–608, 2009.
-  COMSOL Multiphysics, version 5.4. www.comsol.it, 2019.
-  R. F. Verhaart, V. Fortunati, G. M. Verduijn, A. van der Lugt, T. van Walsum, J. F. Veenland, and M. M. Paulides, “The relevance of MRI for patient modeling in head and neck hyperthermia treatment planning: A comparison of CT and CT-MRI based tissue segmentation on simulated temperature,” Med. Phys., vol. 41, no. 12, 2014.
-  S. N. Makarov et al., “Virtual human models for electromagnetic studies and their applications,” IEEE Rev. Biomed. Eng., vol. 10, pp. 95–121, 2017.