Volumetric analyses of ablation dimensions in microwave ablation for colorectal liver metastases

Abstract Background In thermal ablation of malignant liver tumors, ablation dimensions remain poorly predictable. This study aimed to investigate factors influencing volumetric ablation dimensions in patients treated with stereotactic microwave ablation (SMWA) for colorectal liver metastases (CRLM). Methods Ablation volumes from CRLM ≤3 cm treated with SMWA within a prospective European multicentre trial were segmented. Correlations between applied ablation energies and resulting effective ablation volumes (EAV) and ablation volume irregularities (AVI) were investigated. A novel measure for AVI, including minimum enclosing and maximum inscribed ellipsoid ablation volumes, and a surrogate parameter for the expansion of ablation energy (EAV per applied energy), was introduced. Potential influences of tumor and patient-specific factors on EAV per applied energy and AVI were analyzed using multivariable mixed-effects models. Results A total of 116 ablations from 71 patients were included for analyses. Correlations of EAV or AVI and ablation energy were weak to moderate, with a maximum of 25% of the variability in EAV and 13% in AVI explained by the applied ablation energy. On multivariable analysis, ablation expansion (EAV per applied ablation energy) was influenced mainly by the tumor radius (B = −0.03, [CI –0.04, –0.007]). AVI was significantly larger with higher applied ablation energies (B = 0.002 [CI 0.0007, 0.002]]); liver steatosis, KRAS mutation, subcapsular location or proximity to major blood vessels had no influence. Conclusions This study confirmed that factors beyond the applied ablation energy might affect volumetric ablation dimensions, resulting in poor predictability. Further clinical trials including tissue sampling are needed to relate physical tissue properties to ablation expansion.


Introduction
Percutaneous microwave ablation (MWA) is an established treatment technique for malignant liver tumors, including colorectal cancer liver metastases (CRLM) [1]. As a minimally invasive, parenchyma preserving procedure it is associated with low complication rates, short hospital length of stay and enhanced patient recovery [2]. In selected patients, MWA is associated with survival comparable to that after resection [3] and can be performed repeatedly in case of recurrent disease [4,5]. Next to standard image-guidance modalities such as ultrasound, computed tomography (CT) and magnetic resonance, image fusion and stereotactic guidance have been proposed to enhance tumor visibility, targeting accuracy and avoid high radiation exposure [6,7]. This was shown to augment treatment accuracy, safety and targeting success, leading to enhanced efficiency and efficacy when using stereotactic guidance [7][8][9].
Next to targeting accuracy, the effective creation of ablation volumes leading to complete tumor coverage with adequate ablation margins is crucial for treatment success [10][11][12]. The generation of ablation volumes is a function of the delivered energy, controlled via power and duration settings on the ablation device. Standardized energy settings for different ablation sizes in 2 D are provided in the ablation device manufacturers' instructions for use (IFU). These are usually derived experimentally from ex-vivo (non-perfused) or in-vivo healthy animal livers [13]. Concerns regarding a potential lack of correspondence between these predictions and ablation dimensions actually obtained in clinical practice have previously been raised [14][15][16]. Clinical experience suggests that physical tumor and liver parenchyma characteristics related to tissue perfusion (e.g., perivascular location or tumor vascularity) and liver consistency (e.g., fatty liver disease) might influence the expansion of ablation energy and resulting ablation dimensions, beyond what can be predicted based on energy settings [13,16]. Furthermore, irregularities and inconsistencies of created ablation shapes make predictions of ablation dimensions and thus safety and efficacy of ablation treatment difficult.
The aim of this work was to acquire further knowledge on predictors of treatment success in microwave ablation for CRLM, by investigating i) correlations of predicted and effective ablation volumes, ii) ablation volume irregularities, and iii) factors influencing ablation dimensions beyond energy settings, in a patient cohort extracted from a prospective European multicentre study.

