HIFU treatment time reduction in superficial tumours through focal zone path selection

Abstract

Purpose: This study evaluates the hypothesis that optimising the path of a high intensity focused ultrasound (HIFU) treatment's N focal zone heating pulses can significantly reduce treatment time, and identifies the underlying bio-thermal principles.

Materials and methods: Thirty-one scanning paths were investigated using 3D simulations, with a minimum thermal dose delivered to every tumour position. Treatment time was calculated as the sum of the N, individually optimised heating and cooling periods. Tumours were superficial (skin to tumour distance ranging from 1.3 to 2.5 cm), but always deep enough so that the pre-tumour normal tissue was routinely heated to its constraint temperature (range: 42−45°C). Properties were uniform and constant, and a range of blood perfusion and phased array powers were studied.

Results: The best paths significantly reduced treatment times, with the largest gains occurring when (1) temperature superposition inside the tumour was maximised by successively heating the focal zone positions located in a ‘stack’ along the transducer's axis, and (2) the focal zone was moved laterally to an optimised location and another stack was applied. Stacking takes advantage of the focal zone's elongated shape, which produces axial temperature superposition within the tumour. Reduced tumour heating times also reduced energy deposition in the normal tissues, thus reducing or eliminating the need for inter-pulse cooling.

Conclusions: HIFU treatment times can be significantly reduced by taking advantage of axial temperature superposition in tumours. Further reductions are obtained by correct choice of the transverse scan path.

• High intensity focused ultrasound
• MRgHIFU
• path optimisation
• simulations
• treatment time

Introduction

High intensity focused ultrasound (HIFU) has been shown to be a promising treatment method for several types of cancer and for uterine fibroids [3–8]. However, when considering large tumours embedded in critical normal tissues, the requirement of precisely placing and applying a large number of heating pulses can result in very long treatment times [3], [9], [10], with McDannold et al. [10] recommending at least 50–60 s of cooling between pulses and Fennessy and Tempany [3] using 20–40 s of heating followed by 80–90 s of cooling. For example, 100 pulses with 30 s of heating and 60 s of cooling per pulse would take 2.5 h just for delivering the desired dose, not including treatment preparation time, adjustments and complications. Such long times pose a major challenge to HIFU's broad clinical acceptance. To reduce treatment times, clinicians can control: the electrical power delivered to the transducer; the sequence of focal zone positions; the heating pulses’ periods; the inter-pulse cooling periods; and, for a phased array, the size and shape of the focal zone [11].

One way to reduce treatment times is to simultaneously optimise all of those parameters. However, because of the computational intensity involved in simultaneously optimising the large number of possible combinations of parameters, prior research [3], [4], [6], [10], [12–14] has commonly addressed a fixed treatment path and parametrically changed one of the other treatment variables, including the: applied power magnitude [13]; ultrasound frequency [13] and; inter-sonication delay period [10]. Most investigators have used fixed heating pulse and inter-pulse cooling periods to both provide standardised heating conditions and to avoid normal tissue overheating [3], [4], [6], [14–16]. Because of efficacy and safety concerns, such protocols used conservative heating and cooling periods, resulting in very long treatment times [6], [12]. To help reduce those times, recent research has developed an approach that minimises each individual pulse heating and inter-pulse cooling period during a treatment, resulting in significantly reduced treatment times [17]. That study also showed that larger applied powers always shortened treatments for focused transducers, a non-obvious result when normal tissue constraints are active. They employed the strategy of developing separate optimising approaches for individual treatment parameters (e.g. heating and cooling periods), with the goal of eventually combining all individual approaches into a comprehensive optimisation approach. The current study extends the use of that strategy by optimising one more treatment variable, the focal zone scan path.

Only a few studies have investigated the effects of scan path, and none have been concerned primarily with treatment time – except Payne et al. [17] which only compared two paths. Fan and Hynynen [18] considered three two-dimensional planar paths applied sequentially to two planes perpendicular to the transducer's axis. They used fixed heating and cooling periods and thus did not study treatment time; however, they did show that the volume of tissue ablated was path dependent. Liu et al. [19] performed a numerical study of two heuristic paths (least possible distance and greatest possible distance between sequential small rapid scanning volumes) and showed that cumulative heating time differences (cooling times were fixed at 1 min between sonications) of several minutes were possible, depending on which path was used. More recently Kohler et al. [20], in an experimental study that was extended by Enholm et al. [21], and was based on that group's earlier work with nested spirals [22], [23], used a focal zone at a fixed axial depth that was rapidly electronically steered in concentric 2D circles going sequentially from inward to outward. They used fixed heating times for their circles and thus they did not provide any comparative information on optimal treatment times. Finally, in a 2D and 3D (spherical tumour) computational study, Malinen et al. [24] developed an optimisation algorithm that found the scan path that minimised the total dose delivered to the tumour while guaranteeing minimal dose at all locations. Again, treatment time was not the study's focus, but a more uniform dose delivery would, presumably, help shorten treatments – a reasonable belief bolstered by comparison of their results to those from a ‘standard, non-optimised’ path.

In summary, while several investigators have studied the role of separate physical mechanisms in reducing treatment times, no study exists that concentrates on the effects of scan path on treatment time. Thus, although several different paths have been proposed and studied by individual investigators, at present no type of scan path has been shown to be preferable – regardless of differing claims by various investigators. The current research attempts to help shed light on this problem by providing direct comparisons of the treatment times needed by a range of scanning paths to treat the same tumours under simple, standard conditions of fixed, uniform tissue properties. Such results will help guide, and be a reference base for, future studies on more complex clinical situations.

Methods

