Noninvasive Imaging Evaluation Based on Computed Tomography of the Efficacy of Initial Transarterial Chemoembolization to Predict Outcome in Patients with Hepatocellular Carcinoma.

Purpose: This study aims to develop a new model to more comprehensively and accurately predict the survival of patients with HCC after initial TACE. Patients and
Methods: The whole cohort (n = 102) was randomly divided into a training cohort and a validation cohort in the ratio of 8:2. The optimal radiomics signatures were screened using the least absolute shrinkage and selection operator algorithm (LASSO) regression for constructing the radscore to predict overall survival (OS). The C-index (95% confidence interval, CI), calibration curve, and decision curve analysis (DCA) were used to evaluate the performance of the models. The independent risk factors (hazard ratio, HR) for predicting OS were stratified by Kaplan-Meier (K-M) analysis and the Log rank test.
Results: The median OS was 439 days (95% CI: 215.795-662.205) in whole cohort, and in the training cohort and validation cohort, the median OS was 552 days (95% CI: 171.172-932.828), 395 days (95% CI: 309.415-480.585), respectively (P = 0.889). After multivariate cox regression, the combined radscore-clinical model was consisted of radscore (HR: 2.065, 95% CI: 1.285-3.316; P = 0.0029) and post-response (HR: 1.880, 95% CI: 1.310-2.697; P = 0.0007), both of which were independent risk factors for the OS. In the validation cohort, the efficacy of both the radscore (C-index: 0.769, 95% CI: 0.496-1.000) and combined model (C-index: 0.770, 95% CI: 0.581-0.806) were higher than that of the clinical model (C-index: 0.655, 95% CI: 0.508-0.802). The calibration curve of the combined model for predicting OS presented good consistency between observations and predictions in both the training cohort and validation cohort.
Conclusion: Noninvasive imaging has a good prediction performance of survival after initial TACE in patients with HCC. The combined model consisting of post-response and radscore may be able to better predict outcome.

Introduction

Hepatocellular carcinoma (HCC) is one of the most common malignant tumors with high morbidity and mortality.1 Major risk factors for HCC include chronic alcohol consumption, diabetes, and hepatitis B or C virus infection.2 HCC was divided into five stages in Barcelona Clinic Liver Cancer (BCLC) staging system depending on disease extension, liver function, and performance status: very early stage (BCLC 0), early-stage (BCLC A), intermediate stage (BCLC B), advanced stage (BCLC C), and terminal stage (BCLC D).2 In terms of BCLC stage 0 and stage A, surgical therapies are preferred, supplemented by locoregional therapies. In terms of BCLC stage B, locoregional therapies are adopted preferentially, and the combination of other treatments or switching to systemic therapy will be applied when there are ineffective or resistant;3 as for BCLC stage C and stage D, systemic therapy and best supportive care are recommended, respectively. Individualized treatment is often chosen according to BCLC stages and the personal experience of doctors. The American Association for the Study of Liver Diseases (AASLD) guidelines state that transarterial chemoembolization (TACE) is the only recommended standard treatment for patients of HCC with BCLC stage B or intermediate stage.4 TACE can induce tumor ischemic necrosis via embolizing the supply artery of the tumor, and cause cytotoxicity in the tumor tissue at the same time,5 so as to prolong the overall survival of patients to a certain extent.

In recent years, the development and application of radiomics in the liver have become more widespread and in-depth, such as lesion segmentation, microvascular invasion, diagnosis tumor grade, classification, treatment planification, prognosis assessment, et al.6 In particular, radiomics are increasingly used for response assessment and survival prediction after locoregional therapies in patients with HCC, commonly based on computed tomography (CT) and magnetic resonance (MR), to build predictive models using preoperative noninvasive markers and to stratify the risk of disease, which can be of great advantage in guiding clinical treatment and benefiting patients.7–9 However, these previous studies were limited to the analysis of preoperative variables without clinical characteristics after TACE or postoperative radiological response, and the accuracy and completeness of their prediction remain a controversial issue. Han et al10 suggested that the response after the therapy may be helpful in survival prediction.

Through constructing survival prediction models, patients who have the best benefit after TACE can be selected out, and the risk factors can also be identified. If patients without the better benefit predicted by the models, they will be recommended to transfer to other locoregional or systemic therapies. In our study, we aim to develop a new model to more comprehensively and accurately predict the survival of patients with HCC after TACE.

Patients and Methods

Patients