Patient population
Clinical data from a prospective European multicentre cohort trial (Microwave Ablation Versus Resection for Resectable Colorectal Liver Metastases, MAVERRIC, clinicaltrials.gov: NCT02642185) were analyzed. In this trial, 98 patients with a maximum of 5 CRLM lesions of 30 mm diameter, amenable to both CT-guided stereotactic microwave ablation (SMWA) and surgical resection as assessed in local tumor boards, were prospectively included for SMWA treatment in three tertiary HPB Centers (Bern, Switzerland; Groningen, Netherlands; Stockholm, Sweden). While the main study endpoint was overall survival after SMWA as opposed to surgical resection, secondary endpoints included technical analyses of resulting ablation volumes. All clinical and procedure-related data were stored on a secured RedCap database. All imaging data were retrieved from the SMWA systems and stored on secure file storage. Ethical approval was obtained by all respective regional committees.
CT-guided SMWA was performed using the commercially available navigation system CAS-One IR (CAScination AG, Bern, Switzerland) and the Needle Positioning System (NPS; DEMCON Advanced Mechatronics) [8] for enhanced planning of ablation probe trajectories, stereotactic positioning of ablation probes, quantitative validation of ablation probe positions and qualitative validation of ablation success. Detailed setup and procedural workflows of both systems have been described in detail previously [8,[17][18][19]. For patients in whom kidney function allowed a second contrast injection, CT scans for qualitative validation of ablation success were performed immediately after SMWA treatment. Ablation systems routinely used in each center were applied, including Acculis pMTA (Angiodynamics, Latham, NY USA), Solero (Angiodynamics, Latham, NY USA), Amica (HS Medical Inc, Rome, Italy) and Covidien EmprintTM (Covidien/ Medtronic AG, USA, Minneapolis).

Data extraction and analysis
Clinicopathological characteristics (e.g., KRAS mutational status) and ablation data (energy settings) were extracted from the study database, and corresponding imaging data was obtained from each center. Retrospective analysis of ablation dimensions generated after SMWA of CRLM was performed using the immediate post-ablation contrast-enhanced (CE)-CT validation scans. In cases of immediate re-ablation in the same treatment session due to initial incomplete ablation, only the first ablation validation scan was used for ablation volume analysis. We excluded cases in which i) no immediate CE-CT validation scan was available (e.g., due to restricted tolerance of contrast administration), ii) multiple parallel ablation probes were applied to create larger ablation volumes, and iii) merging ablation volumes resulted due to multiple tumors ablated in proximity, resulting in the difficult distinction of individual ablation volume boundaries. In included cases, the ablation volumes were manually segmented using semi-automatic segmentation tools available in the software Amira (ThermoFisher Scientific, USA), and stored as binary masks. All segmentations were performed by a trained technical expert and verified by an experienced interventional radiologist from the respective clinical institutions.

Effective and predicted ablation volumes
Obtained effective ablation volumes (EAV) resulting after SMWA were correlated to predicted ablation volumes (PAV) and applied ablation energies, stratified according to the type of applied ablation device. Subgroups with 10 EAV samples available were excluded from the analysis. Ablation energy was calculated from the applied power and time settings (energy ¼ power Ã time), as set on the ablation device by the treating physician. PAV was extracted from the ablation device manufacturer's guides. For applied ablation energies for which no corresponding PAV dimensions were available from the manufacturer's IFU's, but which were situated within the minimum and maximum energy range provided from the IFU's, PAV was interpolated using linear regression (PAV ¼ b 0 þ energy Â b 1 ). Additionally, the percentage of tumor volume in the ablation volume was calculated as described by Heerink et al. [14].

Ablation volume irregularity
We hypothesize that ablation volume contours are of importance for obtaining complete ablation. Therefore, we introduced ablation volume irregularity (AVI) as a potential risk factor for incomplete tumor coverage by the ablation zone. In the present study, correlations between AVI and applied ablation energies and factors influencing AVI were investigated. AVI was defined as the deviation of the ablation volume from a perfectly ellipsoid shape. To this end, the minimum enclosing ellipsoid volume (MEV) and maximum inscribed ellipsoid volumes (MIV) of obtained ablation volumes were calculated as illustrated in (Figure 1). AVI was calculated as the proportion of the volume between MIV and MEV compared to the total MEV (see equation below). AVI ¼ 1 represents a highly irregular ablation consisting of only spikes and AVI ¼ 0 represents a perfectly ellipsoid shape with MEV ¼ MIV.