The goals of this study were to determine whether significant treatment time gains could be obtained by optimising the scan path when significant normal tissue heating was present, and to identify the thermo-physical mechanisms underlying those gains. To do so, numerous treatment scenarios were simulated using different scanning paths, tissue perfusions, transducer powers, tumour sizes and depths, and normal tissue constraints. In this study we concentrate on superficial tumours, but with all such tumours located at a sufficient depth so that there was sufficient ‘pre-focal’ normal tissue between the skin and tumour so that normal tissue overheating due to scanning-related thermal build-up was a significant factor in differentiating among scan paths.

Treatment time

Each scan path involved a sequence of N discrete focal zone locations, at each of which a heating pulse was applied and a subsequent power off, inter-pulse cooling period was initiated (if needed) to avoid heating the normal tissue past its constraint value (generally 43°C) during the next heating pulse. Each path's treatment time (t_TREAT) was calculated as the sum of its N heating and N-1 inter-pulse cooling periods:

Here, Δt_Heat,_n and Δt_Cool,_n are the pulse heating period and the subsequent inter-pulse cooling period for the nth heating pulse. Each heating pulse was applied at the same total transducer power magnitude at every focal zone position within a tumour. In order to find the path that minimised t_TREAT a scan path was set, the pulse's heating and inter-pulse cooling periods were individually minimised at each of the N sequential focal zone positions, and t_TREAT was calculated and evaluated with respect to all other paths’ treatment times values. To find the minimal heating and cooling times at each focal zone position, the simulations used the approach of Payne et al. [17] in which each pulse's heating period was made just long enough so that the desired thermal dose was delivered to the focal zone (as determined both during that heating pulse and during the immediately ensuing decay of the focal zone temperature); and power was then turned off during an inter-pulse cooling period that was also optimised. This approach ensured that every path studied was standardised to consistently meet the same important conditions; the same applied power magnitude was used for every pulse, and every heating and cooling period was individually optimised. The thermal dose [25] was monitored in the hottest voxel of the full width half maximum region of the focal zone's specific absorption rate (SAR) pattern, and the inter-pulse focal zone spacing was made small enough to ensure that the tumour was treated to at least 240CEM43°C at all locations (as verified after each treatment). Any ‘prior dose’ delivered (i.e. delivered to the current focal zone position during earlier heating of other focal zone positions) was subtracted from the target dose before the optimal heating and cooling periods were calculated. Conversely, no accounting for future dose was implemented. That is, the dose that will be delivered to the currently heated position by later pulses was not anticipated, an approach that results in the (clinically conservative) overdosing of some tumour locations. The metric used for comparing scan paths, the treatment time, is equivalent to a rate of thermal ablation when a fixed tumour volume is divided by the treatment time.

Thermal model

Temperatures were calculated using the Pennes bio-heat transfer equation [26]:where T is temperature (°C), ρ is the tissue density (kg/m³), C is the specific heat of tissue (J/kg°C), k is the thermal conductivity of the tissue (W/m°C), Q_ap is the applied power density deposited by an external applicator (W/m³), W is the Pennes perfusion parameter (kg/m³ s), C_b is the specific heat of blood and T_b is the arterial blood temperature. Although this equation has its limitations, it has been shown to adequately predict the major aspects of in vivo experiments [26–29]. An explicit finite-difference approximation of Equation 2 was implemented with a temporal resolution of one second and a spatial resolution of 3 × 3 × 1 mm. These resolutions are similar to those in Malinen's study [24] which used an isotropic 2-mm finite element nodal spacing, and 0.5-s time steps. Their preliminary studies showed that finer spatial resolutions did not significantly affect their results, a finding similar to ours for the 1-s time step. Also as in that study: all boundary conditions were kept at 37°C, as were the initial conditions; tissue properties were constant, and their computational region size was chosen to include the major features of the scanning process (primarily the heating of the proximal and distal normal tissues, and the tumour), while being mindful of the computationally intense optimisation problems being solved. The present study used a comparable computational region size for the same reasons. We applied a range of normal tissue temperature limits [42°, 43°, and 45°C) on constraint planes located 1 cm behind and in front of the tumour – which is the location where the maximum normal tissue temperatures almost always occurred due to thermal build-up, a phenomenon that is minimal in the transverse direction due to the rapid transverse fall-off of the focal zone's SAR pattern.

Power deposition

The power deposition pattern was determined as in Payne et al. [17] using the HAS (hybrid angular spectrum) method [30] to simulate the SAR field from an existing ultrasound phased array with 256 elements arranged randomly on a spherical surface (Imasonics, Besançon, France), with the focal zone steered electronically in the axial direction and mechanically in the transverse direction (Image Guided Therapies, Bordeaux France). The HAS method is an extension of the traditional angular spectrum method to inhomogeneous tissues and has been found to be both accurate and much faster computationally than other methods (i.e. the Rayleigh-Sommerfeld integral) [30]. The resulting focal zone size, expanded to adjust for scattering and transducer element variability, is ∼2 × 10 mm (FWHM) which is very close to that seen when using this phased array to heat agar phantoms and ex vivo tissues [31]. The simulated treatment geometry is shown in Figure 1.

Figure 1. Geometry of the treatment domain. The origin is located in the centre of the transducer's face. The geometric focus of the transducer is at 13 cm and the water filled distance from the transducer to the skin is 11.4 cm (11.0 for some simulations, as measured from the centre of the face of the transducer). A small (1.1 cm³) and medium (2.6 cm³) volume tumour were studied with this geometry. A deeper medium and a larger (4.5 cm³) tumour were studied with a slightly modified geometry using an expanded simulation region.

Every tumour studied was treated with heating pulses at three focal zone depths along the transducer's axis (z direction). Thus, three electronically steered power deposition distributions were simulated for each tumour, one for each focal zone depth. Each of those distributions was translated in the x–y plane at its fixed depth to obtain the SAR distributions for the tumour's (N/3) transverse positions at that depth. This approach simulates mechanical scanning in the transverse directions. The focal zones’ peak absorbed power density values (W/m³) used were all below the conservative limits set by the FDA [32] to avoid cavitation. For every scan pattern the total applied transducer power was the same at all N locations, and was fixed at a value that gave a maximum absorbed power density (in the hottest voxel) of 0.67 × 10⁷W/m³ when the focal zone was centred on the front plane of the small, superficial tumour studied (Figure 1). This power density value is close to optimal in the sense that lower values yield significantly longer treatment times, while larger values give only small gains in treatment times, a prior result [17] that was confirmed in the current study (Figure 6). The property values used are shown in Table I.