The ethics committee approval of the First Affiliated Hospital of Sun Yat-sen University was obtained for this study, and informed consent was waived due to the retrospective nature of this study. This unicentric retrospective study enrolled 102 patients with HCC who were initially treated with TACE from the First Affiliated Hospital of Sun Yat-sen University between January 2013 to December 2019. All patients were randomly divided into training and validation cohort in the ratio of 8:2. The inclusion criteria were: a) HCC was diagnosed by histological or noninvasive radiologic criteria;11 b) patients were received initial TACE; c) patients were accepted CT enhancement both in pre-TACE and post-TACE; d) patients underwent laboratory tests after initial TACE; e) after the initial TACE, repeat TACE or other local treatment was chosen depending on the patient’s condition and the clinician’s experience. The exclusion criteria were: a) patient has a history of other cancers in addition to HCC; b) patients received other invasive treatments in initial pre-TACE or post-TACE, such as radiofrequency ablation (RFA), surgery, radiotherapy, et al; c) poor quality of CT images that might affect analysis; d) patients had incomplete clinical and laboratory data. The detailed selection flow chart is shown in Figure 1.

Figure 1

Flowchart of patients’ enrollment.

Treatment Procedure and Follow-Up

The treatment included the conventional TACE and TACE with blank or drug-loaded microspheres. All procedures via a traditional femoral artery approach, with the introduction of a 5F catheter and then a 5F YASHIRO catheter. Under the guidance of the guidewire, selective catheterization was performed to the abdominal trunk, and high-pressure syringes were used for angiography. According to the angiographic findings, a 2.7F microcatheter was used to enter the tumor supply branch with super-selective intubation guided by the guidewire followed by injection of embolization drugs and chemotherapy drugs. The embolization drugs include mixed suspension (2 mL - 20 mL) both with lipiodol (2 mL - 20 mL) and epirubicin (20 mg - 40 mg) or microspheres (100 um - 300 um) followed by injection of 350–560um gelatin sponge, which volume of the emulsion was determined based on the tumor volume. The chemotherapy drugs include epirubicin (10 mg - 80 mg), raltitrexed (4 mg - 5 mg), irinotecan (240 mg), nedaplatin (40 mg), and oxaliplatin (200 mg - 400 mg). When the lipiodol deposition can be seen in the re-imaging and the blood supply to the tumor is obviously reduced, the operation is completed. The dose of embolization drugs and chemotherapeutic drugs was based on tumor burden and patients’ characteristics.

Contrast-enhanced CT was performed within 4–6 weeks after the operation to check the embolization of the tumor and determine whether there was complete necrosis or “viable tumor”.12 According to the modified Response Evaluation Criteria in Solid Tumors (mRECIST),13 the postoperative radiological responses of tumors were divided into four categories: complete response (CR), partial response (PR), stable disease (SD), and progressive disease (PD). Radiological responses were independently assessed by two radiologists with 5–10 years of experience and were blinded to clinical information. After the disease stabilized, patients were followed up with contrast-enhanced CT every 8–12 weeks. All included patients were followed up by telephone or medical record system for overall survival (OS). OS was defined as the time interval between initial TACE and death from a cause related to this study. Data from patients who were lost to follow-up or survived at the last follow-up (June 16, 2021) were considered censored data.

Image Acquisition and Segmentation

All patients underwent contrast-enhanced CT both pre-TACE and post-TACE using a 64-slice spiral CT scanner (Aquilion 64; Canon Medical Systems). The tube current is 250mA, the tube voltage is 120Kv; the layer thickness is 5mm, and all the thin layer reconstruction (1mm) is performed. The contrast agent (Ultravist, Bayer, Germany) was injected via a high-pressure syringe through the antecubital vein at a rate of 3.0 mL/s, while 40 mL of saline was injected at the same rate, and CT images of the arterial phase and venous phase were obtained at 35s and 65s, respectively.

Firstly, images of all patients were transmitted into the radiomics platform (Huiying Medical Technology Co., Ltd). Secondly, the outline of the region of interest (ROI) was chosen to be sketched on 1 mm thickness images on the arterial phase independently by two reviewers with more than 5 years of abdominal imaging experience, both of whom were blinded to the clinical and laboratory information. Thirdly, the ROI is outlined as close to the tumor margin as possible, excluding normal liver tissue, bile ducts, and larger blood vessels. The volume of interest (VOI) was attained after continuous outlining of ROI on all layers of the lesion, and the first and last layer sections were excluded due to the presence of partial volume effects. The detailed outline details are shown in Figure 2.

Figure 2

Overview of radiomics analysis in this study. The volume of interest (VOI) was formed after region of interest (ROI) outlining on the arterial phase of CT enhancement. Features were extracted from the VOI and selected using least absolute shrinkage and selection operator (LASSO) regression. A combined radscore-clinical model was developed to predict overall survival (OS) in patients with hepatocellular carcinoma (HCC) after transarterial chemoembolization (TACE), and a Kaplan–Meier (K-M) stratification analysis of risk factors for predicting OS was performed.