AVI ¼ MEVÀMIV MEV
This definition was chosen to account for the fact that MEV-MIV increases with larger ablation volumes and to create an outcome parameter with a fixed range. A convex optimization algorithm as defined by P. Boyd and L. Vandenberghe [20] was applied for computation of MIV and MEV, and the CVXPY library [21] was used for implementation.

Factors influencing ablation dimensions
Factors potentially influencing resulting ablation dimensions were investigated using linear regression analysis. Separate linear regression models were created, using EAV per kJ energy, MEV/kJ, MIV/kJ and AVI as outcome variables. EAV/ energy ( ml kJ ) was chosen as a parameter defining the expansion of ablation energy across the tissue, acting as a possible surrogate variable for tissue properties affecting heat generation and conduction (e.g., heat capacity, microwave absorption). MEV/energy and MIV/energy ( ml kJ ) were chosen as parameters defining the maximum and minimum expansion of ablation energy. All factors available from the study database and thought to potentially influence the spread of ablation energy and thus ablation volume dimensions (EAV/ kJ and AVI) beyond the applied ablation energy, were included in the regression models. An overview of factors thought to potentially influence ablation expansion and thus ablation dimensions are provided in Figure 2. Their description and rationale for inclusion are provided in the following sections.
Device-specific factors. Due to different wavelengths, antenna designs and other technical factors, different ablation devices yield varying ablation volume dimensions, depending on the applied energy [22]. The type of applied ablation device was therefore included as a random effect in both models. The applied ablation energy was included as an independent variable in the model investigating AVI, and as part of the outcome parameter in the EAV/energy model.
Tumor specific factors. Tissue necrosis due to MWA has derived from direct heat production and absorption of microwaves in the vicinity of the ablation probe and by heat conduction further away from the probe [23]. Water and lipid concentrations in cancerous tissue are different from those in the healthy liver parenchyma. This results in variation in physical properties affecting direct heat creation and heat conduction [24]. To capture different effects of microwaves onto cancerous versus non-cancerous liver tissue, tumor radius defined as the distance between tumor center and non-cancerous parenchyma was included as an independent variable. A subcapsular tumor location ( 5 mm distance from the liver capsule) was included as a physical factor limiting ablation expansion and promoting tissue shrinkage. Additionally, differences in vascularity in the periphery of the liver compared to more centrally located parenchyma could influence ablation expansion. Proximity (distance 5 mm) to intrahepatic vessels (!3 mm in diameter) was entered as a variable in the model to account for the potential heat sink effect known to influence ablation dimensions [25].  Patient/liver-specific factors. Liver steatosis was described to affect heat conduction properties during ablation [26]. As a surrogate measure for hepatic steatosis, we compared Hounsfield units (HU) of the liver and spleen on the intraprocedural native CT scan, with hepatic steatosis defined as HU Liver -HU Spleen < À10 [27]. To assess if Kirsten Rat Sarcoma (KRAS) tumor mutation, previously described as a significant influencing factor for ablation site recurrence [28], affects ablation expansion, KRAS mutation status was included as an independent variable.

Statistics
Outcome parameters were reported as medians and interquartile ranges (IQR). Linear correlation analysis was applied to investigate relationships between applied energy and EAV, EAV and PAV and AVI. Correlation <0.4 was considered weak, <0.7 moderate and >0.7 strong [29]. Univariable and multivariable linear regression analyses were performed to investigate factors potentially influencing EAV/kJ, MEV/kJ, MIV/kJ and AVI, with all independent variables included in the multivariable model. Linear mixed-effects modeling was applied to account for clustering of data for patient-specific factors (hepatic steatosis surrogate measure, KRAS mutation) and device-specific factors (device type), with 'patient' and 'ablation device' included as random effects. The effect size of all tested variables was reported as linear regression coefficients with 95% Wald confidence intervals (CI). Statistical significance was set at p < 0.05 level. Statistical analyses were performed using R (R Core Team, 2019) and RStudio (RStudio Inc., USA) and graphs were generated using ggplot2. The R package lme4 was used for mixed-effects modeling. The full reproducible statistical code and data for R are available in the supplementary material 'statistical_analysis.pdf' and 'data.csv'.

