Jin-Ping Li , Sheng Zho , Hui-Jie Jing , , Ho Jing , Lin-Hn Zhng , b , Zhong-Xing Shi ,Ting-Ting Fn , Song Wng

a Department of Radiology, the Second Affiliated Hospital of Harbin Medical University, Harbin 150086, China

b Department of Nuclear Medicine, the First Affiliated Hospital of Harbin Medical University, Harbin 150086, China

c Department of Radiology, Longhua Hospital, Shanghai University of Traditional Chinese Medicine, Shanghai 20 0 032, China

Keywords:Hepatocellular carcinoma Dual-energy Radiofrequency ablation Tumor response Texture analysis

A B S T R A C T

Introduction

In China, liver cancer ranks the second most common gastrointestinal cancer and is the third most common cause of cancerrelated deaths, annually accounting for approximately 369 0 0 0 deaths [1] . Primary hepatocellular carcinoma (HCC), a subtype of liver cancer, is one of the malignant neoplasms with a poor prognosis, and the primary treatments for HCC are surgical resection,systemic therapy, transplantation, and locoregional therapy [2] . Unfortunately, surgical resections are unavailable for most patients,because of tumor stage and liver reserve function; only 5%-15% of patients can undergo surgical resection [3] .

As a locoregional therapy, radiofrequency ablation (RFA) has become one alternative for patients who are not candidates for surgical resection [2] . Sometimes, for an HCC lesion with a diameter less than 2 cm and with unfavorable locations such as adjacent to large vessels, RFA is more cost-effective than surgical resection [4] .Correspondingly, RFA causes less damage and fewer complications than surgical resection, which is why the National Comprehensive Cancer Network (NCCN) recommends RFA as the first-line treatment option for small HCC [5] .

Although RFA is a minimally invasive therapy that can effectively reduce postoperative injury in patients, residual tumor lesions may present in the ablation zone due to its technical limitations. For instance, low ablation temperature, inadequate ablation time, or carelessly selected ablation zone can lead to incomplete tumor ablation. Compared with surgical resection, RFA had worse overall survival and a higher recurrence rate [6] . Therefore, it is still a challenge to effectively assess the post-RFA HCC-tumor response. Dual-energy computed tomography (DECT), also known as spectroscopy CT, may be one solution to this problem. DECT can be used to quantify the material composition of the lesion and thus,improve the possibility of distinguishing residual tumor in the ablation zone from surrounding non-tumor parenchyma [7] .

CT texture analysis (CTTA), a “radiomics” medical image interpretation technique, provides an aim and quantitative assessment of tumor heterogeneity by analyzing the distribution and relationship of pixel or voxel gray levels in images [8] , which can be used to predict pathological types [9] , assess treatment response [10] and estimate the prognosis in a variety of cancers [11] . In contrast, previous studies have demonstrated that clinical variables were also prognostic risk factors for late recurrence after RFA [ 12 , 13 ]. Notably, systemic immune-inflammation index (SII) [14] and albumin-bilirubin grade (ALBI) [15] , two variables based on laboratory tests, can be calculated by peripheral platelet count, neutrophil count, lymphocyte count, bilirubin, and albumin for simple calculations, which can be used to stratify poor prognosis of patients with HCC.

To summarize, numerous factors can be used to predict post-RFA recurrence. Hence, we utilized post-RFA-DECT imaging, combined with CTTA and clinical variables, to develop multiple predictive models to quantify post-RFA HCC-treatment outcome, and attempted to use these models to predict post-RFA tumor progression within 12 months.

Patients and methods

Patients

This was a retrospective study that was approved by the Medical Ethics Committee of the Second Affiliated Hospital of Harbin Medical University (KY2019-217). Due to the nature of this retrospective study, the patients’ written informed consent was waived.Of the 147 patients enrolled from August 2017 to September 2020,63 met the following inclusion criteria: (i) HCC was confirmed by histology or met the noninvasive diagnostic criteria of European Association for the Study of the Liver (EASL) [2]; (ii) patients underwent enhanced DECT with 2-4 weeks after RFA; (iii) the tumor size was less than 3 cm or the sum of two tumors less than 3 cm;(iv) patients received regular post-RFA enhanced CT or MRI followup within 12 months; (v) patients had no extrahepatic metastasis and/or vascular invasion. Exclusion criteria were as follows: (i) patients received pre-RFA anti-tumor therapy; (ii) patients were lost to follow up; (iii) patients had a history of other primary tumors.