Figure 2. Treatment time vs. perfusion for three focal zone paths treating the small superficial tumour with 45 positions with the standard total applied power and a range of perfusions. Also shown are the percentages of the treatment time spent cooling the normal tissue. The lines are for the AS (MBF), XY (Ra) (squares); 3D Kn (circles) and PL (BMF), XY Ra (triangles) scan patterns.

Table I.  Thermal and acoustic properties used in all simulations [1], [2].

Thermal conductivity, k (W/m/°C)	0.5
Specific heat, C (J/kg/°C)	4190 (blood and tissue)
Speed of sound, c (m/s)	1500 (tissue and water)
Attenuation coefficient, α (nP/m/MHz)	4.8
Density, ρ (kg/m^3)	1000
Pennes’ perfusion, W (kg/m^3/s)	0.5, 1.0–10

Tumour models

Four tumour geometries were studied (Table II), labelled as small superficial, medium superficial, medium deeper, and larger. The choice of geometries was made in order to provide (1) a reasonable range of sizes (the largest tumour dimension in the current study was ∼2 cm, a value near the middle of the range for the largest tumour dimension representing approximately 1/3 of breast tumours [33]), (2) tumours at several depths, including depths where significant axial electronic steering was needed, with its concomitant focal zone SAR reduction, and (3) tumour sizes and dimensions similar to those used in other path studies, e.g. Liu et al. and Malinen et al. [19], [24]. The tissue region dimensions for the small and medium superficial tumours were 3.3 × 3.3 × 3.1 cm, while for the medium deeper and larger tumours they were 3.3 × 3.3 × 5.1 cm. These region sizes are close to those used by Malinen et al. [24]. To verify that the region size was large enough to prevent boundary condition induced changes in scan path rankings, multiple additional simulations were performed using a larger, region up to 1.8 × 1.8 × 1.4 cm³. The larger region gave the same rank ordering of scan paths as the smaller region. The focal locations were spaced 3 mm apart in all three directions for all tumours except for the larger tumour, for which the focal zones were 6 mm apart in the z-direction. The 3-mm axial spacing had considerable axial SAR overlap among the three axial focal zone positions within a ‘focal zone stack’ and was thus a conservative spacing ensuring that all voxels within the tumour achieved a minimum dose of 240 CEM. This matches the recommendation of a previous study [18] that some focal zone overlap is present to avoid untreated regions in the tumour. This 3-mm spacing is on the same order as the 2-mm spacing between focal zones used by Fan et al. [34] (which is more conservative than the 5-mm axial and 3-mm lateral spacing used by Hynynen [6]), and is similar to the 4-mm inner annulus diameter used by Kohler et al. [20]. The 6-mm axial spacing used in the larger tumour provided a larger spacing that could treat the tumour in a reasonable amount of treatment time while still providing enough SAR overlap to avoid regions of untreated tumour tissue between focal zones.

Table II.  Simulation configurations. All distances are given in cm from the transducer's centre, the transducer's geometric focus is at 13 cm, and all dimensions are in cm. The coordinate order is X, Y, Z.

	Small	Medium, superficial	Medium, deeper	Larger
Number of FZs	45	108	108	108
Treatment grid	3 × 5 × 3	6 × 6 × 3	6 × 6 × 3	6 × 6 × 3
Skin distance	11.4	11.4	11.0	11.0
Treatment plane distances	12.7, 13.0, 13.3	12.7, 13.0, 13.3	13.5, 13.8, 14.1	12.4, 13.0, 13.6
Constraint plane distances	11.7, 14.3	11.7, 14.3	12.5, 15.1	11.4, 14.6
Tumour dimensions	1.5 × 0.9 × 0.8	1.8 × 1.8 × 0.8	1.8 × 1.8 × 0.8	1.8 × 1.8 × 1.4
NT constraint (°C)	43	42, 42.5, 43	45	45
Distance from tumour to skin	1.3	1.3	2.5	1.4

Finally, for the more difficult to treat ‘medium deeper’ and ‘larger’ tumours a higher normal tissue temperature limit was used (45°C) for the constraint locations, which were deeper in the tissue and thus further from the skin's nerve endings. Post-treatment analysis of the normal tissue temperatures showed that this higher limit (imposed to speed treatments) was approximately equivalent to a 30 CEM dose limit. Even though the tumours studied were superficial, near-field thermal build-up activated the normal-tissue tissue constraints requiring interpulse cooling for approximately 90% of the paths studied.

Scanning paths

Due to the extremely large number of possible paths (approaching N) an exhaustive search was not practical. The 31 paths studied are a combination of ones that are exactly or roughly equivalent to those used previously [13], [14], [18], [20], [24], plus several that appeared promising due to potentially improved tumour heating and/or decreased normal tissue heating. The paths are broadly categorised as: axially stacked (AS), planar (PL), and 3-dimensional (3D). Within each category multiple patterns were studied, as summarised in Table III, while Figures 3 and 4 show example patterns.

Figure 3. The successive focal zone locations for: (Left) the ‘planar (BMF); XY raster’ path, in which the back plane is treated first with an X-Y raster scan, then the middle plane and then the front plane; (Right) the ‘axially stacked (MBF), XY raster scan’ in which the transducer is fixed at one transverse (X-Y) location and the focal zone is electronically scanned through the stack of three axial M,B and F locations. The transducer is then moved through the remaining transverse locations in an XY raster pattern, applying the same axial stack order at each location. The numbers specify the order in which the focal zone positions were scanned.