Radiomics Feature Extraction and Selection

To decrease the variability of radiomics features, the image data is normalized using a z-score. A total of 1409 quantitative radiomics features were extracted from the VOI. The first-order statistics comprised 18 features that reflected the quantitative depiction of the distribution of voxel intensity in CT images. The shape-based features comprised 14 features that reflected the shape and size of the region. The texture features comprised 24 gray level co-occurrence matrix features, 16 gray level run length matrix features, 16 gray level size zone matrix features, 5 neighboring gray tone difference matrix features, and 14 gray level dependent matrix features, which can quantify regional heterogeneity differences. Additionally, we applied higher-order filter features—logarithm, exponential, gradient, square, square root, local binary pattern (LBP), and wavelet filters—to images and then extracted first-order statistics and texture features based on the filtered image (1302 filter features). The interpretations of these radiomics features are displayed in the Pyradiomics file (Ver. 2.2.0) at .

The inter-class correlation coefficient (ICC) was calculated for evaluation of the consistency of the ROI sketching by the two readers. Reader 1 and Reader 2 randomly selected the ROIs of 30 CT images from the whole cohort, and Reader 1 repeated these operations two weeks later. If a patient’s ICC value is greater than 0.8, the manual tracing is considered to be in good compliance, and subsequently, the images depicted by Reader1 will be used in the model construction. The intra-class correlation coefficient (ICC) was calculated to ensure the stability and reproducibility of each radiomic signature. When ICC is greater than 0.8 it indicated good consistency and reproducibility of feature extraction and is retained for further radiomics evaluation, while other features are excluded.

Construction of the Radiomics Signature

To prevent overfitting of the radiomics signature, we obtained the optimal radiomics features using the least absolute shrinkage and selection operator algorithm (LASSO) in the training cohort in Anaconda 3 platform () with scikit-learn () and matplotlib packages (). Additionally, the Spearman correlation was used to evaluate the mutual correlation among the initially selected features. When the value is less than 0.9, it is considered that there is no correlation between the initially selected features.14 Then, the selected features with nonzero coefficients were applied to construct a radiomics signature. Each patient’s radiomics score (radscore) was calculated by a linear combination of selected features, weighted by their respective LASSO coefficients.

Statistical Analysis

Statistical tests were performed using SPSS (version 25.0, IBM) and R statistical software (version 3). Categorical variables were summarized as number (percent), and continuous variables were summarized as mean (standard deviation, Sd) or median (interquartile range, IQR) based on the underlying distribution of the data. The Chi-square test or Fisher’s exact test was used for categorical variables, and the t-test was used for continuous variables. Independent risk factors for predicting OS using univariate and multivariate cox regression analysis in the training cohort, including basic clinical characteristics, CT images features, and radscore. Then, nomograms for individualized prediction of survival probability were established based on the radscore, clinical model, and radscore-clinical (combined) model, respectively.

Harrell’s consistency index (C-index) was used to quantify the discriminatory performance of the models with values ranging from 0.5–1. When the value is 0.5, it means that the model has no predictive effect; when the value is 1, it means that the predicted results of the model are in perfect agreement with the observed results. The calibration curve was conducted to evaluate the consistency (goodness of fit) of the model between the actual and predicted mortality risks. The decision curve analysis (DCA) was performed to quantify the clinical utility of radiomic features from the training cohort and validation cohort to assess the ability of derived radiomic features. Survival curves were estimated using the Kaplan-Meier (K-M) method, and Log rank tests were used to compare. K-M analysis and Log rank test were used to stratify the different risk groups to assess disease progression, and hazard ratio (HR) and 95% confidence interval (CI) were used to measure patient survival differences. A two-sided P < 0.05 was considered that the differences were statistically significant.

Results

Patient Characteristics

A total of 102 patients with HCC were enrolled in our retrospective study (Figure 1). The mean age for the whole cohort was 54.88 ± 11.20 years, and males accounted for the majority (n = 91, 89.2%). Most patients were evaluated as BCLC-stage B (n = 58, 56.9%), and Child-Pugh class A (n = 76, 74.5%). The whole cohort was proportionally divided into a training cohort (n = 81) and a validation cohort (n = 21), and the mean age was 54.78 ± 11.43 years, 55.29 ± 10.85 years, respectively, and males accounted for the majority similarly (n = 71, 87.7%; n = 20, 95.2%, respectively). After statistical analysis, we found that the two cohorts were comparable in terms of basic clinical characteristics and CT images features (P > 0.05), except for etiology (P = 0.021). Details were shown in Table 1.