Results
From a total of 168 CRLM treated with SMWA in 98 patients included in the MAVERRIC trial, 116 tumors were suitable for retrospective volumetric analyses of ablation dimensions and 105 tumors for correlation analysis of EAV to PAV. PAV was interpolated in 74 of the 105 cases included for correlation analyses (70%). Figure 3 illustrates the reasons for the exclusion of treated tumors and the number of included tumors per applied ablation device. In the 116 lesions included for volumetric analyses, 17 (15%) arose in steatotic livers, 49 (42%) had KRAS mutation of the primary tumor, 53 (46%) were located in a subcapsular position and 18 (16%) in proximity to major blood vessels.
Median ablation energies in kilojoules (kJ), tumor radius in millimeters (mm), EAV, PAV, MEV, MIV in milliliters (mL) and AVI as ratios are summarized per ablation device and overall (n ¼ 116) in Table 1.

Effective and predicted ablation volumes
Correlation analysis between EAV and applied energy per ablation device are illustrated in Figure 4. Overall, the correlation between applied energy and EAV was moderate (r ¼ 0.4, R 2 ¼ 0.16, p < 0.001) for the Acculis device, moderate (r ¼ 0.5, R 2 ¼ 0.25, p ¼ 0.01) for the Amica device, and weak (r ¼ 0.38, R 2 ¼ 0.15, p ¼ 0.08) for the Solero device. This represented a 16%, 25% and 15% of the variability of EAV being explained by the applied ablation energy, when using the Acculis, Amica and Solero device, respectively.

Ablation volume irregularity
Correlation analysis between applied ablation energy and ablation volume irregularity (AVI) per ablation device is illustrated in Figure 6. Overall, the correlation between applied energy and AVI showed a weak correlation for the Acculis system (r ¼ 0.28, R 2 ¼ 0.08, p ¼ 0.02) and a moderate correlation for the Amica system (r ¼ 0.36, R 2 ¼ 0.13, p ¼ 0.07) and no correlation for the Solero system (r ¼ À0.01, R 2 ¼ 0.00013, p ¼ 0.96). This represented an 8%, 13% and 0.013% of the variability of EAV being explained by the applied ablation energy when using the Acculis, Amica and Solero device, respectively.

Factors influencing EAV/kJ
Univariable and multivariable linear regression analyses of factors potentially influencing obtained EAV per applied kJ energy are summarized in Table 2. The univariable analysis yielded tumor radius as the only statistically significant factor regarding EAV per applied kJ energy, with decreasing EAV per applied energy for each increasing millimeter (mm) tumor radius (B ¼ À0.03, p ¼ 0.002). In multivariable analysis, a larger distance between tumor center and healthy liver parenchyma (tumor radius) yielded significantly smaller ablation volume per applied kJ energy (B ¼ À0.03, p ¼ 0.01). KRAS mutated tumors yielded on average 14% larger ablation volumes per applied kJ ablation energy (p ¼ 0.09). Other factors did not show a significant influence on ablation volume expansion in our model.

Factors influencing MIV/kJ
Univariable and multivariable linear regression analyses of factors potentially influencing MIV per applied energy are summarized in Table 3. On univariable analysis, only the tumor radius was shown to be a significant factor influencing MIV/kJ (B ¼ À0.02, p < 0.001). On multivariable analysis, larger tumor radius remained leading to smaller MIV per applied energy (B ¼ À0.02, p < 0.001). All other factors did  Values are reported as the median and interquartile range (IQR). Abbreviations: EAV: effective ablation volume; PAV: predicted ablation volume; MEV: minimum enclosing ellipsoid volume; MIV: maximum inscribed ellipsoid volume; AVI: ablation volume irregularity.
not have a statistically significant influence on MIV per energy in our model.

Factors influencing MEV/kJ
Univariable and multivariable linear regression analyses of factors potentially influencing MEV per applied energy are summarized in Table 4. On univariable analysis, the tumor radius (B ¼ À0.07, p ¼ 0.002) and KRAS mutation status (B ¼ 0.44, p ¼ 0.02) were shown to be a significant factor influencing MEV/kJ. On multivariable analysis, larger tumor radius remained leading to smaller (B ¼ À0.07, p ¼ 0.008) and KRAS mutations leading to larger (B ¼ 0.44, p ¼ 0.04) MEV per applied energy. All other factors did not have a statistically significant influence on MEV per energy in our model.