Figure 4. Axial view of four transverse, XY paths used to scan an axial stack: inner middle outer (IMO), concentric annuli (top left), knight (top right), adjacent square annuli (bottom left), and adjacent rectangular annuli (bottom right).

Table III.  Scanning patterns studied: F indicates the front position in an axial stack, or the front plane in a planar scan; M and B indicate the middle and back. MFB indicates that the M position (or plane) was treated first, etc.

Axially stacked (AS) scans		Planar (PL) scans		3-Dimensional (3D) scans
Axially stacking order	Transverse scan pattern of stacks	Order in which planes are treated	Scan pattern within each plane
All six permutations of F, M and B	XY raster	MFB, BFM, and BMF	XY raster	XYZ knight
MBF and MFB	Concentric annuli: all six permutations of I, M and O	LMR	XZ raster	Maximum distance first
MBF	Adjacent annuli: squares and rectangles	FBM	XY knight	Maximum distance last
MBF and MFB	XY knight			Random

For the axially stacked patterns, when the transducer is located at a given transverse (X-Y) location these patterns first treat all three axial focal zone positions in an ordered ‘stack’ (e.g. successively in a straight line along the z-axis in one of the six possible sequences involving the back (B), middle (M), and front (F) treatment depths). That axial stack sequence is then repeated at multiple, sequential, transverse locations which, taken together, cover the complete tumour. These transverse patterns include raster scans, which have been used by [5], [18], Annular scans, which have been applied clinically in hyperthermia treatments and HIFU experiments [20], [22], [35–37], and knight jumps (the L-shaped moves from a chess game), which were studied in an attempt to reduce normal tissue cooling times by using a transverse pattern that had a fixed distance between successive stacks that is close to the largest that can be consistently applied between all pairs of successive stacks throughout the complete treatment.