Table 1

Characteristics of Patients in the Training and Validation Cohorts

Characteristic	Training Cohort, N (%)	Validation Cohort, N (%)	P
Patients	81 (79.4)	21 (20.6)
Gender
Male	71 (87.7)	20 (95.2)	0.318
Female	10 (12.3)	1 (4.8)
Age, year (mean ± Sd)	54.78 ± 11.43	55.29 ± 10.85	0.855
Etiology
HBV	10 (12.3)	7 (33.3)	0.021
Non-HBV	71 (87.7)	14 (66.7)
BCLC‐Stage
B	47 (58.0)	11 (52.4)	0.642
C	34 (42.0)	10 (47.6)
Child‐Pugh class
A	63 (77.8)	13 (61.9)	0.137
B/C	18 (22.2)	8 (38.1)
Longest diameter, mm (median, IQR)	81.96 (43.35, 111.22)	73.64 (41.23, 103.90)	0.585
Tumor number
1	36 (44.4)	10 (47.6)	0.794
≥ 2	45 (55.6)	11 (52.4)
Pre‐vascularization15
Type-1	9 (11.1)	2 (9.5)	0.973
Type-2	6 (7.4)	2 (9.5)
Type-3	26 (32.1)	6 (28.6)
Type-4	40 (49.4)	11 (52.4)
Pre‐pseudocapsule
Complete	46 (56.8)	12 (57.1)	0.977
Incomplete	35 (43.2)	9 (42.9)
Pre‐cirrhosis (CT)37
I/II	67 (82.7)	17 (81.0)	0.850
III/IV	14 (17.3)	4 (19.0)
Pre‐thrombus
Absent	48 (49.3)	11 (52.4)	0.569
Present	33 (40.7)	10 (47.6)
Pre‐AFP, ug/L
≤ 400	40 (49.4)	11 (52.4)	0.807
> 400	41 (50.6)	10 (47.6)
Pre‐ALB, g/L
≤ 35	36 (44.4)	11 (52.4)	0.516
> 35	45 (55.6)	10 (47.6)
Pre-ALT, U/L
≤ 40	43 (53.1)	7 (33.3)	0.107
> 40	38 (46.9)	14 (66.7)
Pre‐AST, U/L
≤ 37	18 (22.2)	5 (23.8)	0.877
> 37	63 (72.8)	16 (76.2)
Pos‐ALB, g/L
≤ 35	36 (44.4)	9 (42.9)	0.896
> 35	45 (55.6)	12 (57.1)
Post‐ALT, U/L
≤ 40	24 (29.6)	6 (28.6)	0.924
> 40	57 (70.4)	15 (71.4)
Post‐AST, U/L
≤ 37	8 (9.9)	1 (4.8)	0.461
> 37	73 (90.1)	20 (95.2)
Post‐response
CR	9 (11.1)	1 (4.8)	0.594
PR	36 (44.4)	12 (57.1)
SD	23 (28.4)	4 (19.0)
PD	13 (16.0)	4 (19.0)

Abbreviations: AFP, alpha-fetoprotein; ALB, albumin; ALT, alanine aminotransferase; AST, aspartate aminotransferase; CR, complete response; HBV, hepatitis B virus; IQR, interquartile range; PD, progressive disease; PR, partial response; SD, stable disease; Sd, standard deviation.

Feature Selection, and Radiomics Signature Building

Noninvasive imaging was used to identify features for characterizing the survival of patients with HCC treated with initial TACE. A total of 1409 features were extracted from the VOI of each patient, the 1248 features (88%) had satisfactory intra-observer and inter-observer reproducibility (ICC > 0.8). Subsequently, using the regularized regression with the penalty (α is denoted as the weight of penalty term), the number of features was reduced to 9 (Figure 3). The result of Spearman correlation showed that coefficients between these selected features were less than 0.9 (). Finally, the radscore was constructed by the 9 wavelet features, including 2 first-order features, 1 GLCM feature, 1 GLRLM feature, and 5 GLSZM features (Table 2). The radscore was calculated using the following formula:

Figure 3

The results of least absolute shrinkage and selection operator (LASSO) regression. (A) Mean squared error (MSE) path. (B) Lasso path, 9 optimal radiomics features were selected to calculate the radscore.

Table 2

Characteristics of Each Radiomics Feature Extracted and Corresponding Coefficients for Predicting Overall Survival (N = 102)

	Filter Type	Feature Class	Statistic	Coefficients
	Intercept	–	–	1.130