Factor influencing AVI
Univariable and multivariable linear regression analyses of factors potentially influencing ablation volume irregularity (AVI) are summarized in Table 5. On univariable analysis, only the applied ablation energy was shown to be a significant factor influencing AVI (B ¼ 0.001, p ¼ 0.001). On multivariable analysis, higher ablation energy lead to increasingly irregularly shaped ablation dimensions (B ¼ 0.002, p ¼ 0.002). KRAS mutant status, steatotic liver parenchyma, a subcapsular tumor location or proximity to blood vessels did not have a statistically significant influence on AVI in our model.

Discussion
This study investigated ablation dimensions using a quantitative approach based on 3 D segmentation of effective ablation volumes, in a prospective multicentre patient cohort treated with SMWA for CRLM. To the best of our knowledge, this is the first work analyzing factors influencing ablation expansion and dimensions in 3 D using clinical data on ablation of malignant liver tumors. The weak to moderate linear correlations between obtained effective ablation volumes and applied ablation energies and predicted ablation volumes confirmed that other factors than ablation energy affect ablation dimensions [13,16]. At most, 25% of the variability in the EAV was explained by the applied energy, highlighting overall poor predictability of created ablation volumes solely by the applied ablation energy. The weak to moderate correlations between EAV and PAV are in concordance with previous works reporting significantly smaller clinical ablation volumes than the ones predicted in the manufacturer's IFU [14,30,31]. This can be explained by the physical properties of tumor and diseased (steatotic, fibrotic or cirrhotic) human liver tissue being different than in ex-vivo or in-vivo healthy animal tissue, where PAV are usually  extracted from [15]. A heat sink effect must further be expected leading to smaller ablation volumes in perfused livers. These findings supported our rationale for investigating factors potentially influencing ablation dimensions using multivariable analyses. We chose EAV per applied ablation energy to investigate effective ablation volumes, to i) include the main exposure parameter (applied energy) in the outcome variable and ii) create a potential surrogate parameter for the expansion of ablation energy across tumor and liver tissue. This is thought to be caused by microwave absorption, heat capacity, heat conductivity and other factors influencing heat generation within the ablated tissue. These parameters are difficult to assess in a clinical setting, especially when considering the heterogeneity between tumors, and would require tissue biopsies of the tumors to be ablated and the surrounding Figure 6. Correlation between energy and AVI for the Acculis, Amica and Solero systems. Lines represent mean correlations, 95% confidence intervals displayed as colored areas.   liver parenchyma. Whether a tendency toward larger EAV per energy in KRAS mutations of the primary tumor and hepatic steatosis is reproducible must be investigated in further studies. Also, KRAS mutation led to potentially larger ablation volumes (B ¼ 0.14, CI [-0.02, 0.3], p ¼ 0.09), which could be perceived counterintuitive as previous studies reported higher local recurrence rates in KRAS mutant tumors after thermal ablation [28]. The presence of satellite lesions, altered lipid metabolism and tissue composition related to genetic tumor properties [32,33], or different chemotherapy regimens applied depending on KRAS mutational status [34], might explain these differences. Lipids having a lower heat capacity than water molecules would result in less heating energy required for ablation expansion [35] and explain larger ablation volumes in steatotic livers or tumors with altered lipid metabolism. In contrast to other studies [26], we found no significant effect of hepatic steatosis surrogate measure on ablation volumes created per applied energy unit. This could be due to the low incident of patients in the hepatic steatosis surrogate measure subgroup (8 patients, 17 lesions). Even though these results were not statistically significant, they suggest that genetic mutations and potentially hepatic steatosis might influence the underlying tissue properties relevant for heat generation and conduction, such as lipid and water concentrations. This must be confirmed in clinical studies objectifying hepatic steatosis in the biopsy specimen.
Another interesting finding was that tumors with a larger radius created significantly smaller ablations per kJ ablation energy (B ¼ À0.03, p ¼ 0.01). Similar to the effect of hepatic steatosis shown in other studies, this could be explained by the higher water and lower lipid content in cancerous tissue as opposed to liver parenchyma [24]. Thus, with liver parenchyma of higher lipid concentration being closer to the positioned ablation probe in the tumor center, higher heat might be created faster at the same energy level. More tumor tissue of higher water content in the vicinity of the ablation probe (direct heat) causes gas bubble formation, which might act as a heat insulator due to the lower thermal conductivity of gases [35]. Enhanced vascularity in cancerous tissue compared to parenchyma might also contribute to this effect. Similar to Heerink et al. [14], the tumor volume represented a low percentage of the total ablation volume (7.7%, IQR ¼ [3.0%, 15.7%]), suggesting that physical characteristics of the surrounding parenchyma might contribute more to the ablation zone expansion than the characteristics of the tumor tissue. This holds especially for smaller tumors where an ablation volume increases to the power of 3 with the added margin. Other factors such as a subcapsular location and vessel proximity did not significantly influence but on average yielded smaller EAV per applied energy. This is probably explained by the inherent limitation of heat expansion at the liver capsule and due to the well-known heat sink effect [25].
To investigate surface irregularities of the generated ablation volumes and potential differences between the ablation volume yielding tumor necrosis and the one yielding spikes and potential collateral damage, a novel measure of approximating AVI was proposed. Hypothesizing that irregularities in the ablation volume surface caused by varying heat profiles around the MWA antenna might lead to increased risk of recurrence or need for re-ablation, we analyzed factors influencing such irregularities. The MIV described the volume where all tumor cells would be covered by the ablation volume with high certainty, the MEV the volume representing danger for inadvertent ablation of surrounding structures. While the perfectly ellipsoidal ablation volumes described by the manufacturer's guides would yield an AVI of 0, large differences between MIV and MEV and AVI heavily differing from 0 were shown in our clinical data, suggesting a significant variability in ablation shapes. Larger irregularity with higher applied energies was confirmed in multivariable analysis, with however a small effect size, corresponding to the maximum 13% variability in correlation analysis. The tendency toward larger ablation irregularities in subcapsular tumors might arise due to the typically wedge-shaped ablations generated in these locations. Even though not statistically significant, these findings and the correlation of AVI with ablation site recurrence must be validated in further clinical trials. Regarding MEV per applied energy, a significant (44%) increase in MEV per kJ of energy in KRAS mutated tumors was found. This could not be observed for MIV per kJ of energy, suggesting that KRAS mutations might have a higher effect on tissue properties toward the periphery of the tumor.
A limitation of this study is that even though prospectively acquired clinical data was used, ablation dimensions were analyzed retrospectively. This led to 29.5% of ablated tumors being excluded from analyses. Additionally, 70% of PAV had to be interpolated since predictions were not available from the manufacturers' IFU. This might have affected the precision of results from correlation analyses. Another notable aspect includes the fact that MWA is known to cause significant tissue shrinkage immediately post-treatment [31,36,37]. This might have led to an underestimation of EAV in this study since CT scans for ablation validation were taken immediately after the SMWA treatment.
In conclusion, this study confirms poor predictability of ablation dimensions solely based on the applied ablation energy. Other clinical factors influencing the physical properties of tumor and liver tissue affect the size and shape of ablation volumes after SMWA for CRLM. Further clinical trials including tissue sampling are needed to confirm these results.

Author contributions
IP, PT: Study planning, data collection, analysis and interpretation, manuscript writing; JF, SR, KPJ, SW: Study planning, analysis and interpretation, and manuscript proofreading. All authors provided approval for the publication of the final version of the manuscript.

Disclosure statement
Stefan Weber and Daniel Candinas are co-founders and shareholders of CAScination, manufacturer of one of the navigation systems applied for stereotactic microwave ablation in this study. Iwan Paolucci, Jacob Freedman, Simeon J. S. Ruiter, Koert P. de Jong and Pascale Tinguely have no conflict of interest to declare.

Funding
This project was partially funded by the Swiss Cancer League (Pascale Tinguely), the Professor Dr. Max Clo€ etta Foundation (Pascale Tinguely) and funds from the Region of Stockholm (ALF) (Jacob Freedman). All funding institutions had no contribution to the design or execution of this project.

Data availability statement
All data used in this study is available in the supplementary material "data.csv". The full statistical analysis and the code are available in the supplementary material "statistical_analysis.pdf".