RFA procedure

Preoperative ultrasound and other imaging data, such as preoperative CT or MRI, were used to determine the puncture point and needle path. After routine disinfection and spreading of the towel, the skin of the puncture site was under local anesthesia with 2% lidocaine. Procedures were performed with a commercial radiofrequency system (Celon Olympus, Berlin, Germany) and internal cooling bipolar electrode needles (Celon ProSurge, Berlin,Germany), and the radiofrequency output power was set at 40-80 W. The ablation time was 15-20 min, and the electrode was carefully guided to the intended location of the target lesion under ultrasound guidance. For lesions with clear tumor feeding vessels,we first placed electrodes through the normal hepatic parenchyma to block the tumor blood supply and then ablated HCC lesions.When performing ablation, the number and distribution of the radiofrequency-electrode needles depended on the size, shape, and location of tumors. Our strategy for complete tumor ablation was to ensure that the ablation zone was at least 5 mm from the tumor margin, with the aim of eliminating peritumor satellite foci and microvascular invasion to achieve complete pathological ablation [16] .

CT scanning protocol

A high-definition CT scanner (Discovery TM CT 750HD; GE Healthcare, Milwaukee, WI, USA) was used for three-phase dynamic DECT scanning. After the patients fasted for 4 h, the Discovery 750HD was used to scan the entire liver in the preset abdominal LIVER + GSI mode, and the arterial phase (AP), portal vein phase (PVP), and delayed phase were routinely acquired.The tube current and tube voltage settings were slightly different from patient to patient. The tube current range was approximately 20 0-40 0 mA and the tube voltage range 100-120 kV (80-140 kV transient high-speed switching). The layer thickness, layer spacing,pitch, and reconstruction interval were 5.0 mm, 5.0 mm, 0.984:1,and 1.25 mm, respectively. After establishing intravenous access for drug administration in the patient’s anterior elbow vein, we injected 60-80 mL of non-ionic iodine contrast (Iobitridol, Guerbet,Villepinte, France) at a rate of 3 mL/s using a double-barrel highpressure syringe (Stellant?, MEDRAD, Oelrich, Germany); we then injected 30 mL of saline at the same rate after the contrast agent injection was completed. The CT value of abdominal aorta at the level of the abdominal trunk was monitored to facilitate the scan(Smart Prep technique) with a monitoring threshold of 150 HU and a 6 s delay in starting the scan after reaching the threshold.

Acquisition of quantitative DECT data

DECT imaging was performed 2-4 weeks after RFA. All monochromatic images were transferred to a GE post-processing workstation (Advantage Windows 4.7) with a GSI viewer software package (GE Medical Systems, Waukesha, WI, USA). The radiologists selected the largest dimension of the ablated lesion diameter in the cross-section and manually outlined regions of interest (ROIs) along the tumor margin to include the ablation zone as much as possible. For each measurement, iodine concentration (IC),water concentration (WC), and effective atomic number (Zeff) were obtained based on the ROI; to ensure consistency in measuring the same ROI on AP and PVP images, the same ROIs were separately placed on AP and PVP images by the copy and paste function. In addition, we also outlined circular ROIs in the abdominal aorta as the background and acquired the slope of the energy spectrum(λAP(40-100keV)) in the energy range of 40-100 keV in AP; the calculation formula wasλAP(40-100keV)= (CT 1 - CT 2 )/60 keV, and the normalized iodine concentration (NIC) and normalized water concentration (NWC) were calculated as follows: NIC = IC lesion /IC aorta ,NWC = WClesion/WC aorta . To reduce measurement error, all data were measured by two independent radiologists who were blinded to patients’ information, while the final measurement result was taken as the average of the two measurements.

CT texture features extraction