R1	wavelet.HHH	First-order	Entropy	0.116
R2	wavelet.LHH	First-order	Mean	0.048
R3	wavelet.HLH	GLCM	DE	0.577
R4	wavelet.HHL	GLRLM	HGLRE	0.028
R5	wavelet.HLH	GLSZM	GLV	0.293
R6	wavelet.HLL	GLSZM	SZNN	0.005
R7	wavelet.HLL	GLSZM	HGLZE	0.273
R8	wavelet.LLH	GLSZM	LGLZE	0.022
R9	wavelet.LLH	GLSZM	SALGLE	0.245

Abbreviations: DE, difference_entropy; GLCM, gray level co-occurrence matrix; GLRLM, gray level run length matrix; GLSZM, gray level size zone matrix; GLV, gray_level_variance; HGLRE, high_gray_level_run_emphasis; HGLZE, high_gray_level_zone_emphasis; LGLZE, low_gray_level_zone_emphasis; SALGLE, small_area_low_gray_level_emphasis; SZNN, size_zone_non-uniformity_normalized.

radscore = 1.130430763 + (wavelet.LHH_First-order_Mean × 0.047667017) + (wavelet.HLL_GLSZM_Size Zone Non-Uniformity Normalized × 0.004595354) + (wavelet.HLL_GLSZM_High Gray Level Zone Emphasis × 0.272802634) + (wavelet.LLH_GLSZM_Low Gray Level Zone Emphasis × 0.022080334) + (wavelet.LLH_GLSZM_Small Area Low Gray Level Emphasis × 0.244658219) + (wavelet.HLH_GLCM_Difference Entropy × 0.577169645) + (wavelet.HLH_GLSZM_Gray Level Variance × 0.292887147) + (wavelet.HHH_First-order_Entropy × 0.116398082) + (wavelet.HHL_GLRLM_High Gray Level Run Emphasis × 0.027502833).

Development and Validating of the Models in Survival Prediction

By univariate cox regression analysis, of all the baseline characteristics and radscore, 8 features were risk factors for predicting OS, including BCLC-stage, Child-Pugh class, longest diameter, pre-vascularization,15 pre-thrombus, pre-AFP, post-response, and radscore. Then, after multivariate cox regression analysis, we found that post-response (HR: 1.880, 95% CI: 1.310–2.697; P = 0.0007) and radscore (HR: 2.065, 95% CI: 1.285–3.316; P = 0.0029) were independent risk factors for predicting OS (Table 3). Subsequently, in the training cohort, according to the two independent risk factors, we built the radscore, clinical model, and combined model, respectively.

Table 3

Univariate and Multivariate Cox Regression Analysis for Predicting Overall Survival in Hepatocellular Carcinoma Treated by Transarterial Chemoembolization in the Training Cohort

Variable	Univariate Analysis		Multivariate Analysis
	HR (95% CI)	P	HR (95% CI)	P
Gender	1.480(0.671–3.264)	0.358		–
Age	0.996(0.973–1.020)	0.751		–
Etiology	1.544(0.701–3.403)	0.258		–
BCLC‐stage	2.225(1.315–3.766)	0.003	2.712(0.294–25.040)	0.382
Child‐Pugh class	2.293(1.293–4.065)	0.008	1.427(0.741–2.748)	0.290
Longest diameter	1.010(1.004–1.017)	0.002	1.004(0.995–1.013)	0.378
Tumor number	1.007(0.602–1.686)	0.978		–
Pre‐vascularization	2.983(1.802–4.977)	< 0.0001	0.940(0.682–1.295)	0.705
Pre‐pseudocapsule	1.656(0.979–2.801)	0.065		–
Pre‐cirrhosis (CT)	1.165(0.570–2.377)	0.682		–
Pre‐thrombus	2.209(1.302–3.746)	0.005	0.489(0.057–4.226)	0.518
Pre‐AFP	2.431(1.447–4.085)	0.001	1.442(0.769–2.704)	0.257
Pre‐ALB	0.6531(0.391–1.092)	0.106		–
Pre‐ALT	0.878(0.526–1.466)	0.620		–
Pre‐AST	1.222(0.672–2.224)	0.507		–
Post‐ALB	0.806(0.484–1.342)	0.410		–
Post‐ALT	0.888(0.515–1.532)	0.673		–
Post‐AST	1.929(0.700–5.321)	0.507		–
Post-response	2.091(1.501–2.915)	< 0.0001	1.880(1.310–2.697)	0.0007*
Radscore	2.305(1.551–3.426)	< 0.0001	2.065(1.285–3.316)	0.0029*

Note: *P<0.05, the difference is statistically significant.

Abbreviations: HR, hazard ratio; CI, confidence interval.