For the planar scans, axial stacking was not used and each of the three transverse planes (one at each axial depth; B, M and F) was scanned completely with a transverse pattern before moving to the next plane. In a very few cases three different planes were scanned (with left (L), middle (M), and right (R) signifying the planes as referenced to a plane containing the transducer's axis).

The 3D patterns (which, like the planar scans, did not use axial stacking) provided the opportunity to utilise larger distances between successive heating pulse locations, including (1) a maximum distance first (3DMaxFirst) path, which maximised the distance between successive FZ locations (i.e. the first two pulses occurred at the ends of the longest diagonal in the treated volume, etc.), resulting in numerous small separations at the end of the treatment, (2) a ‘3DMaxLast’ path (generated by reversing the order of the maximum first path), and (3) a ‘3D knight’ path, generated with 3D knight jumps between successive (x, y, z) focal zone locations, thus providing a relatively large fixed distance between all successive scan locations for the complete treatment. Finally, 100 randomly generated 3D paths were studied for the small superficial tumour.

Elemental reference cases

Four elemental reference cases that are composed of ‘building blocks’ from which more complex paths are constructed were developed and are presented in Appendix A. These cases serve as a reference base against which to compare the treatment time gains seen in the individual paths developed below.

Results

Results are presented in terms of the scan path treatment times. For most cases the ‘standard treatment conditions’ of a fixed applied transducer power level (giving a front treatment plane maximum power density of 0.66 × 10⁷ W/m³ for the superficial tumour), a perfusion magnitude 0.5 kg/m³s (typical for muscle tissue [1]), and a 6°C normal tissue constraint temperature were used. Due to the large number of paths and treatment conditions studied, and the lengthy computations required for the optimisation of the heating pulse and inter-pulse cooling period at every focal zone position along each path, computation times of over a week on a large, multicore computing node were sometimes required to find the treatment time for a single path. Thus, when it was clear that a scan path being studied had treatment times far from optimal (e.g. the optimal treatment time was already significantly exceeded after only half of the tumour was scanned), the treatment times were extrapolated to an estimated final time, and the results indicated by asterisk.

Effects of scan path

To illustrate the large effect of scan path, Figure 5 presents the treatment times for fourteen scan paths with the standard treatment conditions for the medium superficial tumour. Clearly the axially stacked X-Y raster scans with the heating pulses starting in the middle of each stack (the AS (MFB), XY Ra and AS (MBF), XY Ra scans) are the fastest, with treatment times that are a factor of 5.5 better (6.7 vs. 36.7 min) than the worst case simulated, i.e. the PL (BMF), XY raster scan that successively treats the back, middle and front planes with XY raster scans. (They are also a factor of 4.5 better than the case of all N pulses applied independently, i.e. than the ‘independent pulse’ case given in the appendix (6.7 vs. 30 min).)

Figure 5. Treatment time vs. scanning path type for the medium superficial tumour (108 focal zone positions) and the standard treatment conditions. Also shown are the percentages of the treatment time spent cooling the normal tissue. When the cooling time percentages are nearly identical, only one value (located between the paths with similar values) is shown.

The small superficial tumour was also studied with most of the paths in Figure 5, and the resulting relative path rankings were unchanged from those in Figure 5. The small superficial tumour was also treated with 100 randomly generated paths, with their treatment times all lying roughly midway in the relative path rankings.

Effects of applied power and perfusion

Figures 6 and 2 show the effects of applied power and perfusion magnitude on the treatment times and the relative rankings of the scan paths. Results are shown for one of the best (AS (MBF), XY Ra), worst (PL (BMF), XY Ra) and middle (3DKn) scans from Figure 5. (The other scans from Figure 5 that were studied at different powers and perfusions followed the trends in Figures 6 and 2, and are thus not shown.) The trends of decreasing treatment times with increased power in Figure 6 are the same as shown previously [17] and illustrate why the standard power density level of 0.66 × 10⁷ W/m³ was chosen – since larger powers yield only marginal treatment time gains, while lower powers greatly increase the treatment times. Figure 2 shows that for the best cases (e.g. AS (MBF), XY Ra scan) for which no normal tissue constraints are activated even at low perfusions, increasing the perfusion magnitude provides no additional benefit in the normal tissue and makes it harder to heat the tumour, thus increasing the treatment times. Conversely, for the slowest scan patterns (e.g. the PL (BMF), XY raster scan), which already have long treatment times at low perfusions (since they activate the normal tissue constraints so frequently), increasing the perfusion cools the normal tissue, and thus reduces the treatment times. A small number of scans were also performed in a tumour model with heterogeneous perfusions, where the tumour perfusion level was twice the normal tissue perfusion. These scans used a subset of the paths described above, and in all cases the ranking of the paths was unchanged from those in Figure 5.

Figure 6. Treatment time vs. power density for three scan paths treating the small superficial tumour with 45 positions. The standard perfusion and constraint values were used but the total applied power was varied. Also shown are the percentages of the treatment time spent cooling the normal tissue. The lines are for the AS (MBF), XY (Ra) (squares); 3D Kn (circles) and PL (BMF), XY Ra (triangles) scan patterns.

To illustrate why the axially stacked scans have the shortest treatment times and other scans have much longer treatment times, Figure 7 shows the heating and cooling times at each successive position for the AS (MBF), XY Ra and the PL (BMF), XY raster paths for the medium superficial tumour. In Figure 7A, the benefits of axial stacking can be seen by noting that the axially stacked case has one long (initial) heating pulse (at M) followed by two much shorter pulses (at B and F) within each axial stack. The short pulses arise since axial stacking allows for sequential pulses to build up the tumour temperatures (before they decay) in the direction in which the most SAR overlap can be utilised, i.e. axially stacking takes advantage of the elongated nature of the focal zone. In contrast, the PL (BMF); XY raster path has nearly constant long heating times at each focal zone position in every plane (Figure 7A) since it does not take advantage of such SAR overlap. In addition, it treats the back plane first, thus spending the most time (the first pulses) at locations for which the ratio of SAR in the focal zone to that in the normal tissues is the lowest, thus activating the inter-pulse cooling much more frequently. Figure 7B shows how these longer heating times for this planar scan result in much longer inter-pulse cooling periods. A second trend in Figure 7A is the steady reduction in heating times as the treatment progresses for both scan patterns, corresponding to increased thermal build-up in the tumour, an effect present for both scan patterns, but to different degrees.

Figure 7. Pulse heating and inter-pulse cooling times (s) during treatment of the medium, superficial tumour under the standard treatment conditions: (A) The heating times at each sequential focal zone position for both the axially stacked (MBF), XY raster path (squares) and the planar (BMF), XY raster path (x). (B) The cooling times at each sequential focal zone position for the planar (BMF), XY raster path. The cooling times for the axially stacked (MBF), XY raster are not shown because they are zero for almost all positions in the treatment and are only a few seconds for the remaining positions.

The bio-physical mechanism by which the treatment time benefits are gained through axially stacking is illustrated in Figure 8, which shows the temperature vs. time curves in the centre of the first focal zone treated (M) for a heating pulse applied to the middle of the larger tumour. It also shows the corresponding curves for the voxels that are immediately proximal (F), distal (B) and transversely adjacent to the focal zone. Here one can see that if after heating the middle position, the focal zone is rapidly moved to the front position for the second heating pulse, that the temperature elevations resulting from the second pulse will be superimposed on the already significantly elevated tumour temperatures at the front position, and thus the front pulse will need a much shorter heating period. Since this superposition occurs at high temperatures it also takes maximum advantage of the gains to be obtained from the non-linear temperature–thermal dose relationship, thus further reducing the needed heating times. Similar but smaller heating time gains would occur if one immediately moved to the back position when applying the second pulse. However, if one moved off-axis by a transverse motion, then not only would there only be a small amount of temperature gain due to superposition, but also the temperatures in the front and back axial locations would be allowed to decay, and thus when they are heated at a later time, long heating periods will be required to bring them back up from basal conditions to the desired target temperature. Such longer heating times will increase the energy deposited in the normal tissues, thus necessitating additional inter-pulse cooling times (as can be seen in Figure 7B for the PL (BMF), XY Ra scan), thus further lengthening the treatment times. Finally, the largest time gains arise from heating the middle position first since this generally maximises the total power (and energy) deposition in the tumour during the first (and longest) axial heating pulse, thus producing a larger ‘pre-heated’ temperature base at the other two axial focal zone locations.

Figure 8. Temperature vs. time for the middle, front, back, and transversely adjacent (to either side of the middle position in the plane perpendicular to the transducer axis) locations with respect to the FZ location when the beam is focused on the middle (M) position. Results are shown for the first heating pulse for the larger tumour under standard treatment conditions.

Transverse path patterns of axially stacked scans

Since the axially stacked scans starting in the middle position gave the shortest treatment times for a wide range of conditions when studied with the XY transverse raster path treating the superficial tumours, additional transverse patterns were studied with the AS (MBF) axial stacking sequence to see if further treatment time gains could be realised. These treatments were performed on the medium superficial tumour with standard treatment conditions. In addition, since a primary safety concern in treatments is thermal build-up in normal tissues, in order to study more difficult treatment conditions the more stringent normal tissue temperature limit of 5°C was also studied. Figure 9 shows (ranked in order of the treatment times for the 5°C limit) the resulting treatment times for the following axially stacked AS (MBF) transverse scan patterns: AS (MBF) XY raster; the AS (MBF) adjacent annuli; the six permutations of the AS (MBF) O, M and I annuli for the concentric annuli scans; and the AS (MBF) XY knight jumps. In the 6°C constraint case there was essentially no inter-pulse cooling needed for any scan patterns (the adjacent square annuli case used 3% of the time for cooling; all others were zero or <1%), and thus normal tissue considerations were not significant. For the 6°C constraint case, all of the transverse scan patterns had very similar treatment times, and as expected since tumour heating was the dominant consideration, the fastest scan paths were generally those with the most adjacency among focal zone stacks (the adjacent square annuli had the shortest treatment time – 418 s) while the longest treatment times had the least adjacency (the AS (MBF) XY knight pattern – 469 s). However, once the normal tissue constraint conditions become very active (i.e. the 5°C limit), then path choice became significant, and the ordering of the paths essentially reversed – with the path with the least focal zone adjacency (knight) becoming the fastest and the path with the most adjacency (squares) becoming the slowest. This change clearly illustrates the problem dependency of the optimal scan paths when normal tissue constraints are very active. It should also be noted for this 5°C normal tissue constraint that while all transverse patterns (except the adjacent square annuli) still gain from tumour temperature superposition (their treatment times are less than the 1860 s independent pulses treatment time – see Appendix A), some paths lose part of that advantage due to activation of the inter-pulse cooling periods. This is evidenced by the monotonically increasing cooling times as the treatment time increases. Finally, for this stricter normal tissue constraint limit, the shortest treatment times are for the AS (MBF) transverse patterns of concentric annuli starting at the tumour centre (IMO and IOM); the XY knight jumps; and the XY raster scan. Note that at the 5°C limit these four patterns have treatment times that are less than the treatment times for all of the Figure 4 cases at the 6°C limit (except the AS (MBF) and AS (MFB) raster scans), so no other patterns from Figure 4 were worth running since their treatment time would only increase if run with the 5°C limit. These four transverse patterns are an interesting combination of quite different transverse scans which all balance the goal of optimising the gains from increased adjacency inside the tumour, with avoiding overheating the normal tissues due to SAR overlap in the near field.

Figure 9. Treatment time vs. transverse scan path for axially stacked (MBF) stacks treating the medium, superficial tumour for 6°C and 5°C normal tissue temperature limits under the standard treatment conditions. The scans are rank ordered based on their times for the 5°C constraint case. The path abbreviations are: Kn for the AS (MBF) XY knight jumps; inner I, middle M, outer O for the AS (MBF) concentric annuli; Ra for the AS (MBF) XY raster; Rec for the AS (MBF) adjacent rectangular annuli; and Sq for the adjacent square annuli. Paths where the treatment times were extrapolated have an asterisk. The percentages of treatment time spent cooling the normal tissue for each path with the 5°C limit are shown. The percentages for the 6°C limit are not shown because they were all less than 3% and were not a significant factor in determining total treatment time.

Deeper and larger tumours: SAR overlap

To study the effects of path selection in other, more difficult to treat situations, treatments were also performed on the medium tumour moved to a deeper location and on an axially larger tumour. Since the 6°C constraint limit resulted in large amounts of normal tissue thermal build-up and thus impractically long treatment times for these tumours, the NT constraint temperature was set at 8°C for these studies. (This 8°C temperature constraint is approximately equivalent to a thermal dose limit of 30CEM43°C since in retrospectively examining the results for these cases, a NT dose of greater than 30 CEM was observed for at most one NT voxel in every treatment.) Figure 10A and B show the resulting treatment times for the medium deeper tumour and the larger tumour for the four fastest transverse scanning patterns from Figure 9, i.e. the AS (MFB): XY knight jumps; IOM and IMO concentric annuli: and XY raster. The results for these medium and deeper tumours differentiate these four transverse paths even further, showing that the AS (MFB) XY raster and AS (MFB) IMO concentric annuli paths are the most consistently fast scans for all cases examined. Finally, scan path treatments were also performed on both this medium deeper and this larger tumour using all other five permutations of (M, F, B) for the axially stacked XY raster path. In all cases, the treatment times were slower than for the corresponding MFB scans, thus reinforcing the optimal nature of the MFB axially stacking sequence.

Figure 10. Treatment time (s) vs. transverse path for (A) the medium deep (top) and (B) larger tumour (bottom) for the standard treatment conditions and a NT temperature limit of 8°C.The independent single pulse treatment times for these tumours are 3852 and 2556 s, respectively (Appendix). The percentage of treatment time spent cooling is shown above each path (or between points in cases where the percentages were essentially identical).

In order to better understand the effects of transverse path on treatment times, the heating times for each axial stack during the treatment (i.e. for all 36 stacks) for each of the five paths of Figure 10 are shown in Figure 11 (A–E), where stacks which were preceded by a previous stack requiring a long cooling time (>30 s) are marked with a tilde. The heating time of the initial stack for all paths was 39 s and decreased for later stacks in all paths as focal zone adjacency came into play and thermal build-up in the tumour accumulated. Also note the ‘spikes’ in the heating time that occurred either after a long cooling period or because of a large change in the path position that reduced the adjacency between successive focal zone stacks.

Figure 11. The stack heating time vs. stack number for the paths in Figure 10. The graphs, in order from shortest to longest treatment time, are for the AS (MFB); (A) raster, (B) IMO concentric annuli, (C) IOM concentric annuli, (D) OMI concentric annuli, and (E) knight. Cooling times of greater than 30 s in a preceding stack are marked with a tilde.

Discussion

The results presented confirm this study's primary hypothesis, that large treatment time reductions can be realised by judicious scan path selection, at least for the tumours and conditions studied herein. They also identify and quantify the effects of the underlying thermo-physical principles that make scan paths ‘fast’, thus providing guidance for future studies with more complicated conditions (e.g. changing blood flow and attenuation coefficients). The most important effect is ‘axially stacking’ of successive focal zones in the tumour, which arises from the combination of axial SAR overlap and temperature superposition due to the focal zone's elongated shape. More specifically, the MFB axially stacked paths studied are generally fastest since they apply their longest heating pulse (M, the initial pulse of each stack) when the focal zone is depositing the most power inside the tumour, and then move to the tumour location that was preheated the most (F) by the initial pulse in order to take maximum advantage of temperature superposition in the tumour. This effect is shown to hold for a wide range of perfusion and transducer power levels, and multiple tumour conditions that activate the normal tissue constraints. In summary, while most previous scan-related studies have concentrated on the problem of avoiding thermal build-up in the normal tissues by having successive focal zone locations spread far apart, the current study shows how in many cases the completely opposite approach of placing successive pulses close together (in an axial stack and with maximum stack adjacency) can shorten treatment times considerably.

The second most important thermo, physical factor, is adjacency of successive transverse stacks, an effect important in treatment paths used in other research [20]. The effect of adjacency can best be seen by comparing the heating times for the AS (MFB) XY knight path (Figure 11E) which has little adjacency, to the other AS (MFB) paths of Figure 11 which all have considerable adjacency. This reduced adjacency leads to much longer average stack heating times for the AS (MFB) XY knight pattern at the start of the treatment, and thus creates longer treatment times both directly because of those longer heating periods and indirectly because the normal tissue is also being heated for longer times. This results in higher normal tissue temperatures and an increased need for inter-pulse cooling following later heating pulses.

The current study concentrates on superficial tumours with fixed properties, and extends previous treatment time optimisation studies by showing that the time gains from path selection add to those obtained by optimising the individual pulse heating and cooling times [17]. The results shed light on the best approach for treating the presented tumour configurations, and provide a basis for further studies with more complicated conditions. An interesting topic for such future studies is temporally changing tissue properties, both changing tumour and normal tissue perfusion, and increased tumour attenuation. Attenuation changes have been observed both experimentally [38–40] and when matching simulation data to previous experimental results [41].The scan-related convention that has been followed by several previous investigators is to treat from the back to the front in a tumour (e.g. using a planar BMF, XY raster scan), with the idea being that if the front of the tumour is treated first that any resulting increased attenuation would block the ultrasound from reaching the deeper positions. The current study indicates that if indeed such a PL (BMF), XYraster scan strategy is adopted, that treatment times will be very long. Thus, it is important to study this question in much more depth in the future to quantify the treatment time trade-offs involved. While several studies have begun investigating tissue attenuation changes during thermal therapy treatments, much remains to be determined experimentally in vivo before definitive conclusions can be drawn about the role of attenuation change in path optimisation. In particular, the speed with which attenuation changes occur during an in vivo HIFU treatment is not precisely known (i.e. whether such changes occur within seconds of ablation or on a scale of minutes, etc.). Previous ex vivo studies have shown that the level and rate of tissue attenuation changes due to thermal ablation depend on the temperature to which the tissue is heated [42–46], the frequency of the ultrasound beam [40], [43], [44], [46], and both the rate of dose delivery and the maximum dose delivered [42], with tissue attenuation changing more after it has been heated and then cooled, vs. when it is kept in a heated state [45]. Further study is needed to determine what affect tissue attenuation changes have on the optimal path choice, and whether an axial treatment which starts in the middle of the tumour could be sufficiently ‘fast’ to treat an axial stack before the tumour's attenuation changes significantly. However, the evidence presented in Figure 4 shows that, if further research shows that tissue attenuation changes are both fast enough and large enough in magnitude to make a back-to-front focal zone path necessary for some treatments, starting with a back treatment plane and axially stacking the focal zones is still faster than a transverse raster scan path that does not use axial stacking. Similar timing and magnitude considerations are present for possible changes in tumour and normal tissue perfusion changes during treatments.

Overdosing (tumour positions above 240CEM) was also examined but was found, both through post-treatment analysis and exploratory simulations on selected paths, to not be a significant enough effect to change the rank ordering of the paths presented for the treatment conditions studied herein – although reductions in treatment times would be possible on all paths through the reduction of overdosing. Avoiding large amounts of overdosing in practice could only be done by anticipating all future doses to be delivered to a given position, a difficult task that would rely on accurate modelling. The problems of overdosing and attenuation changes are interrelated; that is, attenuation changes will likely be very significant in treatments that go to very high temperatures and thus deliver very high tumour doses (with concomitant overdosing for treatment efficacy ‘insurance’).

One of the main results from the simulations conducted in this paper was an examination of the physics present in treatment paths that have significant levels of axial stacking. Even though the temperature increase due to this stacking may vary for a particular treatment, the results of these simulations lend credence to the belief that the underlying physics of axial stacking and stack adjacency will be a useful tool for treatment acceleration in many, and more complicated, clinical situations. In summary, these results provide useful guidelines for reducing treatment times through path selection in many tumours, and a reference base for future studies of more complex situations, e.g. for treating larger, deeper tumours with heterogeneous perfusions, irregular geometries and temporally changing properties. Given the very large treatment time gains realised in the current study, it is reasonable to expect that by applying the basic bio-thermal approaches identified in this study, those gains will translate, in some proportion, to other cases. Given the immense scope of an exhaustive search over all tumour anatomies and physiologies, transducer configurations, focal zone sizes and shapes and spacing, and N! paths for any tumour, it is unlikely that a truly universally optimal scan path will ever be found, nor be needed since it is likely that some ‘close to optimal’ patterns will give practical treatment times.

Conclusion

Optimisation of the focal zone scan path in HIFU treatments yields significant treatment time savings. The optimised treatment paths use axially stacked (along the transducer axis) focal zones, with each such stack starting with the focal zone in the tumour centre. These time savings occur primarily because the optimised axial focal zone stacks take advantage of SAR overlap and temperature superposition in the tumour in successive axial focal zone heating pulses. The tumour is generally best treated by multiple such stacks scanned through the tumour in an optimised transverse pattern that maximises the amount of ‘stacking adjacency’, with the best path being dependent on the normal tissue constraints. These time gains are present over a range of transducer powers and tissue perfusion levels, and for several different tumour configurations. Such treatment time reductions will be very important clinically, thus motivating further path studies and searches for even better thermophysics-based heuristic guidelines for a wider range of expected clinical conditions.

Acknowledgements

We greatly appreciate the help of Dennis Parker, Urvi Vyas, and Doug Christensen, and use of the UCAIR facilities for this work. An allocation of computer time from the Center for High Performance Computing at the University of Utah is gratefully acknowledged.

Declaration of interest: This work was partially supported by grants from NIH (R01-CA134599), Siemens Medical Solutions, the Focused Ultrasound Foundation, a University of Utah Synergy Grant and the Ben B. and Iris M. Margolis Foundation. The authors alone are responsible for the content and writing of the paper.

Appendix A

Four cases that represent the most elemental ‘building blocks’ used to construct all 31 of the scans studied are presented to investigate the thermo-physical mechanisms underlying the scans’ differences in treatment times. These idealised cases include ‘independently optimised’: (1) pulses; (2) single axial stacks; (3) pairs of adjacent axial stacks; and (4) corner triplets of adjacent axial stacks. All reference stacks used MFB axially stacking since this stacking order always gives the shortest treatment times of all six stacking orders. All four of these cases are ‘independent’ in the sense that for the independent pulses, each pulse is applied under the idealised conditions that it is not affected by, i.e. is independent of, all other pulses. Such pulse ‘independence’ could arise if there was a high blood perfusion everywhere, and the sequential heating pulses were widely separated in space so they did not affect each other. Similarly for independent stacks, each stack is applied so that it is not affected by any other stacks, and similarly for the pairs and triplets of stacks. Each case is ‘optimised’ since the pulse heating and inter-pulse cooling periods are optimised separately for each case. That is, the independent single pulse case uses one heating pulse period that is optimised for the independent heating of one location, the independently optimised single axial stack has three heating periods and two cooling periods that are optimised for the treatment of the three axial positions within one independent stack, etc.

Independently optimised single pulses

If each heating pulse was independent of all other pulses then it would begin its heating pulse with the tumour at basal conditions and treat its focal zone location to the needed dose, but would have no effect on any other location. This means that there would be no treatment time gains at that location from 3D SAR overlap among successive pulses, nor from long term thermal build-up in the tumour, nor from the delivery of prior or future dose. Similarly, there would be no treatment time delays from inter-pulse cooling since there would be no thermal build-up in the normal tissues (i.e. Δt_Cool,_n = 0 for all n) other than from the single pulse, which will not require inter-pulse cooling by definition. The resulting treatment time is thus N times the independent pulse's optimised heating time, which is the same for all pulses due to their independence.

Independently optimised single axial stacks

Next, a series of ‘independent’ MFB stacks, whose three pulse heating and two inter-pulse cooling times were optimised within that stack, was studied. This approach takes advantage of SAR overlap and temperature superposition within the tumour to reduce the second and third pulse heating times (when compared to the first pulse). In the tumours studied it also resulted in some normal tissue thermal build-up due to the three heating pulses, but never enough to reach the constraint temperature, so normal tissue thermal build-up was not a consideration (Δt_Cool,_n = 0 for all n). There is no interaction with other stacks, and thus there is no long-term thermal build-up in the tumour or in the normal tissues. The difference between this case's treatment times and the independent single pulse case's treatment times is thus a measure of the heating time gained by using axial stacking. The resulting total treatment time is N/3 times the optimised independent axial stack heating time, which is the same for all stacks due to their independence.

Independently optimised pairs of adjacent axial stacks

Next, since in most transverse scan patterns studied, successive stacks are placed adjacent to each other (in x or y, but not diagonally), those pairs gain from their adjacency (the second stack will have lower heating times since it will have been preheated by SAR overlap heating and transverse conduction from the first stack). If such pairs are independent, they will have interactions between themselves, but not with any other stacks. Comparing the treatment time for the optimised independent pairs of adjacent stacks to that for independent single axial stacks gives a measure of the time gains obtained solely from the ‘two stack adjacency’. The resulting treatment time is N/6 times the heating time needed to heat one such pair of stacks, which is the same for all such pairs due to their independence. Again, for this case, the normal tissue constraint temperature was never reached during the six pulses used by a pair of stacks, and thus Δt_Cool,_n = 0 for all n.

Independently optimised corner triplets of axial stacks

Finally, the sequence of three adjacent stacks that form corners (e.g. stacks 1, 2 and 3 at the beginning of the AS (MFB), IMO pattern of Figure 3) not only gain from the immediate adjacency of each pair of successive stacks, but also from the diagonal heating between stacks 1 and 3. When compared to the paired stack times they give a measure of the time gains realised solely from adding a third stack to form a ‘corner triplet’. Their treatment time is N/9 times the time required for one such triplet. Again, in these cases there is little thermal build-up due to the nine heating pulses and thus normal tissue constraints were never activated (Δt_Cool,_n = 0 for all n).

The ‘reference case’ results for the standard treatment conditions are presented (Figure 12) to both show the gains attainable from the different ‘building blocks’ scanning approaches, and as a baseline for comparison with the results for the 31 paths studied. Also shown are the results from those paths (out of all of the 31 paths studied) with the shortest treatment times. Note that those ‘fastest’ paths have the additional benefits (not present in any of the ‘independent’ reference cases) of both (1) thermal build-up in the tumour, a benefit that results in further shortening of the treatment times, and (2) accounting for ‘prior dose’ to reduce treatment times. These results clearly illustrate the gains achievable through optimising the use of the thermal physics (primarily axial stacking and secondarily stack adjacency) during treatment path selection.

Figure 12. Treatment time vs. reference case for the medium superficial, medium deeper, and larger tumour treated with ‘independent’ pulses and stacks under the standard treatment conditions. Also shown are the times for the fastest scans obtained in this study for the same treatment conditions. The fastest case for each of the tumours was; AS (MFB) XY raster for both the medium superficial tumour and for the larger tumour, and AS, MFB IMO concentric annuli for the medium deeper tumour.