Subsequently, the CT images of all patients were transferred into the free and open-source software 3D slicer for texture analysis. The ROI was manually traced around the tumor border on the slice of the ablation zone’s maximum diameter in AP.Each ROI was quantified by an extension ( https://www.slicer.org/wiki/Documentation/Nightly/Extensions/Radiomics ) to extract texture features, and a total of 56 features were extracted, including first-order texture features, a gray level co-occurrence matrix(GLCM), and a grey level dependence matrix (GLDM). The definitions of these features have been described in detail [ 17 , 18 ]; and no spatial resampling was performed.

The least absolute shrinkage and selection operator (LASSO) regression algorithm was performed using the R language “glmnet”package, and 10-fold cross-validation was used to select the optimal CTTA features. When the smallest optimal tuning parameter(Lambda) reached the minimum, the optimal CTTA features that may reflect the tumor response of patients within one year after the RFA were selected, and the data of each feature were refit with reference to the previous report [19] , which was calculated as follows:

CTTA-score =β0+β1x 1 +β2x 2 +β3x 3 + ... +βnx n

Where x n denotes the CTTA features selected by the LASSO regression algorithm;β0is the constant term of the regression equation;βn is the regression coefficient corresponding to the corresponding CTTA feature in the regression model. CTTA-score is the sum of the linear integrals of the constant terms, features, and their corresponding weight coefficients, and the CTTA-score can be calculated for each patient according to the formula.

Follow-up

DECT imaging was performed 2-4 weeks after RFA. The patients were usually monitored with dynamic hepatic CT or MRI or/and tumor markers in follow-up and every 2-3 months thereafter for 12 months, which was decided by the attending doctor in charge of the patient. Moreover, we reviewed the data on the baseline variables recorded by the hospital information system, including the patients’ age, sex, etiology of liver disease, bilirubin,albumin, ALBI ratio, platelets, neutrophils, lymphocyte, SII, alphafetoprotein (AFP), and carcinoembryonic antigen (CEA). ALBI was calculated with the formula ALBI = 0.66 × log10bilirubin (μmol/L)- 0.085 × albumin (g/L); SII was calculated with the formula SII = platelets × neutrophils/lymphocyte.

Two radiologists, who were unaware of the patients’ clinical information, performed the liver imaging reporting and data system (LI-RADS) treatment response (LR-TR) categories (“nonviable”,“equivocal”, and “viable”) on each follow-up image to assess the tumor response. According to LR-TR categories, we defined tumor progression as LR-TR viable reported; otherwise, we defined progression-free as LR-TR nonviable or LR-TR equivocal reported within 1 year of follow-up; the final treatment response was based on the worst LT-TR categories during follow-up. The follow-up ended until the patients’ ablation zone showed LR-TR viable on dynamic hepatic CT or MRI, or the follow-up time was exceeded 12 months, whichever came first. Patients with LR-TR viable underwent multidisciplinary discussions to reach a consensus on treatment, including but not limited to re-RFA, transcatheter arterial chemoembolization (TACE), surgical resection, and transplantation.

Statistical analysis

All statistical analyses were performed using R software (version 4.0.4, http://www.r-project.org/ ). Continuous variables were expressed as mean ± standard deviation (SD) or median (interquartile range, IQR) and compared using Student’st-test or Mann-WhitneyUtest. Categorical variables were expressed as frequency (percentage) and compared using Fisher’s exact test or Chi-square test. The univariate logistic regression model was used to determine clinical variables, DECT quantitative parameters, and a CTTA-score that could be included in the multivariate model.For variables withP< 0.05 in the univariate logistic regression model, they were included into the multivariate logistic regression model. Based on clinical variables, DECT quantitative parameters, and a CTTA-score, a multivariate logistic regression analysis was performed to construct a combined model via backward stepwise model. Taking LR-TR as the diagnostic criteria, the “pROC”package was used to draw the receiver operating characteristic(ROC) curves of the CTTA-score, clinical model, and clinical-texture model, and the DeLong test was used to compare the area under the curve (AUC) of each ROC curve. Furthermore, decision curve analysis (DCA) was used to evaluate the net benefits for clinical applications under different probability thresholds. Bootstrap resampling (20 0 0 times) was used for the internal validation of the clinical model and clinical-texture model. The effectiveness of prediction model diagnosis was evaluated by calculating the consistency index (C-index), and the calibration curve was drawn to explore the consistency between the frequency of the observed results and the predicted probability of the model. RMS package was used to draw a nomogram to visualize clinical-texture model. AP< 0.05 was considered statistically significant.

Results

Patient characteristics

Of the 63 patients enrolled in the study, 53 were males (median age 56 years, 53.0-60.0 years) and 10 were females (median age 57 years, 57.75-59.50 years). As shown in Fig. 1 , other patients were excluded for the following reasons: excessive time between the first DECT and RFA (n= 25), a tumor diameter greater than 3 cm or over 2 lesions (n= 42), patient was lost to follow up(n= 12), or the patient had a presence of extrahepatic metastases before RFA (n= 5). At the end of the 12-month follow-up, 63 patients were divided into two groups, 20 patients (31.7%) were classified as tumor progression, and the remaining 43 patients (68.3%)were classified as tumor progression-free. Among all baseline clinical variables, AFP (P= 0.007), ALBI (P< 0.001), and SII (P= 0.023)were significantly different between patients with tumor progression and tumor progression-free ( Table 1 ).

Feature extraction and selection

The LASSO algorithm and 10-fold cross-validation method were used to extract the optimal texture feature subset ( Fig. 2 A and B).Finally, 6 texture features were extracted based on 63 patient images to construct the CTTA-score ( Fig. 2 B), including (i) first-order feature, (ii) GLCM features, and (iii) GLDM features. The specific CTTA features and weight coefficients are shown in Table 2 , and the CTTA-score for each patient was obtained according to the following equation:

CTTA-score = -10.16831 -0.0 0 083?firstorder_Median -3.73924 ?glcm_Correlation + 0.059682?glcm_SumSquares + 0.546096?gldm_LDE + 0.001644?gldm_LDHGLE + 10.71753?gldm_LDLGLE

Fig. 2. The LASSO algorithm and 10-fold cross-validation were used to extract the best subset of texture features. A: When the optimal tuning parameter was minimum,log(Lambda) = -2.77634; 10-fold cross-validation was performed to select non-zero coefficient features; B: LASSO coefficient distribution of 56 texture features, including 6 non-zero coefficient features into the construction of CTTA-score; C: the bar graph of CTTA-score. The red and green bars indicate tumor progression and tumor progressionfree, respectively. D: The violin chart of CTTA-score between tumor progression and progressive-free patients. ??: P < 0.001.

Fig. 1. Flow chart for enrolling patients. HCC: hepatocellular carcinoma; RFA: radiofrequency ablation; DECT: dual-energy computed tomography.

Table 3 shows the difference between the median and IQR of the 6 optimal texture features and the calculated CTTA-score between patients with tumor progression and tumor progressionfree. Among 63 patients, the CTTA-score and the 6 optimal features were significantly different between patients with tumor progression and tumor progression-free (P< 0.05). Specifically, the CTTAscore for tumor progression was significantly higher than that for tumor progression-free (1.093 vs. -1.407,P< 0.001). The CTTAscores of the two groups were displayed as bar graphs and violin graphs ( Fig. 2 C and D).

Development of models and evaluation

To build a clinical model, the IC, WC, and Zeff were extracted from AP and PVP images;λAP(40-100keV), ALBI, and SII were entered into univariate logistic regression ( Table 4 ), and variables withP< 0.05 were entered to a multivariate logistic regressionmodel using backward selection. ALBI [odds ratio (OR) = 2.77,95% confidence interval (CI): 1.35-6.65,P= 0.010),λAP(40-100keV)(OR = 3.21, 95% CI: 3.16-5.65,P= 0.045), and IC AP (OR = 1.25, 95%CI: 1.01-1.62,P= 0.028) were risk factors of tumor progression (Cindex = 0.917).

Table 1 Comparison of baseline patient characteristics based on LR-TR.

Table 2 CTTA features selected by LASSO algorithm and their weight coefficients.

In addition, we developed a clinical-texture model based on a clinical model, and the CTTA-score was additionally entered in the multivariate logistic regression analysis; the results of the clinical-texture model demonstrated that ALBI (OR = 2.40, 95%CI: 1.19-5.68,P= 0.024),λAP(40-100keV)(OR = 1.43, 95% CI: 1.10-2.07,P= 0.019), and CTTA-score (OR = 2.98, 95% CI: 1.68-6.66,P= 0.001) were independent risk factors of tumor progression (Cindex = 0.957).

Furthermore, we plotted ROC curves to predict the accuracy of the clinical model, clinical-texture model, and CTTA-score, which were used to predict tumor progression in patients with HCC(within 12 months) after RFA; the AUCs of each ROC curve were compared to assess their diagnostic accuracy ( Table 5 ). As shown in Fig. 3 , the clinical model, clinical-texture model, and CTTA-score performed well in predicting tumor progression in patients with

Table 4 The univariate logistic regression of clinical features and CTTA-score between patients with tumor progression and tumor progressionfree.

HCC (within 12 months) after RFA (AUCs = 0.917, 0.962, and 0.906,respectively,P< 0.001). The DeLong test showed that the AUC of the clinical-texture model was superior to that of the clinical model (P= 0.040); although the AUC of the clinical model and clinical-texture model was greater than that of the CTTA-score,there was no significant difference (P= 0.800 and 0.116, respectively). The DCAs of the clinical model, clinical-texture model, and CTTA-score are shown in Fig. 3 B. DCAs showed that there were more net benefits than “treating all patients” or “not treating any patients” when using the clinical model, clinical-texture model, or CTTA-score to predict the probability of tumor progression within 12 months. Moreover, the clinical-texture model added more benefits than using only the CTTA-score or clinical model at any given threshold probability.

Nomogram construction and validation

Based on a clinical-texture model constructed by multivariate logistic regression analysis, we obtained 3 independent risk factors (ALBI,λAP(40-100keV), and CTTA-score) for tumor progression in patients with HCC (within 12 months) after RFA; we further developed a nomogram. As shown in Fig. 4 , a prediction score was assigned to each risk factor, and the sum of the three scores wason the total score axis, which corresponded to the prediction of the probability of tumor progression.

Table 3 Comparison of 6 selected texture features and CTTA-score between tumor progression and tumor progression-free.

Table 5 Prediction performance abilities of three models.

Fig. 3. Clinical model, clinical-texture model, or CTTA-score for clinical usefulness evaluation. A: Receiver operating characteristic (ROC) curve for clinical model, clinicaltexture model, or CTTA-score; B: decision curve analysis (DCA) for models and CTTA-score; the Y-axis represents net benefit, and the X-axis represents the threshold probability; compared with the clinical-texture model, or CTTA-score, the all-treatment strategy (gray line), no-treatment strategy (black dashed line), the models, and CTTA-score achieved more net benefits within most of the threshold probability range; moreover, at any given threshold probability, the clinical-texture model added more benefits than using only the CTTA-score or clinical models.

Fig. 4. Nomogram based on clinical-texture model incorporated ALBI, λAP (40-100 keV) , and CTTA-score. ALBI: albumin-bilirubin grade; CTTA: CT texture analysis.

The calibration curve indicated the degree to which the predicted probability and actual probability of the clinical model and clinical-texture model ( Fig. 5 ). To evaluate whether the prediction model was optimally calibrated, unreliability statistics were calculated to reflect the reliability of the calibration curve. The Hosmer-Lemeshow test showed that there was no significant difference between calibration curves and ideal curves in the clinical model and clinical-texture model (P= 0.792 and 0.716, respectively) and therefore, there was good agreement between the predicted values and the observed values.

Discussion

Fig. 5. Calibration curve of clinical model ( A ) and clinical-texture model ( B ). The horizontal axis represents the predicted probability, and the vertical axis represents the actual probability. The diagonal dashed line is the reference line, which represents the “ideal” situation where the predicted value is equal to the actual value. The solid blue line represents the performance of the model, the closer to the reference line, the better the predictions.

As a minimally invasive therapy, RFA can improve patients’postoperative quality of life with less postoperative injury and theoretically, better protect liver function. However, hepatectomy is superior to RFA in terms of overall survival and disease-free survival [20] ; the recurrence and progression-free survival for post-RFA patients with HCC ranged between 59.8% and 63.1% and 16.5 and 22 months, respectively [21–23] . Thus, routine imaging surveillance was crucial. LI-RADS is one of the most common imaging-assessment criteria used to effectively collect, interpret,report, and evaluate data in patients with HCC [24] . LR-TR categories can accurately predict viable tumors after locoregional therapy, and these categories correlate with overall survival [ 25 , 26 ].

We can quickly intervene to treat and benefit patients as soon as we discover the risk of HCC recurrence. An experimental animal study showed that DECT quantitative parameters can be effective on distinguishing tumor marginal areas from non-tumor hepatic parenchyma in VX2 rabbits at 14 days after RFA [27] . Therefore, we selected DECT images of patients about two weeks after RFA to attain the DECT quantitative parameters. Yue et al.found that the arterial iodine fraction calculated by dividing ICAPby IC PVP showed good performance in differentiating tumor remnants from adjacent hepatic parenchyma in HCC lesions after TACE treatment [28] , and the slope of the energy spectrum curve can be used to effectively distinguish the residual tumor foci from the tumor necrosis area, which was similar to our results. Contrastingly, our results showed thatλAP(40-100keV)was an independent risk factor for tumor progression in both clinical and clinicaltexture models (OR = 3.21 and 1.43, respectively). The potential reason was that the slope of the energy spectrum curve reflected the average attenuation characteristics of a given material; hightemperature radiofrequency can generate both coagulative necrosis and partial liquefied necrosis [29] , and if there was residual tumor tissue in the ablation zone, the attenuation characteristics may differ from the completely necrotic area, resulting in different slopes of the curve.

Similarly, texture analysis is also an effective tool for extracting image information. It is used to extract and process the pixel distribution within ROIs. Both CT/MRI and PET/CT imaging are suitable for analysis, which has been successfully used in distinguishing the pathological features of HCC [ 30 , 31 ]. Interestingly, based on CTTA and tumor size, Vosshenrich et al. constructed a multiparameter decision tree model and found that this model had an AUC of 0.96 in predicting the tumor response of HCC to TACE [32] .Although their treatment methods and model significantly differed from our study, our data supported the possibility of CTTA modeling as a predictor in tumor progression of patients with HCC. In contrast, we did not include tumor size as one of variables in this study. Because the edema band mostly disappeared two weeks after the RFA, it was still difficult to calculate the tumor size, and we manually outlined the ROI precisely along the tumor edge before performing texture feature extraction, which contained information on tumor size within the texture features. Considering the above, including tumor size, one of our variables may overfit the model. As a method of simply evaluating liver function by numerical calculation, the ALBI correlated with the prognosis of patients with HCC after various antitumor therapies [ 33 , 34 ]. In the present work, we found that ALBI was an independent risk factor for tumor progression in patients with HCC after RFA. Our data corroborated the findings of the previous work by Oh and colleagues; they demonstrated that ALBI was an independent risk factor for 5-year overall survival after RFA in patients with early-stage HCC [35] .

Our study has several limitations. Firstly, the sample size of our study was small, which may be insufficient to detect the differences of some indicators between tumor progression and progression-free. Secondly, we manually segmented lesions, but two radiologists consistently interpreted the DECT images to minimize this error. Finally, because this study was a single-center retrospective study, the long case collection time and lack of external validation may have affected the robustness of our conclusions.Therefore, further prospective studies are necessary to validate the performance of the proposed prediction model.

In conclusion, DECT quantitative parameters, CTTA, and clinical variables are valuable in predicting post-RFA tumor progression of HCC. The constructed clinical prediction model can be used to provide early risk warnings of potential tumor progression in patients after RFA, process postoperative efficacy evaluation, and formulate clinical personalized treatment plans.

Acknowledgments

None.

CRediT authorship contribution statement

Jin-Ping Li: Data curation, Formal analysis, Writing - original draft. Sheng Zhao: Methodology, Writing - original draft.Hui-Jie Jiang: Conceptualization, Funding acquisition, Supervision,Writing - review & editing. Hao Jiang: Investigation, Supervision.Lin-Han Zhang: Resources, Visualization. Zhong-Xing Shi: Project administration, Resources. Ting-Ting Fan: Software. Song Wang:Funding acquisition, Writing - review & editing.

Funding

This study was supported by grants from the National Key Research and Development Program of China ( 2019YFC0118100 ),the National Natural Science Foundation of China ( 81671760 ,81873910 and 62171167 ), the Natural Science Foundation of Shanghai ( 19ZR1457800 ).

Ethical approval

The study was approved by the Medical Ethics Committee of the Second Affiliated Hospital of Harbin Medical University(KY2019-217). All study procedures were carried out in accordance with theDeclarationofHelsinki.

Competing interest

No benefits in any form have been received or will be received from a commercial party related directly or indirectly to the subject of this article.