In the radscore and clinical model, the C-index was 0.834 (95% CI: 0.711–0.958), 0.694 (95% CI: 0.675–0.755), respectively, with the latter being significantly less effective than the radscore (). The efficacy of the combined model still showed a better result with a C-index of 0.806 (95% CI: 0.697–0.953) (Figure 4A), which was slightly lower than that of radscore. In the validation cohort, the results displayed that the efficacy of the radscore was also significantly higher than that of the clinical model, with a C-index of 0.769 (95% CI: 0.496–1.000), 0.655 (95% CI: 0.508–0.802), respectively. Notably, the efficacy of the combined model was higher than that of the radscore, with a C-index of 0.770 (95% CI: 0.581–0.806) (Table 4). The calibration curve of the combined model for predicting OS presented good consistency between observations and predictions in both the training and validation cohorts (Figure 4B and C). The DCA in the training cohort showed that the radscore and combined model had higher overall net benefit than the clinical model (Figure 5).

Figure 4

Development and validation of the model. (A) The nomogram of combined model for predicting overall survival (OS) was consisting of post-response (1= complete response, 2 = partial response, 3 = stable disease, 4 = progressive disease) and radscore, with a C-index of 0.806 (95% CI: 0.697–0.953) in the training cohort. The calibration curve demonstrating predictions from the model to the actual observed probability in the training cohort (B) and validation cohort (C).

Table 4

Predictive Performance of the Survival Models

Prediction Model	Training Cohort		Validation Cohort
	C-Index	95% CI	C-Index	95% CI
Radscore	0.834	0.711–0.958	0.769	0.496–1.000
Clinical model	0.694	0.675–0.755	0.655	0.508–0.802
Combined model	0.806	0.697–0.953	0.770	0.581–0.806

Figure 5

Decision curve analysis of the models. The net benefit of both the radscore and combined model are higher than that of clinical model.

Survival Prediction

A total of 74 patients in our whole cohort had been completed follow-up for OS. In addition, left-censored data occurred in 12 cases due to transfer to another hospital, death from other causes, and loss to follow-up, and right-censored data occurred in 16 cases due to survival at the end of follow-up. The median OS was 439 days (95% CI: 215.795–662.205) in the whole cohort, and in the training cohort and validation cohort, the median OS was 552 days (95% CI: 171.172–932.828), 395 days (95% CI: 309.415–480.585), respectively (P = 0.889) (Figure 6A).

Figure 6

Kaplan–Meier analysis of overall survival (OS). (A) Kaplan–Meier analysis of OS between training cohort and validation cohort (P = 0.889). (B) Kaplan–Meier analysis of OS of radscore (low-score and high-score) divided by the cut-off value (1.36) in the training cohort (P = 0.0001). (C) Kaplan–Meier analysis of OS of post-response evaluated by mRECIST (CR, PR, SD, and PD) in the training cohort (P < 0.0001).

Furthermore, we stratified two independent risk factors for predicting OS selected in the training cohort. To facilitate the analysis of clinical practice, we determined the cut-off value (1.36) of radscore by using maximally selected log-rank statistics and classified it as the low-score subgroup and high-score subgroup (). In the training cohort, the median OS was 885 days (95% CI: 568.839–1201.161) in the low-score subgroup, and the median OS was 220 days (95% CI: 155.588–284.412) in the high-score subgroup, which was significantly lower than that in the low-score subgroup (HR: 3.296, 95% CI: 1.825–5.954; P = 0.0001). In addition, the mortality rate of the high-score subgroup is explicitly higher than that of the low-score subgroup, with values of 93.55% and 60.00%, respectively. In terms of post-response, we divided the four subgroups: CR, PR, SD, and PD. The result demonstrated that the median OS in the four subgroups was 944 days (95% CI: 602.734–1285.266), 741 days (95% CI: 468.055–1013.945), 190 days (95% CI: 175.061–204.939), and 185 days (95% CI: 57.845–312.155), respectively, with statistically significant differences among the four subgroups (P < 0.0001). Firstly, taking CR as the reference, the HR of OS in PR, SD, and PD were 1.308 (95% CI: 0.631–2.712), 2.881 (95% CI: 1.235–6.724), and 5.647 (95% CI: 1.597–19.976), respectively. Secondly, taking PR as the reference, the HR of OS both in SD and PD were 2.202 (95% CI: 1.157–4.1912) and 4.316 (95% CI: 1.385–13.453), correspondingly. Thirdly, in the PD subgroup, there is still a higher risk of death than that in the SD subgroup (HR:1.960, 95% CI: 0.581–6.614). Furthermore, the mortality rates of CR and PR were significantly lower than those of SD and PD, with values of 55.56%, 66.67%, 86.96%, and 76.92%, respectively. Survival curves according to radscore and post-response were shown in (Figure 6B and C).

Discussion

Although noninvasive imaging has been proved to be important in survival prediction, there is still no better model to accurately predict survival in patients with HCC undergoing TACE for the first time. In our study, we developed the new, more comprehensive combined model to predict OS (C-index = 0.806), including radscore consisted of nine radiomics signatures and a clinical model composed of one characteristics of post-response. Furthermore, we discovered that in addition to radscore in the pre-TACE, radiological response in the post-TACE was also proven to be an independent risk factor for predicting OS. Moreover, to avoid overfitting and optimize the models, we combined the radscore with the clinical model. The result showed that the efficacy of the combined model was slightly higher than that of the radscore in the validation cohort, both of whom were significantly better than the clinical model.

Recently, most studies have attempted to use artificial intelligence (AI) to assess the efficacy and predict the survival of patients with HCC treated with different modalities in order to help patients choose the optimal treatment modality for the best clinical benefit. Machine learning, as a branch of AI, includes support vector machines (SVMs), artificial neural networks (ANNs), regression trees, and deep learning, which is a new machine learning model using technologies such as convolutional neural networks (CNNs) to develop multilayer neural network algorithms.6,16 Although both the European Association for the Study of the Liver (EASL) and AASLD guidelines recommend that the treatments of HCC are selected according to the BCLC staging system, the actual clinical treatment has not been adopted strictly in accordance with the recommendations, which may be related to the experience of doctors and the systemic condition of patients. Choi et al17 developed a machine learning-based clinical decision support system for recommending initial treatment options and predicting OS for HCC, aiming to overcome the complexity of treatment decisions and obtain the survival prediction curves for each patient based on the chosen survival models. This finding could provide personalized reference options for clinicians to develop treatment plans, which may be a direction for further research in the future. Zhang et al18 built models based on machine learning to predict OS of HCC after surgical resection and indicated that radiomic signature could predict OS well (C-index = 0.88), and the prediction performance can be further improved by combining clinical-radiological features (C-index = 0.92). Similarly, Wu et al19 and Kim et al20 explored the importance of radiomic signature and clinical characteristics in predicting survival of patients with HCC, who were treated with stereotactic body radiation therapy or TACE, and both of them showed that the combined model outperformed the clinical characteristics model alone in predicting survival. It is worth acknowledging that machine learning has considerable value in the field of diagnosis and treatment in patients with HCC, especially in the efficacy assessment and survival prediction of TACE because of the relatively large amount of TACE data available in the clinic.

In previous studies,9,20,21 for predicting survival of HCC treated with TACE, the categories and number of features in the models were diversified, mainly including tumor size, tumor number, BCLC-stage, Child-Pugh class, AFP, irregular tumor margin, and rad score. Similarly, our model differs from the models developed in previous studies, but there has some overlapping features. Firstly, in terms of radscore, which is essential feature for predicting survival, the number and type vary in each radscore model as reported in previous studies. In the study of Niu et al,21 radscore consisted of eight features (including two volumetric texture features and six wavelet texture features) that were screened using LASSO regression analysis. In our study, 9 wavelet transform features were filtered to form the radscore by the same method. It has been suggested that wavelet-based features can be used for disease diagnosis and treatment response prediction22,23. GLCM, GLRLM, and GLSZM are all matrix-based features: GLCM describes the pairwise arrangement of voxels with the same gray value and is used to highlight local heterogeneity information; GLRLM_HGLRE is used to measure the distribution of high gray values, and the HGRE value is expected to be larger for images with higher gray values.24 A study25 comparing the differences in radiomics features in normal liver, cirrhosis, and hepatocellular carcinoma parenchyma pointed out that first-order, GLCM, GLRLM, and GLSZM had low variation, and the values of two radiomics signatures (first-order entropy and GLCM-sum entropy) were consistently larger, and the value of GLCM inverse variance consistently decreased from healthy liver to cirrhosis and HCC. Entropy reflects the degree of randomness of a region, and the larger its value, the higher the image complexity. Furthermore, we performed a high-low risk subgroup analysis of the radscore based on the cut-off value and found a statistically significant difference in median OS (885 days vs 220 days, P = 0.0001). Our radscore is similar to that of Niu et al, including three GLSZM features (gray level variance, size zone non-uniformity normalized, and low gray level zone emphasis). The C-index of radscore was 0.834 and 0.769 in the training and validation cohort, respectively, which the prediction performance may be better than that of Niu et al (0.723, 0.734, correspondingly). Meng et al9 filtered six features (one Log.sigma.1.0.mm and five wavelet texture features) to construct the radscore, with a C-index of 0.68 in the training cohort, which is lower than our model. However, there were no overlapping features with our radscore. Moreover, Kim et al20 screened five radiomics signatures to build the radscore, which had no features overlapping with ours as well. Since Kim et al’s study did not use the C-index as a quantitative indicator of the model, it was not possible to compare the predictive performance.

From the above discussion, we summarize several reasons that can be reasonably explained. First, both in our study and in Niu et al’s study, the timing of image acquisition was obtained in the arterial phase and venous phase at 35s, 65s, respectively. Second, in terms of scanning devices and parameters, both ours and Niu et al's study involved only 1 device (64-row CT) and the same tube voltage of 120Kv, whereas the other two studies involved 2 and more scanning devices and the tube voltage was not constant. Nonetheless, a recent systematic review26 noted that no two studies can be directly compared with each other because proprietary feature extraction tools generate thousands of quantitative variables that have no meaning outside the context of individual studies, and suggested that standardization of protocols and outcome measures, sharing of algorithms and analysis methods, and external validation are necessary to avoid model heterogeneity.

However, radscore is characteristic based on pre-TACE, which may not be comprehensive enough for survival prediction. Subsequently, we analyzed the characteristics in the post-TACE as well and found that post-response was also an independent predictor of OS (P = 0.0007). Previously, the importance of post-response for predicting survival has rarely been discussed in previous studies. Before that, Han et al10 developed a pre-TACE model consisting of baseline clinical characteristics, and a post-TACE model including response after TACE. Then, they compared their models with the current hepatocellular artery embolization prognosis (HAP) score and its modified version (including tumor size, tumor number, AFP, bilirubin, and albumin), and emphasized that the post-TACE model had the best predictive performance (C-index = 0. 723). Likewise, the post-response after initial TACE was an independent risk factor for predicting the survival in our study. A study27 pointed that GLCM moment based on pre-TACE dynamic CT texture is one of the significant predictors of CR after TACE. Currently, mRECIST is an accurate, and convenient criteria for assessing radiological response. Previous studies28,29 pointed out that the mRECIST and EASL criteria were more accurate than the RECIST 1.1 criteria in predicting the rate of response and overall survival, and the mRECIST and EASL criteria are highly consistent and predictable. Furthermore, the liver function indicators of pre-TACE and post-TACE included in our study were not found to be predictive of survival. Similarly, in most previous studies laboratory indicators of liver function besides AFP were not found to be significant. In terms of the preoperative pattern of tumor vascularization, Hasdemir et al30 indicated that the pattern of vascularization of the largest tumor could predict survival after TACE (P = 0.012). Nevertheless, in our study, this result was not be found. The possible reason is that the study by Hasdemir et al included all intermediate stage HCC and the efficacy assessment criteria was RECIST1.1, whereas we included patients with BCLC-stage B and C and used the mRECIST criteria.

Some studies suggested31,32 that TACE is a good option for unresectable HCC before moving to the next step of treatment, such as surgical resection and liver transplantation. Because TACE can reduce tumor volume, control cancer thrombus as well as improve tumor microenvironment to some extent, providing transformation opportunities and enhancing overall survival.33,34 Furthermore, TACE treated before RFA can provide better clinical benefit and allow RFA to cause a larger area of necrosis due to embolization of the blood vessels surrounding the tumor to reduce arterial and venous blood flow.35

There are some limitations in our study. Firstly, only CT images were analyzed, and iodine oil deposition after TACE may be able to interfere with the assessment of post-response. However, CT images are more commonly and inexpensively compared with MR images and can be provided sufficient samples for analysis. Secondly, in terms of treatment modality, we did not distinguish whether conventional TACE and drug-loaded TACE would have an impact on survival due to the small sample size. However, the two therapies had similar results in the treatment of HCC and can be used interchangeably.36 Thirdly, although many studies (including ours) have now shown that radiomic models can predict survival very well in patients with HCC after TACE, these models are not completely consistent, which may be related to the standardization of imaging data, the way of features extraction, and the underlying conditions of patients and ethnic differences. There is still no better explanation for that, and comprehensive analyses such as systematic evaluation or meta-analysis are urgently needed to explain the reasons and try to find the best model. Finally, due to the small amount of data meeting the inclusion criteria, this study was no external validation cohort and was a single-center study. This will be addressed in the future with further expansion of data collection.

Conclusion

Noninvasive imaging methods offer significant advantages in predicting survival of patients with HCC after initial TACE. By combining clinical characteristics including pre-TACE and post-TACE, we constructed a novel, reliable model to predict outcome, including radscore and post-response.